simple2/person_network.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "slide_end",
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "# 参加者間のネットワークについて"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "subslide_end",
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "例2の結論でも述べたように、すべての発言者が$[0,1]$の間の一様な乱数を発言として発するという設定では、参加者の違いや距離などは全く寄与しない。\n",
    "したがって、次に考えることは個人間のリンクにおいて成り立つ関係はどのようになっているかということである。\n",
    "ここでも、まずは先行研究としてなされていることも多いと思われるので、そのまとめと一般化、それから簡単なオリジナルの場合を考えることにする。\n",
    "最終的な結論に行き着くために、これら個人間のリンクについての性質をまとめておき、その後にそれを反映した評価関数や閾値などを設定することにする。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "slide_end",
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "時刻$k$に$i$が発言を行い、その後時刻$k+1$に$j$が発言$x_{k+1}^{j}$を行う確率$p_{k}(i,j)$を\n",
    "\n",
    "$$p_{k}(i,j) = \\frac{g_{k}(d_{ij})}{\\sum_{j} g_{k}(d_{ij})}$$\n",
    "\n",
    "として考えていた。\n",
    "\n",
    "これは適応度モデルと違い、(初期のやりとりを新たなノードの追加と見なさなければ)新たなノードを追加することはない。\n",
    "また、はじめにスケールフリー性と呼ばれる性質を導くネットワークを作るものとなったWSモデルとも、そのモデルでは枝のつなぎかえの操作であるのに対して、こちらは最初からリンクの数を増やしていく点で異なる。修正WSモデルでは、少数の新たな枝を追加するようにしているが、そのオーダーなどから考えても、現在考えているモデルとは異なることに留意すべきである。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "subslide_end",
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "~~さて、現段階(2014/12/09)で良さそうな先行研究なども見つかりそうになかったので~~\n",
    "(複雑ネットワークは基本的にノードを増やしていく操作を行っている。これはむしろ意見の作るネットワークに対して適用すべき事柄である。)\n",
    "、上のような場合を考えるために、設定を別の視点から見直して、シミュレーションしやすい設定に変える。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_type": "subslide"
    },
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "$n$人の参加者がいるときに、発言者$i$を、アルファベット$A=\\{0, \\dots , n-1\\}$の一点と考え、確率変数$x\\in A$で表すことにする。次にこうしてラベリングされた参加者を円周上に並べる(図)。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end"
    },
    "slide_helper": "subslide_end",
    "slideshow": {
     "slide_type": "-"
    }
   },
   "source": [
    "<img src=circle_walk1.jpg width=200 />"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "slide_end",
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "このとき点$i$から時計回りに$a$目盛進む確率を$p_{i}(a)$と書くことにすると、進む目盛を確率変数$X$で表した時、$X$に対する確率密度関数は$p_{i}(X)$と書ける。\n",
    "また、ここで$X\\in \\{ 0, \\dots , n\\} = A$である。\n",
    "\n",
    "ここで、時刻$k$の粒子の位置を$x_{k}$として、その確率過程を式で表すと\n",
    "\n",
    "$$x_{k+1} = x_{k} + X \\mod n$$\n",
    "\n",
    "となる。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true,
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "slide_end",
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import operator\n",
    "\n",
    "def accumulate(iterable, func=operator.add):\n",
    "    \"\"\"Return running totals\n",
    "    \n",
    "    Usage:\n",
    "    accumulate([1,2,3,4,5]) --> 1 3 6 10 15\n",
    "    accumulate([1,2,3,4,5], operator.mul) --> 1 2 6 24 120\n",
    "    \"\"\"\n",
    "    it = iter(iterable)\n",
    "    total = next(it)\n",
    "    yield total\n",
    "    for element in it:\n",
    "        total = func(total, element)\n",
    "        yield total\n",
    "\n",
    "def circle_walk(P, n=6, i=0, K=100):\n",
    "    '''\n",
    "    alternative view of malkov chain\n",
    "    \n",
    "    [Keyword arguments] = [default value]:\n",
    "    n = 6   -- a number of participants\n",
    "    i = 0   -- initial condition\n",
    "    K = 100 -- trial\n",
    "    '''\n",
    "    k = 0\n",
    "    speaker = [] # time series of participants\n",
    "    speaker.append(i)\n",
    "\n",
    "    def gonext(i):\n",
    "        w = list(accumulate(P[i])) # create weight list\n",
    "        rn = np.random.rand()\n",
    "        for j, _w in enumerate(w):\n",
    "            if rn <= _w:\n",
    "                return j\n",
    "\n",
    "    while True:\n",
    "        if k+1 > K:\n",
    "            break\n",
    "        k += 1\n",
    "        i = (i + gonext(i))%n\n",
    "        speaker.append(i)\n",
    "    return speaker"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false,
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "slide_end",
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAERCAYAAACU1LsdAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXuUZVdd5z9VXV3V1c+8Oo/utKZJwktWeEVEECkgIDpR\nwKxRQIUhS2dhYABFhLBAOktnzeiooI6oIIaIgGaCwwTkIRkoYBQjhICQh4QkkKS6093pdHdVd996\n9p0/du26u3btvc/e531v/T5r1aqqe8/ZZ599ztm//f39fnsfEARBEARBEARBEARBEARBEARBEARB\nEARBEARBEARBEIS+5gzgJuAu4E7gmc1WRxAEQWiSG4Crl/8eAXY0WBdBEAShQXYA9zVdCUEQBKHH\ncIPH3gscBq4Hvg68H9jcYH0EQRDWPU0ahRHgacB7l3+fBN7WYH0EQRDWPSMNHvuh5Z+vLv9/E5ZR\nuPjii7v33ntv3fUSBEHod+4FLsmzY5NK4WHgQeCxy/9fAdxhbnDvvffS7Xblp9vlXe96V+N1aMuP\ntIW0hbRF+Ae4OG/H3KRSAPgvwIeBUZRle02z1REEQVjfNG0Uvgn8cMN1EARBEJZp0n0kJDAxMdF0\nFVqDtEUPaYse0hblMNR0BTLoLvvHBEEQhEiGhoYgZ/8uSkEQBEFYQYyCIAiCsIIYBUEQBGEFMQqC\nIAjCCmIUBEEQhBXEKAiCIAgriFEQBEEQVhCjIAiCIKwgRkEQBEFYQYyCIAiCsIIYBUEQBGEFMQqC\nIAjCCmIUBEEQhBXEKAiCIAgriFEQBEEQVhCjIAiCIKwgRkEQBEFYQYyCIAiCsIIYBUEQBGEFMQqC\nIAjCCmIUBEEQhBVGGj7+94BpYAlYAJ7RaG0EQRDWOU0bhS4wATzacD0EQRAE2uE+Gmq6AoIgCIKi\naaPQBW4Bvgb8SsN1qZ1uFxYXez/dbtM1EgRhvdO0++jZwAFgJ/A54G7gy+YG+/btW/l7YmKCiYmJ\n2ir3yCPw0z8NX/lK2n6vfCW8/vXwrGeFt7v6arjhBhgehtOn4TWvgQ98ILzPP/0TvO99ar82cuQI\nvPjF8NWvVn+sW26Bm26CP//z6o9VJ1ddBe96F1x2WdM1aRd/+Zdw6BC8/e1N16R9TE5OMjk5WUpZ\nTRuFA8u/DwP/GxVo9hqFujl2DO66K32/++5T+2UZhcOH4eab4cor4dOfhve8J7vsBx5QP23l2DG4\n/XZYWoING6o91v79qq0HjYcegocfFqNg8/3vw8GDTdeindgD5uuuuy53WU26jzYD25b/3gK8CPhW\nc9VZy/w8zMyku3VmZlSHFVP+xo3q7+3b1X4xZXc6afWpk/l5ZRDqeHg7nbg26zfm59t9jZtiZmYw\nr3fbaNIonIdSBd8AbgU+Cfxjg/VZw8KCcuucOpW23/Q0TE3FlT86qv7etk3tF1N2mzuMhQX1O8Yo\nFqXTiWuzfmNhod3XuCmmpwfzereNJt1H9wNPafD4megObmYGtmyJ329mJt4oaKWwbdtgKAXdZlNT\ncPnl1R5rUJWCGAU3ohTqoenso1YzP69+p4xOut14ozCo7iOIO/+iDKpREPeRG1EK9SBGIYCpFGI5\ndUq5nPK6j7LiF/3iPqrLKMS0Wb8hSsGNKIV6EKMQQI96U27EmRk45xw4erS3f6h8rRRGR2FkBGZn\ns8tvc4dRt1I4fbrd7ZEHUQpuxCjUgxiFAHrUmyJZp6fhjDPg3HPhwIHwtqZSgLhgs1YKbR0dLyyo\nVNS6jAIMnktBlIIbcR/VgxiFAHncRzMzKj6we3d2Bo4ZaIa4uIJOkc1SIU2xsAB79tRrFAZt9ChG\nwc3MDMzNtffeHxTEKATIE2ienlYj/t27sztG030E8UoB2ttpzM/D3r2iFIog7qO1dLvqOm/dOniD\ngLYhRiFAUaWQ1THa7qNYpQDt7TQWFpTrbH4eTpyo9liDqBS6XVEKLmZnVcztrLMG63q3ETEKAfLG\nFIoohRijsGVLezsNbehi3GdF6XRUWwxSJ7G0pH639fo2hX6uYid5CvkRoxAgb/ZRilLI4z4677z2\ndhra0MWcf1E6HdUWg9RJ6Huurde3KfRzFTufR8iPGIUAetRblVJIdR/NzcHQEOzY0d5Ow1QKdRiF\nc88drE5Cq9O2Xt+mmJkRpVAXYhQCzM/D2WenK4UYo6B9xylKQRuc8fH2dhqiFIohSsHN9LQohboQ\noxBgYSGfUdi+HXbtUj5133yCxUUVOBsy3juXdcPrsttsFLShE6WQD1EKbkylMEjXu42IUQiwsKCy\nHfK4j7ZtU5O4jh/3l22qBBgMpaDdR7t2iVLIgxgFN6ZSGKTr3UbEKATI6z7avl39HRot25lHEKcU\n2m4U6nYfDZpSmJ9X6rGt17cpRCnUhxiFANp9lEcpQLhjtIPMkH3D94v7qI5Ac7erAu+DZhQWFtQ1\nbuv1bQpRCvUhRiFAkZgC5FMKg+A+2rhRuY8OHuzl3ZfN7KwyPjt2DFYnIUbBjSiF+hCjEMB0H8Uu\nQJeiFFwxhVilkLWaalNoYzc6qhYGPHSomuPMzqp2GLROYn5eGbq2Xt+mkMlr9SFGIcDCgup4Nm6M\nH7nFKgWf+2gQlII+rypnNXc6PaMwSJ3EwoJa32dhoTqV1Y/I5LX6EKMQYH5edXApHY+tFHydYt5A\nc9tjCuZ5VRlX0EZh0DoJfc9t2iRqwUQmr9WHGIUA2sUT2/F0u2oRuCKB5kFQCnUahUHrJHT7tfka\nN4FMXquPkaYr0Gb0Axrrt+50em9Qg3Cufiim0O2untSm0aOlNq+iaRq7KucqmEYh1Gb9hm4/MQqr\nkUBzfYhSCKClfGwanOk6AjWx6siR3oQku2zbKIyNwfCwSrV0Ie6jHp2OcrHEvsa0X9Dt1+Zr3ASS\nklofYhQCpCoFM8gMqrPauRMefthdtu0+gvBN3y/uIzPQXLVSgMHqKMR95EaUQn20wShsAG4HPtF0\nRWx0B5dXKYC/Y3QpBQjf9P2gFOqOKcBgdRTiPnJjDogWFtzqWyiHNhiFNwJ3Aq17Fb3uuPMqBfB3\njK6YAvS/Uqg7+wgGSymI+2gt3W5PKQwNDdYgoI00bRQuBH4K+EugdWFC031UtlLwuY8GQSno8zrr\nLBUfOXmy/OMMslIQo7Ca2VkVaxsbU/8PWsZZ22jaKLwbeAtwuuF6ODEDzXans7i4NiDsUwquuQoh\n91GWUti0KbvD6HbjO5VTp+K2i8E8r6Gh3hLiRXDVzzYKdpvNzcVN/ordznV8FyHXRsz10PdcUaMQ\nu2/Z21VRhv1cuZ7HTidu1YHQdkXbotuF738f7r9f/Rw4EFeeyfy86ltSKXMA0aRRuBI4hIoneFXC\nvn37Vn4mJyfrqhsQDjS///3wlres/kxLXJPzzlNrALnK9gWaQ0oh1n30L/8CL31peBuA22+HK67I\n3i4W2y12/vnu84/l5En4gR9Y+yDb7iO7zV77WvjYx7LLf/Ob4YYb0ur07W/DC17g/u4P/gB+//fd\n3z3veXDnneGyy1AKp07Bnj1x2/7gD/Ze7BPisY/1LwMfw513wvOfn29f+7lyPY8/9VNw663ZZb3i\nFeDqRhYWVJtlGZZuV92Prmy3z38envAEdZ7Pf77aLrXNfvu34U//NG2fBx6AH/qhyVV9ZRGanKfw\nLOBnUO6jTcB24K+BV5kbFT3BIpiT1+yR6L33wn33rf7M5T7yrVOUGmiem1M35NhYXIdx6BAcOxbe\nBtR5HD6cvV0strHbvLmYEnnoIZXWOzenFJImy310773x5793b1qdQm179Ki/k425JmUEmnWb+eJW\nmqUlde314oIhDh9W57ZjR746HTuWfx0snY6q8T2PMfex7xrMzvbaLNQWs7PwyCPq+Ob9COrzK6+E\nG29U/194obovU9rs0KHePKdYjh2DubkJ9u2bWPnsuuuuSyvEoEml8HZgD7AXeDnweSyD0DTmMhd2\npzM1tTZW4HIfjY66O4nUQLMZaIvpMKan40aA+/eX65+1jV1RN4h2Pdl1zAo0T03Fnf/UVHr9Qm07\nPx92/2XVqYxAs74vs/bXdcmqU7cbPq8YiuyfpRROn1aumti4n+t8Y1+Dqo/hOpY9KMxzDWOfW5P5\n+cFxH9m0LvsopBT271/rK3cphbEx92S01ECzaXBibraZGf8kOJOpqXKDtPZ5FTUKuoNz+ZB9SqHb\nVfvFnP/+/en1C7Xt3FzY/ZdVpzLcR/q+zNpf1yWrTouLvQygvITaJQv7ubJjSIcPqzrGZgi6zld/\nFvNcmb/t78xBYZ5rGPvcmszNDaZR+CLKldQqdAfnUwqHDq226i6l4DMKPveRTymYD0asUog1Cp1O\nvuCWC1sBlWUUUpTC0aPq3LPOX7sM8ozmQkbBdf3m5tQ1jzEKRd1HsUoh1ijo74soBX09UkfBkB1o\n9t0jLnzXLtYo1KEU8hiF2dn45f2zaItRaCXmPAXzJtAj0TPPXJ1h4FIKo6N+pZASUzAfjE2bVN1C\nWTMzM/HuE719GZTtPsqjFPQ+WecfO6K2CbXt/Lz/+sXUqY3uI/19kXukSBlZSiH2HtZqp4j7KEsp\nmPWMyRJ0lZHHfQTlLfUiRiGAb5XUo0eVArj00tVxBZ9S8MUUUpa5MB+MoaHspZVT3Ee6/DKoyn2U\nohT0PlnnH9t52mQphZBRqMN91FalkLeMspRCp6MGUmUoharcR3mVApTnQhKjEMD3PoWpKTX/wJ6Y\nlhJTSM0+skchWTdczM2lFc+ePe1WCq76xSiFGKMwPJzf7+uS6z73kf4s1n2UZ5SpiT2vVKNQNKaQ\nt4wYpRBzD4cMc2pMoSr3Ud6YAohRqAV7noLuBHxGwTVPoYpAM2TfcDE31/HjqvPYvbs8o1BFTOHx\nj083Cjt2xBmFiy7K9+B2u+44TFGlUJb7KOa8+sUoxCgF1z3iKsesi6t+/RpoBjEKtaA7OHt55pBR\niE1JLRJohjilMD8fDj7p8yhz2YAy3UdLS2ri2+Mel+4+esxjsn2zU1NwySX5JD74r+v09Np2D+1j\nUjTQfPq0WpX3MY8pP6ZQNCU1bxlZKalTU2rSWFbZWdcNmg00LyyoPiZvTEGMQsXokaDuuM2OZ//+\nnlEw01JTU1LzBpohTin4RrMafR5lvc2q2y1XKRw+DGecAeeck64UHvOYOKWQxyhkjTh9S6D49jEp\nqhQOHVIq6YwzBkcp2M+VPQjYv79epbBjRzVKIfYesRGlUBMLC0od6Ld5mR2PSymYKzmapLqPylQK\nEL7BylYKi4ur2yymniFC9ctSCnv3xhmFiy/OrxRCnYtd35SYQhGjoA19zP79GmguqhSKBpp3765G\nKcTeIzZiFGrC7rTNjmlqSi30Zr5ustNRD7M9+velpGYFmm33Qx6lAHFGoSyl4FI/RQOmvvq5lIIZ\n86nSfRRKL/WlXsampBZ1H+k2i2n3OlNSqwo0dzpqGZWLLiqmFFJSUn0xuLKUgriPWordaZsdk60U\ntEqw4wkQTkl1GYWxMTXStm/c1FHIzIz/2JqylYLL0NWhFEZHYcOG3gSpY8fUujMhg9jtqjkmeZWC\nTwHOzanvXErBt49JUfeRHrCUrRRc55TC/Hz+MkKBZn2+Me/U0NfA9UykKIVdu6pTCjH3iI0ohZqw\nO21bKejOamREdUKueAL0LrI98g8tvOWKK6SkpGojdc458e6jspSCfU5lGIUspQC963PggFqZdnw8\nfO5HjqjF+s48M59S8LXt3Jw7BhJzPaC4+0i3WdlGwXVOKRQpI6QU7Hs4lFiRdd0gv1LQk0nNRfLy\nKIWYe8RGjEJN6DkKGt0xzc+ryWvnnqs+12rBpxSGh9Uo1g74+txH5rFMXNLUN3nt1ClljDZvrtd9\nVLdSMB9AfQ7ap+5z29llp9av24UTJ+Dss6szCtp9lGeGapVGoWhMIW8Z9r2/ebO61xYWeuerXbdZ\ncbaqjIKuY5F4mr5HUt1HYhRqwqcU9Eh0wwb1uTYKPqUA7rTU0LLGrk4wRZrqbWPcR7t2lec+cp1T\nUaOg62c+hIuLalTmivnoTiLWdZZav5MnlTHSHZPN/Dzs3Om+fjt3Vr/MRcp5pcQUdu4svsxFnjJc\nCRzmKzn1+UK24p2Z8V+DlJRUV6DZ9fzncR/t3JmuFCSmUBO2K0SPRM2bEFYrBZ9RcPkJbSVi4nMf\nxQax9LYh/+Tiolr//fzzyw00l+k+8qXMzs6qcs1Rmd1JZPlmze20kYkhq23LVAriPurVTb+KU+N6\nHrPu4yz30ZYtxZWCSZ3uo5i6xyJGwYPtCtGdju6oNHqugs99BO4OJKQUXEGzvErBd4M9/LAalYyM\n9F+g2Y4nQK/NUo2CXkcqto5Zbetzk4RcFya2Ukhd+bIqo7BjR/hVo1nkdR/5nivzedy1q/dZqPws\n91HM3A4z0Gxem7KUQl6jEFP3WMQoeLA7bbvT0cS6j1xGoWqlEPKrp4ywYilTKZw6pRTBWWetDSK6\njIKtFHwzyTXaNZVaR7Nt7fJPn1aK46yzigeaR0ZUPCqlEz51Sp3H2WfHGwVXpptru7GxYgkJeZWC\n77my3YVQTCnMz2d3rGYCx/Dw6nLKVAp5YgpnnCGrpFaOa55ClvsopBTsC50VaDZHPDqzwZTQWUYh\nK6Zg+2LbFlPQnfbQ0No0XZ9SyOM+0nWMfaBCSkG7BH1KL+aBN++71LbTo2b9dr6sc5qfh61b42IK\nY2NxaZ+hMlyxlix8z5XLfRSjFHwxBa2GQm02O6uMgesal60UUhTi/LyquyiFinG5j/TIRI8woTeB\nLaQUUt1H9ohMd/KxmQ0x7qMqlEKZ7iPb+Jp19CmFPO6j1DqaBtcuPzSijlUKZhumtl2q+pmbU+1a\nl1I4++xylcLx4yrxQ59zjFIIZY1lKQXTQLmeUZdSSBm9z8yoOgwPp730StxHNeFyHxVRCinuI/vm\nTpWmqe6jMpWCfU6pgVxX/ew6ZsUUdu3KNgpmbCil852e9rft3FxvFFnUfZRaL0g3dClGwaeAYik7\nprB9O9x3n7ovdGpy1TEF00D5Bm4meZSCzzUZQoxCTdjZQS4fJqj01CNH1E9KSmrIfWTf3KnSNCYl\n1TwPPZ+h6Cs5XeeUGsh11Q/ilMIDD6i23ro1/GDNzfVcCZBfKbiuqR5Ru9x/O3ZU6z5KNQrz86qu\nse6jIkpB++wXF9M6vJBSuPtu/z3iQqek+mIKWS4Y00BV4T4KqdAQ4j6qidhA88iImsj2ne+Um31U\nhlII3VzmSFnnfZ844d42Ft855ZmIZWd5xSiFu+7q7aONgss3u3+/SsUdHu7VLzXQHHIf+a5fyjIX\nqfWCfEph27Z491FRpeBqmyx8qd729YawUuh2w3GdWPdRSCmUEWiOvU9S656CGAUPrkDzgw+qB3br\n1tXb7t6tRi1lzlMoSymE3EdmbKQMF5LPJZYnruCqX5ZSMEeOw8PqWrk6ANuwp7qPYmIKrusX87DX\nqRTqjinkKSOUkmorhVDZs7NqALd1azH3UYpSSFXIKfdJat1TEKPgwbUgnp7NaKNnOFYZaLYfjNAN\nlxpTgHKCzT6XWF6jkCfQbO7jcyHZBqeslFQz+8jlb85aekOX0UaloF9LW8Qo+OItIbJSUmPv4axn\nIiYltS6lIDGFluJa5gL8RgH87iPfMhehQHMRpZCVkqon3thr1LdNKaS6j2D1Pr4RVxNKQbsFQg97\nt6tiDyMj6fWC6mMKRVNSXW2TRSjQDPHuI33dfG5FnZIaG2iOeUZ1/xE71ySvUhikmMIm4FbgG8Cd\nwH9rsC5rsDu40VH1EzIKqe6j1JRUkyz3Ucg3ac7m1ZShFEIxhZQbVr9S0hzNxygFqN4o5I0pxDzs\nuv30dUmpl91mbXUflakUIF0paLei3VFXkZIK8ddQzxYfH1/f7qNZ4HnAU4DLlv/+sQbrswpXp719\nez6l0ESgOeSusDtFKEcplOU+Ony41/G66tdmpTA62qu3/j42gGi3X0q9HnlkdXrmpk3Zk6BSjUIZ\ngeaylUJsXMw0Lr5rt327uganT2eXETNwg/hraM5FylKUNnUbhWHg58o5lJNTy79HgQ3AoxUeKwlX\np71t2+qbUGOuveLC1Tk3GWh2GYWylEIZ7qOs+oWUgnl9QjGFIkrBdEOYaBeJro++hrbrwod9z6UE\nKu04iZ55G8r6SnEflRVTSC2jTKWg9/E9j5s2qevna7OslNQiSiGrfiHqdh+dBt5azqG8x/8GcBD4\nAsqNtIalJZVGWBUPPbR2ROXqtENKYdMm/8jftvynT6sfvfy2jb7hvvQl9XPPPflTUmM6RUgfwek3\nzpmUpRSy6td29xGs7qD0PiMjvfWRXNhGNaVerjbL2r+KlNSHHlr7mT7fkZF0tRFSCqOjKsVUk6UU\ndDmu50IbrZjBlj5WmUrBrl+T7qORiG0+B/wG8HfASePzMkb1p1Huox3AZ4EJYNLcYN++fXzve/DF\nL8L1108wMTFRwmFX8xM/AR/9KFx2We8zl1K46ip4+tPX7r93L/zKr/jLty+y7Tt2bX/llfCOd/Q+\ne+pTV29TRCkcPgx79qz+LHUEd+WV8MEPwpOf3PusrJjC4cO9lxi56ucyCmNj8NrXqsmE5me+8zfL\nHx9XL06KISbQrOtrKwXtGpibUxMGbYq4j/IahRT3EWTfI5dfDrfdtroudruk3GdHj6oOz2b3brjm\nGvfy6S7MTjt07WIGWxAXaIZ8SiHFKHS7attvfnOSEycm2bcvbr8QMUbh5UAXeJ31+d7ih1/hOPAP\nwOU4jMInPwl33AEV2ANAzYS1byaXK+Sd73Tvv2kT/PEf+8u3L3LIdaT5+78Pf28urWw+GPrNYKGY\nQqeztlPavh0eTTDzjz4a12ZmXWPx1S+kFAD+7M9W/+9z19jlj4/HK9FQaqMeber6mkpBGyu9n8so\ntFkp6ESL0Ch/dhYOHlw7CdJul4MHw8czsd1imi1b4N3vXv2ZNsT2MwGrR+K+axdjFFxKYX5ezdR2\n3ZMpRsGsX2xMYXFRuQpf9KIJhoYmeMc7lCK77rrr4gpwEBNovghlAOyfopwD6DHAOPBC4HbXhp1O\nedLIZnpa3cR2+aHsoFTsixwKMseilYadRdHp9F5N6HMfhfL8Y5mZiW+zPEbBtzS273sXvgA/5BuR\n66WTt24NL3Oh62sqhZDrwty/TqWQJyU1NMrXhtV1X+RRCidPqut31llx24+Oqg7RFRPIGonrOobi\nOL7Ja64FKzUp7qM8SkHXW6+MW0Y/GWMUtgDvBN6//P+lwJXFD80FwOdRMYVbgU8A/9e1YZVGYWqq\ndwyTMjpujct9lKUUYnDdBDG+SV/2TuzDqpcMcLVZWUohVL9OZ/X7mX24zt9Vdmz9Tp1Sxx0ZSY8p\nxDzw9j3XNvdRVofue5Z87ZKFuXx6LL7yzZF43piCTymEFsPMoxRSjIKpwsoyCjHuo+uB24BnLf+/\nH7gJ+GTBY38LeFrMhlUaBd/opqyOG9zuozIMjr4JduzofRaTxVBUKczOquBhrCHNYxTspUSyAs0u\nyjYKMWmN+p4JKYWQUcjrPrLXioLs7KWyU1J9z5Iv1pKFy9Blocu3Y1Ix1y7LfWQHmu2YkYu891as\n+8hs2zqVwsXA7wK6micD21ZCE0qhbPeRK9BclKaUgtk5m1TpPspKSXXh8s0WMQpZfl/TTeJTCqF0\nwyYCzfolO1nzGfQoOvRKzpBScMVasshjFGKUgn0Nut3etUsJNJetFPKkpJpGIc9qxC5ijMIcyuev\nuXj5s9oYBPeR2YHEBJpjcN1wtquiipiC2TmbVOk+6hel4PKdZ7kuNHndR52Oioudffbqz2NiCtod\nFlqKwfRbh1xIoQFW3UrBJjQSX1xU6eHDw/420y5TswxQ19yXjgrVp6SabVunUtgHfAa4EPgIKg5Q\n5dyFNXQ6qpF8Mw2LMDWlHkTXjVyV+6hqpZDXfdRmpaA7pbm55oxCTFqjz30UG1PI4z7avx8uuKC3\nFHjM/nq+hE5ICHVCsSmlvmcpb0qqyyWWRUgp+K5BjF9+bk61rznLXrvTfBPXQuWl1C9EFTGFGKNw\nG3AV8BqUUbgc+F7xQ8ejT7SsF1ObTE2peQZ1B5qrMgoxASufe6YMpVCVUYBeR5tiFKp0H+UJNFfh\nPvKNqEP7687EnDvhQufBu4ydqx6uZ8l2HzWlFHzXLsYv74obaANXtlJISUltKqbwCWABFVj+JLCT\n4kHmJPSJVuFCmpqCSy6pNtBsX+Q63UdVKIUm3EdmHVNiClW6j8pOSc2rFPIYBfMeDHVCS0vKcOiV\nW0MxAd+zlFcp5DUKqUohxgXjihuYSqGoUcirFJpyH/1XlGHYCjwd+F/ALxQ/dDxNGIUyA81VKgVb\nPcV0QK6Uzi1bellFWdThPnKlnOqHsKj7yC67rLRBl1LQE5v0MUMPvEspxKjjvEpB1zVUJ3M78I/E\nu13l7rn44nBMIStYbeKbuBbCp0RCz0WMC8alBkylsN5SUv8BtWDd51CG4WeBfy9+6HiqMgqLi2rJ\ng7171fpCJv06TyFvSurQkMpEmZlxLytgH0OXY1KXUpidbW+g2Xaz6Idd59pXMU+hTqPgUwpHjqjj\nnX12WCmYwerQpDTX8ukxuJSCnnQYiinkcR+ZSsFOgdXUnZJahos9pBT+xPh5PrAduB94PRBY1KF8\nqjIKhw6pm3jbturdR1XOUzDJm5IK8emC09Nw5pn1xxRSlUIbUlLtEWYoplAk0FylUbBfS+u6R3Qd\nXMe0y4i5zw4dUvNvTIMUg6vsuTmVXWS6y0LuI9+MaFsNVKUU8qak1qEUbkOteQQwZPw/ZHxeC1UZ\nBS1PXfm9ZbuPyl7mAvKlpC4uqlGY6/ix6YIzM2pkFNtmqfnTIaVw5Ijyb/tWmDUZG1u7npOrbF0/\n15o5JtPTsHNnr+yYlFS7w6himYuQUjh82H+smJiC2WGC3z2T9SzFuKDs8lLjCb6y7VF+llI4dGht\nuXUEmtsSUwgZhQ8WL74cOh31sFZhFHyjmyrdR1UGmrMW/tKdoqvzS1EKLqNQtftIL6YWoxIg3n2k\njUyWOjRT7fqXAAAgAElEQVQ7eB14XVzs/e3KssnqkEyKBJpdbpaQO6FITMF1j4SepVgXlEmedFRf\n2VmGOTamEAo0+5RC7ICo31JSH4ta1uJOlPvofuC+4oeOp9Mp9yUSmtCN3K/zFFw3lzlbNeR6SVEK\n551X7zwFXb9Dh9KMQoz7KLaO9ogwNOL0jSJTYgqbNoXfBAa9AG/TgeYUo9B2pVB3SurioqqDXjm3\nH5a5uB74c2AR9b6DG4APFz90PJ2O24ddlLqUQhWrpEL2PAU9S3Nxsfd9llEoqhR8RiElAFaWUohN\nSdV1zLq/7BGhK9XY7Gi7XfWazNA+JrZR1XMIQm33yCMqcyz1nFLcRzHxgCyjkBpTyGsUYpRCVkyh\nTKUQc+/bq6ymxBSaSkkdB25BxRK+j5rh/B+KHzqeTkdlKtRtFKpcEK9K95E9KjIf+JBRiHUf+ZRC\nHdlHqUqhTKOQohR0ls3UVH73UUy9Qp1nW5RCG2MKZaekFlEKLvdW291Hs6j3J38XlXn0s6jltGuj\nCaPQL/MUQu4jWDvqKMt9lBJo1vUMLbqmyQqEVxFTMOsYIsYNYWfqTE3FP/Cu9mubUcirFGLiEq7y\nUnGVnWLM86akFjEKWUYrRFPuozcBm4E3oJa4+EXg1cUPHc/sbDNKod8Dza5jl6EUpqfTlMLIiHJj\nxUxYygqEpyiF2JRUyDeiC7mPdH33709LSW2jUbANXZ6Ygu0+qkopuMp2uf3KcB+VlZKa5d4K0ZRR\n+FdgBjgC/CeUUviX4oeOpyqlYOZW236/Kpe5qEop2JN0IM0opCoFV5v5ziv2hs2qX5Puo5SApa6v\nSynELnMRU69Qlk5ZMYUspTA3B8ePq3uiTKWQOnHNLNtUpVUphW3bVMrzwoL/nizjvgrRVEzhWajM\nIz2L+cnAe4sfOo5uV3U+Z55Z7oJ4J04oV8WOHdW7jzZsUCNfHfAtyyjY6W6dTu9VnBr7gS9LKaS4\nj6Aco7B9u/9duC7KNArmu6995btGxCkxhX5wH7k69AMH4Pzz/UtPxxgWk1OnVBn2UuAxjI6q5808\nnypTUrUS9M1vGdSYwnuAFwOPLP//TeC5xQ8dx9ycelC2bClXKZiv+qvafQSrL3RV7iNXwKtspaDV\nSMo8BVddfWTVT5cVQ5kpqadO9d4BbJafpRQOHCgWaM7Kcw+NqOsKNJt1cKnu1JTU/fvTX8NpYpdf\nZUqqfX1tBjWmAPCA9f+ic6sK0A9xWSesMUdYmzatzecvc54CrPYTVuU+co1myo4pzM6qkdj27Wtz\n6OtQCrqsGMpMSXW1bUxMYXExPiW1zpiCnT5bJCXVrEMZKal54wm+8qtMSbWvr42rb3GVHRt3smnK\nffQA8Ozlv0eB3wDuKn7oOOowCsPD6kKYI5wqlIJ+8KoyCq7RjP3AhxaTi1EKZk61nUNfR0wB3Cuo\nunCNgH3nH2MUslSYa0Rs/vbVSVO3+8i3FpC9nXlOmzerz8y5L1lGITUltahRyKMUslwwPqVg/nbh\n6ltcZfeb++hXgdcBu4Ep4KnL/9dCHUYB1pZfZqAZ6nMfuUazKe6jLKVg3sApbRZ7/UJGK1UplBlT\ncE1QMss/fXqtUTTf5xuqkyY10Dw7u3o9Jpsq3EeuV3LaRmF2dvXIOHWZi6qVgiumEKMUXIOCjRvD\nSiFUpq/sfnAfHQZeCZyLesHOL6AykWqhCaOwtKRu6phF12IxH7y6lUKZ7iPzBjaPr1/v6GuzMpSC\nfgibSEl1dQqmwdWG3vSDu0aSKW9ey6qX7zWcmo0b1XVZdDh786akwtqRuK26R0b8gVzX/jZVK4XQ\nm9d0+5vp03pb07BBz0CGlAKkuyb7ISX1YtRLdh5BGYj/Azym+KHjaMIolO06gtUdVD8Hmn1GQbdZ\nkSyMrPrpOjalFEKuOdtFousKaSmpKUYhq/P0JVHo+uZJSYW1g4cs1Z2aklqGUbCVQsh9ZJ+jb7Dl\nurfLMAr2vWUutphFUzGFjwA3AhcAu1BvXvto8UOzB/gCcAfwbdTkuDVUZRTs/G67gyvTdQT1KYWs\nYGhRpeBzH2W1WVlGYfv2ZoxCVhDfNaLW21e1zEXMSqK+7KW87iNYO3gIPUsQv/y2r7xU7PJdPvuQ\nkrHbLDQ5bfv2ctxHWQkiPpqKKYwDH0K9p3kB+BsgMtQXZAH4NeCHgGei4hRPsDdqQimUOUdBU4dR\nKKoUdNpv6JWcPqWQ1WaDqhRMo+DqPHXufqhOmlT3UcyI2rd/UaOgBw/d7tq02CylMD6ungHfSLhu\npeCqX5Zb1jxW2UpB1zEmrmC7j+x4Th5ijMKngWuBi5Z/3rr82VnLP3l5GPjG8t8nUBlNazKuqzAK\nS0tqDZ3zz+99Vof7qOxAs05302mhRVNS9Ss5T5zwH7OflIJWSfoh6XZVW7iyl4qmpPrcLOarOPU+\nZSmFmFm/ZRgFlwLSne7Ro+r7rVv9x7TLcAWrNadPq9z/Cy4In1cIs36uWf6hlFRX/ZtQCrFxBfM6\nDg+rfis2HuEj5h3NP49609p/9nxeRnzhIlRW0632F2UZhQcegC9/Wf09M6NmSPteaFL2HAWoJiV1\neFiVe8MNqr5f+xq88IVrjxtrFKDnGtixw/19PykFHfTUne3srPrtCswWTUn1jahTRoC+mMKXvgQf\ndixWf+ut8Ku/6q+z3r/smMK2bfC5z6lzdrl6spSCLuNDH1o7a/nkSfVdbNqxi23b4JZbVJstLvbS\npzVZM9HbohRiOne7rypj8BxjFH4T+AwwDfwWqvP+HdTrOctgK+olPm9EKYZV/O3f7uP734e/+As4\nenQC9UqHdP7oj+CLX4THPU79/6Y3rf6+aqVQhfsI4A1vUA8AKEP3XGuueUpMAdQaU0eOwJ497u9N\npWD6XrPOqSyj8OpXww//cHY5Gt3uo6PhsmMe3L173WWDe0R96aXw5jf797FxGdYf/3Fl7D/1qbXb\n790Lz3722s9NqnAfvexlcOONvTrZhikrpqD3uXXNEFDxa7/m/jyWiQm47bZe/exnPZSSCmvrf+SI\nerZcvPKV6jqHyLq3XOWnuo8mJyeZnJxkcRF+53ey9wsRYxTeiQo0/xjwfOD3UWsf/UixQwOwEfgY\nKk7xcdcGL3rRPu64A669Fj7ykfwHeugh+PVfVxfRRZ2B5jKVyO/9XvxxIbvT3bVLjf6e8hT39zMz\nPYNRlfvonHP837/85dllmJgPVxGjEBNotju+rVvhjW/072PjasOnP92tEmKpwii85CXqJ/aYrjKu\nvTa77nl52tPCbWa6FYeGso1CKPD90pdm1yd0b3U6Sh3Z93xKoHlsDCYmJpiYmOCv/xquuQb+8A+v\ny97ZQ0xMQYcdrwTeD3wSNbO5KEPAB1CL7b3Ht1FZ7qOs4FXVgWZzxF6FEvGRahR271Zt5cOUuk24\nj1IxfbNFlUJqTMFXn9g3r5VBGe6j1AFMVkyhaYaHe+/khuyYQtHAd1YGmWudp9iYQlbd8xBjFKaA\n96FiCP+AyjyKXTMpxLNR72Z4HnD78s+L7Y2aMApVu4+qiFn4cLmPQv7aLKNgjphT2izFKBTxJ9uY\n7R4qu4yYQsw1zVIKdRmFIkoh9ZixBrNOQteubKMQWtTQV3aelFToZSAVIcZ99HOozvp/AMdQ8xXe\nUuywAPw/IoxLp6NSJc2ZhqkPTrerMhpSjEI/zFOIPa6Zsx2jFL72Nf/3oclrdWQfpWIbhSJKITXQ\nnFUfmyruuzYYhTxlVI3pVsxyH1WpFEJGITUlFbJX1Y0hZsR/EuX3v2f5/wPAPxY7bDzm6C6vWjhy\nRC3kFepszMasyn3UlFEo233kUgptdh9VFVPI4z4aGVEp0a65IHW6j4qkpGZhPktLSyrNtMwlY8rA\nThG3jYI52m7CKORJSc06VixluIEqxXyQ855waj53Ve6jspe5iD1umUZBlIK77NjRsE6PdI0C61QK\nRVJSU46p2yXvuxGqItZ9tLSk3vZXZN5EXqXQVErqujEKKTM/q5qn0IRSSE1JHTSlUJZRyLPMRUyd\nTAYxptDGeAKEDbpZ/4MH185pSqUq95GeEW6++GndGYW8/rJUo1BHoLmt7qNzzlEdvy9YlbUgno8m\njUKs+8h3zq5ZsbrsVKVg72dSlfvIdV51GYU2xhMgPqZQdB0muzwbnxcjRim42nZdGAVzff26lEIV\nMt5OSW3CfaTfdx3qdIeHlVTev3/td92uPyU1xn0UkxXRxpTUTmftqzh12akxBXs/k7a6j4qkpLYt\nHVWTFVPQ9S8aT7DLsykSU3Ddc+vCKDTlPuqXGc1ZmA/83Jw6rm/9fc3u3W6jMDfXW1oDBst9pEeO\n5utFNS6VYJct7iP3MfvFfeTzy1dpFHRWZBGlkLJWVixiFOiVrUeygzRPIbZTNPHFFexga0qbtd0o\nDA0p96RLzbgmrtllF3UfVfFiJ4g3CubCgb7t8hyz7e4jvUhiSClkJalk4VPJjzyiZr277smYmMK6\ndR81oRQGbZ5CWUbBDrb2Q/ZRbEoq+OvoUwpF3Ee2UajKpRhjFOwZvvZ2g+w+Wlrqnb+mLqUQKjvG\nfSRGgfyz9WKCRVW7j9qwzEUZRsFWCoPiPgJ/HX2rZBZRCvYosKrkA9c56Vd0msfzxRXKSkltG/ra\nZblgmjIKMe4jl8dhXRqF1BOem4Pjx/0vN9fUlX3U7dabfZQyUtaE3Ed5lUJM5lhMIDyVFKPgq6Nv\nPf0yYwp1KgXX+6R9ndCgxxSyRttNGgVxH3koahT271cv08kKrtY1T8ElV6ukLUph0yZ/IFcTGwhP\nITYlFepVCi6jUJdScNXVVSdXHnzqMdusFObnszvWKlNSiyoFMQrkO+FYS1+1UtAj9jpVAuQzCrt2\nxQeaY9tMz+QNuf/Kdh1BfEoqhGMKLqVQZkyhTveRS9W46pQ3HmC6edseUwildZ48qb73vUshlpBR\n8AWxJaYQoCmjUJVSqHOOAqy+uWJdMzol1c5GKRJotrd3UYVRMI1i1vmnBpqLuI9s10Dd7iOXUnDF\nOfKM8vtFKWTFFHSnXXSJDokplMziYu/E6zIKVc5TqDPIrI+bGlPQiwceObL685BSiGmzpo1CEfdR\n1SmpbXQf5e3QByWmUEY8wSzPRmIKORkf71nqfnYf6ZuwH9xH4I4rZCmFokah7CAzlBNTCCmFPO4j\nVwdc1X3hCp7XaRTarBS0O7etRkHcRx7Mh7hKo7BxYy9Vr8plLup2HxUxCvasZlspmO+4aKv7qIyY\ngi/QvGGDumeWltI6P1f6ZxPZR1l1yptw0U/zFLLcR2UYBbNv0czOwokT/lfPxi5zsS7dR0WNQmz2\nwNBQr/xBch/lSUmFOKUAaW3Wr+4jX6BZB899nUtMnTR1B5rFfRTnPio6mxlW9y2a/fvVGmO+TDtx\nHwWoSymY5VftPqpz5KTTCRcXyzEK9ojZbLM2KoUqU1J1+b7OJVSnOucpzM6uThqo2iiMjqr7LVVB\n1UkoJVUvd1KWUoC191ZW2eI+ClDEKHS7vRdjxx5Lj3qrcB81oRSgd4OldLqutFRX59h2pVBlSqou\n3+ebjqmTpqr7YsMGNTAwj5diFPI8B+Y6Um13H7mu2/Cw+v6++6o1CqF+SVJSAxQxCo8+qhpty5b4\nY1WpFJqYpwC9G6xK91EZgea2uo9ilUKRlNQq7wv7vFJiCnlH+fqYbVYKoes2Pg7f/W55RsEO+Jeh\nFCSmQPoJp8q/FFdIKk3NU9DHnp8vbhRCSqHN7qMqlULb3Uew9ryqdh+Zx+zHmAKo+p84UU5MQZeX\nahQkpuBh06be33UZhUEKNJvHXq+B5qpSUnX5ZRmFupRCnUahzUrBF1MAVf+dO8sz1HXFFHQ8pwit\nNwpNKYWyH1C91lGn0x9G4dxz1UKC5o3pUgpaFrfVfaRdZ4uLKi0wdSVX+21zrvLzxBT6xX2Ut1M0\njULbYwo+91FZriNdXopRyJuSqjOditC0Ufgr4CDwLd8Gg+I+AtVpnDhR/0OSJ6YwPKwWEtRzFXzv\nKe4X95EuO7Rkgat+s7Oqs/Z12GWlpA6q+6jNSiHLfdSkUcjrPtLHKkLTRuF64MWhDYoYhdQVDqt0\nH0HPKDShFHRMwXTHZWFOYNMdhs+HWZb7KKV+MaScu6t+IZWgyy/DfVSnUpCYQrNGQb+Gswr3kT5W\nEZo2Cl8GjoY2qFMppLhC8tCUUsjjPoLVcQWtEuyRdtuVQopKctUvFGTW5ZeVktpGpZC3TvpZarP7\nKHTdxsfLCzLr8vQ1OHKkt75YqH5NGYXEldLrxzxB1zouX/86fPaz7n2/+c20C1vlPAVQZTahFPK4\nj0C13d/8Ddx7Lxw96u4cU5XCww/7v6/DfZRVvzxKYXY2baTvcg3UGWj2xRRcy2GI+6gcxsfh4x9X\ng6xDh7L7JZeR/tSn4IoretcuFA8pQuuNwq237mPfPvX3j//4BLOzE3S7vRHre9+rpNhll63d9+qr\n4UlPij9WlYFmaIf7KOWGefWr4aabVMe4YQO8611rt0lpszPPhNtu83/ftFE480w1t8UkSymY1zT2\n5UB1u4/s8/IphWPHVn+2XtxHruv7hjfA4x9f3vFe9SrVqU9Pq8Htb/1Wdv3sgcMv/7IyLM94hvrf\nvD6Tk5NMTk4Cypj927/lr2vrjcJP/uQ+3va23v+jo6tX05yagte9Dq68svixBjXQnNd99PSnq58Q\nKW3mWmTPpMqU1JiyYxYBdJU/PZ3W8dUdaLbPa24OduzIrtOgB5pDKakveEG5x3vOc9RPLCMjKvaw\ntKQGZIuLcPDg6jRxs+4TExNMTEysfDc0dF3uujYdU8jEfpBTo/ipx6oy0NyU+yivUYghpc18r/nU\nNB1T2LlzbRpuTExhZiat46tzmQtY2+4uP3+Zb16D/kpJbaPRgtXX5OGHVUq1zyiUSdNG4aPAPwOP\nBR4EXmNvULdRmJ0dPPdR3phCDClt5nvNp6Zp99HwMJx3nnJHakIT13T5MzNpHV/dy1zYRqHqN69B\nfyiF1FTiujGvib5+9nWsou5Nu49ekbVByCh0OuE1yVOpw300MwMXXlh+2VnHzRNTiCGlzbZv700G\nc42+m3YfQa8Dvegi9X+d7qPQcYpgG2NJSc2XSlw35jVxGYVBVQqZhIzCgQPhNcnzHGtQ5ym0wX00\nNBR2IVVRv+Fh5ZOdmUkzCpqYQHOq+6juQLPLfSST19JTievGvCb798PFF6+NDYlRoOeugHJdR7rs\nKrOPdEyhiRnNp06tft91WaSqq7qNAqh6HTuWzyhkKQUdU0hp1zrfvAbKJXbkiDoGyDIXEH7zWlsw\nYwpTUyrrqA73UV8aBa0Uynozkl12VY3dpFI4fjx7mYc8pKqrJoyCTresSimU5T6q6r4YGVFrWelY\niSiF7OyjNmDHFGyjIErB+N80Cv2kFJo0CrGdYiqpbdZvRqGumELV79mwlyxZ7zEF7VY8caKd9YO1\nMYUnPlGlqM7MqM/EKBj/96tRaMp9ZCqFsukHozA6Gn/+qUqhzJTUKu8L87wkJVWR59rVie0+2r17\n7XUUo0D1RmFmRo0gygpemzSZklqlUjhxQrml9PLgIUJpqW1QCrt2rQ7mxSiFMlJSq37PhtmZ1JWS\nql/H2dZON8+1qxPbfWQbBYkpGP9XaRSmp6t7OMfGstf0r+q4VSqFlDbzzWpeWqomEA5p568fOv2i\n+ypiChs3qvM9fbr3WR3uo9AIswr30cmT6hxHmk5895Dn2tWJvibT06odd+wQpQDUbxSOH69u5KAv\n4KC5j1LazOc+0kuXlB0Ih7Tz37JFbX90ee3emMlrqR3L0NBad03d7qM6jMLx42r/Kq5pGfSLUdD9\nnJ3SLUbB+L/TUSO51PclxByrSqWgH/pBcx+ltNn558Phw2tfGViV6wjSz9988GJSUpeW0jt0OwW0\naqVguu3qSkk9dqy9rhnIf+3qQg8czH5OK+2lJdUHxrhsU2m9UbBfjKKNwpEjalRXZkcyPl5tYKxJ\npVClUUhps40b4eyz1eJeJlUahdTzN11IMUrB/J1Sp7qVQt3ZR8eOtXcUDvmvXV3omIKZeq+Nu742\nVaiw1hsF2x+pjULZriNdNlQbU6iy/NBxq3QfQdo5uVxIVRuFlPPX9ZudVSOxUGddplGoI6bQ7dbv\nPmor/WAUTPcR9K5jlQH81hsFG9MolDlxTZcNg+c+GhtTnUEVna5Wcm02CqOjaeev65cVZNZlQ/oD\nascUqnYfbd+uRpXT0/WlpHa77XXNQP5rVxf6mohRyKBKpTA8rC7EoLmP9PGq6HSHh9V5pZyTKy21\naqUA8eXr+mXFE8yyU6+pnQJatfsIeh1KKCVVZ11B8ZRUXW5baep5jMV0H+m+7oIL1JvbOp3q6i1G\nwVH+ILqPoLpON7XNmnAfQTVKoV/cRxAeZeq5OWYCQFH3EbTfKAwPtztl1lYKGzfCWWfBgw+KUlhB\njEL+44pRiNteB2VTlEJRo1C1+wjUeT30kF+V2HUqYhQ2bFDn03aj0Pb62UYB1N/33SdGYYU6jEJV\nskyXO0juI11uyjm5JrBVHVOAfEohyyjkvaZ1z1MAdV733686a1fWil2nopl4VT5LZVClq7gM9OrG\nhw4pt5FGGwVxHy0jSiH/cUUpxG1/7rlq8tqRI9W6j+qcpwA9o+Crq1knnQdfxLUyPt7+kXjb6/fg\ngyqF27w3sq5jUfrWKJQ9cc0sf9ACzfq49pyPsti0KV0puIxCVfVLPf8NG9Q7CL7znfrcR3UohV27\nso2CrlMZefBiFIoxNqaul93PiVGw2LRJTYo5fly9aL2K8kUppJGqFHbsWL0EMLRLKYB68O66q76U\n1LoCzSFftG0Uihqp1MFC3YyOttsojI6q62Wn3u/aJTGFVYyPKyt5/vnVrGRapfuoyWUuoD1GwfVa\nzqpjChs3pi0JsHs33H13fSmpdbmPDhzw19Vc6qKM9yD0g1Jos9EaG1PXy6UUQtexKH1pFDqdalxH\nuvxBdR+1JdAMa+cqVK0UUsvetQvuuafelNSq7ws9kEpxHxWhH4xC2+sHbqNgfl82fWkUoFqjIO6j\nNPK0WZ1KIY9R2L1b5exnKYWREaV8ihiFbrce99HIiIqViFFQiFFwI0bBUb64j9Jou1EYHc1nFCBb\nKehlsIukpC4uKtdWHUtM794ddh+VGVOQlNRi6LrZfd0ZZ1Tbtk0bhRcDdwP3AG+N2aEOozBo7iM9\nmm2T+6gflAJkKwVdfpGU1DpcR5rdu+NSUtdLTKHt9YO1fZ2OyQ2iUtgA/E+UYXgi8ArgCVk7jY72\nGqUKBtF9pN0bbVMK5gS29WgU9Ki8jiCzJssoiPuoPfiMgv5sEI3CM4DvAt8DFoC/BV6StZMe8faz\nUqjbKOhji1KIJ9Z9pMsvYhTapBTEKLQHfd+eccba76o0Ck0uBbUbeND4/yHgR2J23Ly5OqOwefPa\nN1CVxdhYz5VTN1Uahc2b1bsHUtDzAK65Rv1/221w9dXl1w3ynfvWrcogxCqFPCmpN94IjzyiljKo\nUyn46jo2Bu97H0xOqpm0RQ3V5s3t9tn3Q0qqfg2nze7d1d0zTRqFbvYmsG/fvpW/JyYmmJiY4Kab\n4JJLqqnUq15VnVHYvBluuaWasrO48UbYs6easn/pl9LbbM8e+JM/US93B7jsMnjOc8qvG8Bzn7t6\n7ZhYbr4ZLrooe7sPfhAe97i0sl/xitUq5Kqr0vbPy0teAk9+svu7N78ZvvIV9feTngRPeUqxY11z\nTXWz1Mvgyivh8subroWfyy6D6693f3fNNauXOZ+cnGRycrKU4zb5Su1nAvtQMQWAa4HTwO8a23S7\n3SjbIQiCICwzpORFrv69yZjC14BLgYuAUeDngZsbrI8gCMK6p0n30SLweuCzqEykDwB3NVgfQRCE\ndU+T7qMYxH0kCIKQSL+6jwRBEISWIUZBEARBWEGMgiAIgrCCGAVBEARhBTEKgiAIwgpiFARBEIQV\nxCgIgiAIK4hREARBEFYQoyAIgiCsIEZBEARBWEGMgiAIgrCCGAVBEARhBTEKgiAIwgpiFARBEIQV\nxCgIgiAIK4hREARBEFYQoyAIgiCsIEZBEARBWEGMgiAIgrCCGAVBEARhBTEKgiAIwgpiFARBEIQV\nmjIK/xG4A1gCntZQHQRBEASLpozCt4CXAV9q6Ph9x+TkZNNVaA3SFj2kLXpIW5RDU0bhbuA7DR27\nL5Ebvoe0RQ9pix7SFuUgMQVBEARhhZEKy/4ccL7j87cDn6jwuIIgCEJOhho+/heANwNf93z/XeDi\n+qojCIIwENwLXJJnxyqVQiwhw5TrpARBEIT+4mXAg0AHeBj4dLPVEQRBEARBEAShL3gxKm31HuCt\nDdelbvagYi13AN8G3rD8+Vmo4P13gH8Ezmikds2wAbidXoLCem2LM4CbgLuAO4EfYf22xbWoZ+Rb\nwEeAMdZPW/wVcBB17prQuV+L6kvvBl5UUx1LZQMqwHwRsBH4BvCEJitUM+cDT1n+eyvw76jz/z3g\nN5c/fyvw3+uvWmP8OvBh4Obl/9drW9wAXL389wiwg/XZFhcB96EMAcDfAa9m/bTFc4Cnstoo+M79\niag+dCOq3b5LH05F+FHgM8b/b1v+Wa98HLgCZeXPW/7s/OX/1wMXArcAz6OnFNZjW+xAdYQ267Et\nzkINls5EGcdPAC9kfbXFRaw2Cr5zv5bV3pbPAM8MFdxGi7EbFYTWPLT82XrkItSI4FbUBT+4/PlB\nejfAoPNu4C3AaeOz9dgWe4HDwPWoFO73A1tYn23xKPAHwAPAfuAYynWyHttC4zv3Xag+VJPZn7bR\nKHSbrkBL2Ap8DHgjMGN912V9tNOVwCFUPMGXurxe2mIEtXjke5d/n2Stgl4vbXEx8CbUoGkX6ln5\nRWub9dIWLrLOPdgubTQKU6hgq2YPqy3demAjyiB8COU+AmX99QzxC1Cd5aDzLOBngPuBjwLPR7XJ\neuA7fk8AAAHISURBVGyLh5Z/vrr8/00o4/Aw668tLgf+GTgCLAJ/j3I7r8e20PieCbs/vXD5My9t\nNApfAy5FjQJGgZ+nF2BcDwwBH0Bll7zH+PxmVDCN5d8fZ/B5O+qG3gu8HPg88Eusz7Z4GOVWfezy\n/1egsm8+wfpri7tRfvFx1PNyBep5WY9tofE9Ezejnp1R1HN0KfCvtdeuBH4SFUj6LipQsp74MZT/\n/Bsot8ntqBTds1AB10FPt/PxXHqDg/XaFk9GKYVvokbHO1i/bfGb9FJSb0Cp6/XSFh9FxVLmUQOF\n1xA+97ej+tK7gZ+otaaCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAhC\nPn4G+FLTlRCEKtjQdAUEoQ8ZBs5BLSsgCANFGxfEE4S286OoVToFYeAQoyAI6TwDtVDhzwK3NVwX\nQRAEoWG+BFy1/Pd4kxURhLKRmIIgpLEV+AXUG7+6rH5PriD0PWIUBCGNZ6NehXkD6uUu48A9jdZI\nEEpEYgqCkMbjgS+gXo05DhxvtjqCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiCIAiC\nIAiCIAgt5/8DVXrnzQiw7kgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f73d8209490>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "n = 6 # a number of participants\n",
    "K = 100 # trial\n",
    "\n",
    "p = np.random.rand(n) # random probability P_{X}: prob X\n",
    "s = np.sum(p)\n",
    "pp = p/s\n",
    "\n",
    "P = []\n",
    "for i in range(n):\n",
    "    P.append(np.roll(pp, i))\n",
    "P = np.reshape(P, (n,n))\n",
    "\n",
    "ans = circle_walk(P, n=n, i=0, K=K)\n",
    "plt.plot(range(K+1), ans)\n",
    "plt.xlabel(r'$k$')\n",
    "plt.ylabel('speaker')\n",
    "plt.ylim(-1, n)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_type": "subslide"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "-----\n",
    "(2014/12/10)上に述べたようなモデルはマルコフ連鎖で表すことができた。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "(引用)マルコフ連鎖のシミュレーション(山本浩・森隆一・藤曲哲郎, シミュレーションによる確率論, 日本評論社, 137-143, 1993.)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "## マルコフ連鎖のシミュレーション"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "時間と共にランダムな変化をする現象を考える。このような例として、バクテリアの増殖、天候の変化、人口の増減などが挙げられる。これらの現象をモデル化するときには、適当な単位時間を設定し、この何倍かの時刻$n$での状態$X_{n}$を調べる。単位時間としては1時間、1日、1年、1世代などが考えられ、もんっ第二より適切なものを選ぶ。確率論論としては確率変数列$X_{0}, X_{1}, \\cdots ,X_{n}, \\cdots$を取り扱うことになる。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "各確率変数$X_{n}$は問題により人数、位置などをあらわす。$X_{n}$のとり得る値全体の集合を状態空間といい、$S$であらわす。6章では$S$は自然数全体、または、その有限の部分集合として定式化できる問題を取り扱う。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "$\\{X_{n}\\}=\\{X_{n}: n=0, 1, \\cdots\\}$に対して、時刻$n+1$の状態$X_{n+1}$の状態$X_{n+1}$は時刻$n-1$までの状態には関係せず、時刻$n$での状態$X_{n}$だけで決まると仮定する。この仮定はマルコフ性とよばれ、現在が与えられたとき過去と現愛は独立であると表現できる。これにより、現在の状態がわかれば未来(次の時刻)の状態を過去の状態に依存せずに確率的に求めることができる。したがって、計算機を用い、乱数の値により一定の操作を繰り返すことにより$\\{X_{n}\\}$のシミュレーションが容易に実行できる。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "一般に確率変数列$\\{X_{n}\\}$に対して、$X_{n}$または$\\displaystyle M_{n} = \\frac{X_{0} + X_{1} + \\cdots + X_{n-1}}{n}$のどちらかの確率変数の$n\\rightarrow \\infty$としたときの状態が問題となる。$X_{n}$が収束すれば$M_{n}$は同じ極限に収束するが、逆は成り立たない。したがって、モデルによりどちらかの確率変数を取り上げる。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "6.1節ではマルコフ連鎖$X_{n}$の説明と記号の準備を行う。6.2節では、現金引き出し機の順番を待つ人の数を扱ったモデルを取り上げ、$M_{n}$の収束を問題とする。6.3節ではランダムウォークを取り扱い、$X_{n}$の$n\\rightarrow \\infty$の状態を考える。6.4節では個体が出生死滅していく集団の個体数を扱った人口モデルを取り上げ、集団が死滅するか個体数が爆発するかを問題とする。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "### 6．1 マルコフ連鎖"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "この章の所で述べた$\\{X_{n}\\}$に対するマルコフ性という仮定を、条件付き確率を用いて数学的に定式化しよう。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**定義 1** 確率変数$\\{X_{n} :n=0,1,2, \\cdots\\}$が次の性質(マルコフ性)をもつとき、**マルコフ連鎖**とよぶ:\n",
    "\n",
    ">任意の$n$個の状態$i_{0}, i_{1}, \\cdots ,i_{n-1}$と2個の状態$i,j$に対して、\n",
    ">$$P(X_{n+1}=j|X_{0}=i_{0}, \\cdots, X_{n-1} = i_{n-1}, X_{n} = i) = P(X_{n+1}=j|X_{n}=i)$$\n",
    ">このとき、上の条件付き確率を$p(i,j;n)$とあらわし、**時刻$n$での遷移確率**とよぶ。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "定義より$p(i,j;n)$は時刻$n$で状態$i$にいるという条件の下で、次の時刻$n+1$に状態$j$へと移る確率である。$p(i,j;n)$を$(i,j)$成分とする行列$P(n) = (p(i,j;n))_{i,j}$が定義でき、この行列$P(n)$を**と呼ぶ時刻$n$での遷移行列**という。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "$P(n)=(p(i,j;n))_{i,j}$が$n$に依存しないとき、$\\{X_{n}\\}$を時間的に一様なマルコフ連鎖という。この場合はパラメータ$$nの必要性がないので単に遷移行列、遷移確率とよび、$P,p(i,j)$とあらわすことにする。6章では、この時間的に一様なマルコフ連鎖を取り扱う。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "$\\pi(j) = P(X_{0} = j), j\\in S,$で与えられる状態空間$S$上の確率分布$\\pi$をマルコフ連鎖の初期分布という。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "時刻$1$での分布$P(X_{1}=j)$を考えよう。任意の$j$に対して、\n",
    "$$\\{X_{1} = j\\} = \\bigcup_{k}\\{X_{0} = k, X_{1}=j\\}$$\n",
    "右辺の各集合は$k$に対して互いに背反であるから、\n",
    "$$\\begin{align}\n",
    "P(X_{1} = j) &= \\sum_{k}P(X_{0}=k, X_{1}=j)\\\\\n",
    "&= \\sum_{k} P(X_{0} = k)P(X_{1} =j|X_{0}=k)\\\\\n",
    "&= \\sum_{k} \\pi(k)p(k,j).\n",
    "\\end{align}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "$\\pi$を行ベクトルと考えれば、最後の行は行ベクトル$\\pi P$の$j$成分であり、次の補題が成立する。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**補題 1** $P(X_{1} = j) = (\\pi P)(j) = \\sum_{k}\\pi(k)p(k,j)\\ \\ (j\\in S).$ "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**例 1** 3.3節例題6は時間的に一様なマルコフ連鎖である。この例では単位時間として親から子供への1世代をとる。親から2ステップ後の結果とは、孫が大学に入るとき理工系か文科系かをあらわす。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "状態空間$S$は$\\{0=\\text{理工系}, 1=\\text{文科系}\\}$、遷移確率は\n",
    "\n",
    "$$p(0,0) = 0.5,\\ p(0,1) = 0.5,\\ p(1,0) = 0.2,\\ p(1,1) = 0.8$$\n",
    "\n",
    "遷移行列は、\n",
    "\n",
    "$$P = \\left(\\begin{array}{cc} 0.5 & 0.5 \\\\ 0.2 & 0.8 \\end{array}\\right) .$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "各世代が大学で専攻する割合を次のように行ベクトルで表示する。\n",
    "$$(\\text{理工系の確率}\\ \\ \\text{文科系の確率})$$\n",
    "親が理工系であることは、初期分布が行ベクトル$(1\\ 0)$であると言い換えられ、補題1より、子供が専攻する確率は$(1\\ 0)P = (0.5\\ 0.5)$である。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "定義より$p(i,j)$は非負であり、任意の$i$に対して次の等式が成り立つ。\n",
    "$$\\sum_{j}p(i,j) = 1.$$\n",
    "ここで、遷移行列という概念を行列の性質として(再)定義しよう。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**定義 2** すべての成分が非負で、任意の行列に関する和が1となる行列を**確率行列**という。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "マルコフ連鎖が与えられたとき、遷移行列と初期分布が定まる。1．3節で適当な標本空間$\\Omega$での事象に対して加算加法性をもつ集合関数として確率$P(\\cdot)$を定義した。マルコフ連鎖を取り扱う場合、$X_{n}$の各時刻で$n$での値$x_{n}(n=0,1,\\cdots)$を順に並べた$(x_{0},x_{1},\\cdots , x_{n}, \\cdots)$を標本点と考える。したがって、標本空間$\\Omega$を状態空間$S$の元の無限列全体とし、$\\omega\\in \\Omega$に対して$X_{n}(\\omega)$は$\\omega$のf第$n$座標とする($n$は$0$より始める)。この場合、$\\Omega$は一定のものであるから、定義1でのマルコフ性は確率$P$に対する仮定となる。このとき、次の定理が成り立つ。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**定理 1** 確率行列$P$と分布$\\pi$が与えられたとき、任意の$n$と状態空間の$n+1$個の点$j_{0},j_{1},\\cdots ,j_{n}$に対して、\n",
    "$$P(X_{0} = j_{0}, X_{1} = j_{1}, \\cdots , X_{n} =j_{n}) = \\pi(j_{0})p(j_{0},j_{1})\\cdots p(j_{n-1}, j_{n})$$\n",
    "をみたす$\\Omega$上の確率$P$がただ1つ存在する。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "この定理はマルコフ連鎖を数学的に扱う場合に基本となる定理である。証明には測度論的な考えが必要となり、ここでは省略する。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**系** 定理1で定まる$P$に対して、$\\{X_{n}\\}$は$\\pi$を初期分布、$P$を遷移行列とするマルコフ連鎖である。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**証明** 定理1の条件式で$n=0$とすることにより、$\\pi$は初期分布となる。遷移行列を計算する。\n",
    "$$\\begin{align}P(X_{n+1} &= j | X_{0} = i_{0}, \\cdots, X_{n-1} = i_{n-1}, X_{n}=i)\\\\\n",
    "&= \\frac{P(X_{0} = i_{0}, \\cdots , X_{n-1} = i_{n-1}, X_{n}=i, X_{n+1} = j)}{P(X_{0} = i_{0}, \\cdots, X_{n-1} = i_{n-1}, X_{n}=i)}\\\\\n",
    "&= \\frac{\\pi(i_{0})p(i_{0},i_{1})\\cdots p(i_{n-1},i_{n})p(i,j)}{\\pi(i_{0})p(i_{0},i_{1})\\cdots p(i_{n-1}, i)}\\\\\n",
    "&= p(i,j)\\ \\ \\ \\blacksquare\n",
    "\\end{align}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**補題 2** $P = (p(i,j))$を遷移行列とするとき、$P^{n}$も遷移行列となる。 "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**証明** 行列$P^{n}$の$(i,j)$成分を$p^{(n)}(i,j)$とすれば、示すべきことは、すべての$n$に対し$\\sum_{j}p^{(n)}(i,j) = 1$であり、これを$n$についての数学的帰納法で証明する。$n=1$の場合は、遷移行列の定義式であり、明らかである。\n",
    "\n",
    "いま$n$まで等式が成立していると仮定する。$P^{n+1} = P\\cdot P^{n}$より、\n",
    "\n",
    "$$\\begin{align}\\sum_{j} p^{(n+1)}(i,j) &= \\sum_{j}\\sum_{k}p(i,k)p^{(n)}(k,j)\\\\\n",
    "&= \\sum_{k}p(i,k)\\sum_{j}p^{(n)}(k,j) = \\sum_{k}p(i,k) = 1.\\ \\ \\ \\blacksquare \\end{align}$$\n",
    "\n",
    "状態空間の$n+1$個の点$j_{0}, \\cdots ,j_{n}$に対して、事象$A_{j_{0}, \\cdots ,j_{n}}$を次で定義する。\n",
    "\n",
    "$$A_{j_{0}\\cdots ,j_{n}} = \\{ X_{0} = j_{0}, \\cdots , X_{n}=j_{n}\\}.$$\n",
    "\n",
    "$\\pi(i) = P(X_{0} = i)>0$と仮定する。条件付き確率の定義より、\n",
    "\n",
    "$$P(A_{j_{0},\\cdots , j_{n}} | X_{0} = i) = \\frac{P(A_{j_{0}, \\cdots , j_{n}} \\cap \\{ X_{0} = i \\} )}{ P(X_{0} = i) }$$\n",
    "\n",
    "このとき、"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "(1) $j_{0} = i$の場合、\n",
    "\n",
    "$$\\begin{align}\n",
    "A_{ j_{0}, \\cdots , j_{n} } \\cap \\{ X_{0} = i \\} &= A_{ j_{0}, \\cdots , j_{n} },\\\\\n",
    "P(A_{ j_{0}, \\cdots , j_{n} } | X_{0} = i) &= p(j_{0}, j_{1}) \\cdots p(j_{n-1}, j_{n}).\n",
    "\\end{align}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "(2) $j_{0} \\neq i$の場合、$A_{j_{o}\\cdots, j_{n}} \\cap \\{X_{0} = i\\} = \\phi,\\ P(A_{j_{0},\\cdots ,j_{n}} | X_{0} =i) = 0.$上の結果を定理としてまとめる。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**定理 2** $\\pi(i) >0$とするとき、\n",
    "$$P(A_{j_{0}\\cdots ,j_{n}} | X_{0} = i) = \\left\\{ \\begin{array}{ll}\n",
    "p(j_{0},j_{1})\\cdots p(j_{n-1},j_{n}) & i=j_{0}\\text{の場合},\\\\\n",
    "0 & i\\neq j_{0}\\text{の場合}.\n",
    "\\end{array}\\right.$$\n",
    "左辺の条件付き確率は状態$i$より出発した場合の条件付き確率であるから、$\\pi(i) = P(X_{0} =i) =0$の場合には左辺は意味をもたないが、右辺は$\\pi(i)$の値に関係なく意味をもつ。右辺を条件付き確率$P(A_{j_{0}\\cdots ,j_{n}} | X_{0} = i) = P_{\\{X_{0} =i\\}}(A_{j_{0}\\cdots ,j_{n}})$と区別するために$P_{i}$と書きあらわす。これを定義としてまとめよう。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**定義 3**\n",
    "$$P_{i}(A_{j_{0}\\cdots ,j_{n}}) = \\left\\{ \\begin{array}{ll}\n",
    "p(j_{0},j_{1})\\cdots p(j_{n-1},j_{n}) & i=j_{0}\\text{の場合},\\\\\n",
    "0 & i\\neq j_{0}\\text{の場合}.\n",
    "\\end{array}\\right.$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**系** 事象$A$に対して、$P(A) = \\sum_{i}\\pi(i)P_{i}(A).$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "上で定義した$P_{i}$に関して成り立つことのいくつかの性質を挙げる。\n",
    "\n",
    "(1) $\\pi(i)> 0$の場合、$P_{i}(A) = P(A|X_{0}=i).$\n",
    "\n",
    "事象$A$に対して確率$P_{i}(A)$を$i$より出発した場合のAの確率と言いあらわす。実際の用いる場合、系のように$\\pi(i)$をかける場合が多い。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "(2) $P_{i}(X_{1} = j) = \\sum_{k}P_{i}(X_{0}=k,X_{1} = j) = \\sum_{k}P_{i}(A_{k,j}) = P_{i}(A_{i,j}) = p(i,j).$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "(3) 状態空間の$n$個の点$j_{1},\\cdots ,j_{n}$に対して$B_{j_{1},\\cdots ,j_{n}} = \\{X_{1} = j_{1}, \\cdots , X_{n} =j_{n}\\}$とするとき、\n",
    "$$\\begin{align}P_{i}(B_{j_{1},\\cdots ,j_{n}}) &= \\sum_{j}P_{i}(\\{X_{0} = j\\} \\cap B_{j_{1},\\cdots ,j_{n}})\\\\\n",
    "&= \\sum_{j}P_{i}(A_{j, j_{1},\\cdots ,j_{n}}) = P_{i}(A_{i, j_{1},\\cdots ,j_{n}})\\\\\n",
    "&= p(i,j_{1}) \\cdots p(j_{n-1},j_{n}).\n",
    "\\end{align}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "(4) $p(i,j)>0$と仮定する。\n",
    "$$\\begin{align}P_{i}(X_{2}= k|X_{1} = j) &= \\frac{P_{i}(X_{1} = j, X_{2} = k)}{P_{i}(X_{1} = j)} = \\frac{P_{i}(B_{j,k})}{P_{i}(B_{j})}\\\\\n",
    "&= \\frac{p(i,j)p(j,k)}{p(i,j)} = p(j,k)\\\\\n",
    "&= P_{i}(X_{1=k})\\end{align}$$\n",
    "\n",
    "定理2を一般化した次の定理も同様にして証明できる。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**定理 3** 任意の$m,n$と$n+1$個の状態$j_{0},j_{1},\\cdots ,j_{n}$に対して、\n",
    "$$P(X_{m} = j_{0}, \\cdots ,X_{m+n} = j_{n}|X_{m} = i) = P_{i}(A_{j_{0},\\cdots ,j_{n}}).$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "定理3の応用として2ステップの変化を考えてみよう。時刻$m$から時刻$m+2$への変化の確率$P(X_{m+2}=j|X_{m}=i)$を計算する。\n",
    "$$\\begin{align}P(X_{m+2} =j|X_{m}=i) &= \\sum_{k}P(X_{m+2} = j, X_{m+1} = k, X_{m}=i) | X_{m}=i)\\\\\n",
    "&= \\sum_{k}P_{i}(A_{i,k,j}) = \\sum_{k}p(i,k)p(k,j).\n",
    "\\end{align}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "最後の和は$P^{2} = P\\times P$の$(i,j)$成分をあらわし、$m$に無関係である。同様にして、次の定理も証明できる。\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**定理 4** $n$ステップの変化は遷移行列の$n$乗であらわされる。任意の$m\\ge 0, n\\ge 1$と$i,j$に対して、\n",
    "$$P(X_{n+m} = j|X_{m} = i) = p^{(n)}(i,j)$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "補題1と同様にして、次の系が一般に成り立つ。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**系** $P(X_{n} = j) = (\\pi P^{n})(j) = \\sum_{k}p^{(n)}(k,j).$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**例2** この系を用いれば3.3節、例題6は簡単に解くことができる。孫が専攻する確率は$(1\\ 0)\\cdot P^{2} = (0.35\\ 0.65)$となるから、孫が理工系である確率は0.35である。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "補題2と定理4より、マルコフ連鎖の変化を$n$ステップごとにみることは$P^{n}$を1つの遷移行列とするマルコフ連鎖をみることと同じである。状態空間$S$上の関数$f_{j}(s)(j\\in S)$を\n",
    "\n",
    "$$f_{j}(s) = \\left\\{ \\begin{array}{ll}1 & s=j\\text{の場合}\\\\\n",
    "0 & s \\neq j\\text{の場合}\\end{array}\\right.$$\n",
    "\n",
    "と定義し、任意の状態$j$に対して、確率変数列$\\{V_{jn}\\}$を\n",
    "\n",
    "$$\\begin{align}V_{j0} &= 0\\\\\n",
    "V_{jn} &= f_{j}(X_{0}) + \\cdots + f_{j}(X_{n-1})\\ \\ \\ \\ (n\\ge 1)\\end{align}$$\n",
    "\n",
    "と定義すると、$V_{jn}$は時刻$n-1$までに$j$を訪れた回数になる。(条件付き)確率$P_{i}(A)$による確率変数$X$の期待値を$E_{i}[X]$とあらわし、条件付き期待値とよぶ。これは4．1節での期待値の定義で、確率$P(A)$を(条件付き)確率$P_{i}(A)$に置き換えたものである。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**定義 4**\n",
    "$$E_{i}[X] = \\sum_{k}P_{i}(X=k).$$\n",
    "\n",
    "これを用いれば、\n",
    "$$E_{i}[f_{i}(X_{k})] = P_{i}(X_{k} = j) = p^{(k)}(i,j).$$\n",
    "\n",
    "$k$に関する和をとれば、$E_{i}[V_{jn}]$は行列の和$I + P + \\cdots + P^{n-1}$の$(i,j)$成分である(ただし、$I$は単位行列とし、$P^{0} = i$とする)。この結果をまとめると次の定理5となる。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**定理 5** 状態$i$より出て、時刻$n-1$までの間に状態$j$を訪れた回数の平均は行列$I + P + \\cdots + P^{n-1}$の$(i,j)$成分で与えられる。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end"
    },
    "slide_helper": "slide_end",
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "**例 3** ここで、3.3節例題6をもう一度考える。この場合、$V_{jn}$は$n$世代までの$n$人のうち理工系の人数の平均をあらわす。\n",
    "\n",
    "\n",
    "$$\\begin{array}{ll}\n",
    "P = \\left( \\begin{array}{cc}0.5 & 0.5 \\\\ 0.2 & 0.8\\end{array}\\right) &\n",
    "P^{2} =\\left( \\begin{array}{cc}0.35 & 0.65 \\\\ 0.26 & 0.74\\end{array}\\right) \\\\\n",
    "P^{3} =\\left( \\begin{array}{cc}0.305 & 0.695 \\\\ 0.278 & 0.722\\end{array}\\right) &\n",
    "P^{4} =\\left( \\begin{array}{cc}0.2915 & 0.7085 \\\\ 0.2924 & 0.7166\\end{array}\\right) \\\\\n",
    "I + P + \\cdots + P^{4} =\\left( \\begin{array}{cc}2.4465 & 2.5535\\\\ 1.0214 & 3.9786\\end{array}\\right) & \\ \\\\\n",
    "\\end{array}$$\n",
    "\n",
    "これより、親が理工系ならば5世代のうち理工系は約半分、親が文系ならば5世代のうち約1人が理工系であることが期待される。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "subslide_end",
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "### 適用"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "slide_end",
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "以上の一般のマルコフ連鎖についての議論を踏まえて、現在考えている系に適用すると、確率変数は図でラベル付けされた数字であり、遷移行列はその時刻に依存しない。また、各状態$i$において、遷移確率はその点からの位置関係のみで決まっている。\n",
    "\n",
    "点$i$から時計回りに$a$目盛進む確率を$p_{i}(a)$と書くことにしており、しかもこれが$i$に依らないとき($p_{i}(a) = p(a),\\  ^{\\forall} i$)を考えていたので、このマルコフ連鎖の遷移行列$P$は\n",
    "\n",
    "$$P(n) = P = \\left( \\begin{array}{cccccc}\n",
    "p(0) & p(1) & p(2) & p(3) & p(4) & p(5) \\\\\n",
    "p(5) & p(0) & p(1) & p(2) & p(3) & p(4) \\\\\n",
    "p(4) & p(5) & p(0) & p(1) & p(2) & p(3) \\\\\n",
    "p(3) & p(4) & p(5) & p(0) & p(1) & p(2) \\\\\n",
    "p(2) & p(3) & p(4) & p(5) & p(0) & p(1) \\\\\n",
    "p(1) & p(2) & p(3) & p(4) & p(5) & p(0) \n",
    "\\end{array} \\right)$$\n",
    "\n",
    "のようにかける。\n",
    "\n",
    "定理5から、状態$i$より出て、時刻$n-1$までの間に状態$j$を訪れた回数の平均は行列$I+P + P^{n-1}$の$(i,j)$成分で与えられるから、\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false,
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "subslide_end",
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "def cal_f(n, P, K):  \n",
    "    f = np.identity(n) + P\n",
    "    if K == 1:\n",
    "        return f\n",
    "    for k in range(2, K+1):\n",
    "        _k = 2\n",
    "        tmp = P\n",
    "        while _k < k+1:\n",
    "            tmp = np.dot(tmp, P)\n",
    "            _k += 1\n",
    "        f += tmp\n",
    "        \n",
    "    return f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false,
    "internals": {
     "slide_type": "subslide"
    },
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div style=\"max-height:1000px;max-width:1500px;overflow:auto;\">\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>0</th>\n",
       "      <th>1</th>\n",
       "      <th>2</th>\n",
       "      <th>3</th>\n",
       "      <th>4</th>\n",
       "      <th>5</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td> 17.767557</td>\n",
       "      <td> 16.464353</td>\n",
       "      <td> 16.534830</td>\n",
       "      <td> 16.594450</td>\n",
       "      <td> 16.918806</td>\n",
       "      <td> 16.720004</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td> 16.720004</td>\n",
       "      <td> 17.767557</td>\n",
       "      <td> 16.464353</td>\n",
       "      <td> 16.534830</td>\n",
       "      <td> 16.594450</td>\n",
       "      <td> 16.918806</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td> 16.918806</td>\n",
       "      <td> 16.720004</td>\n",
       "      <td> 17.767557</td>\n",
       "      <td> 16.464353</td>\n",
       "      <td> 16.534830</td>\n",
       "      <td> 16.594450</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td> 16.594450</td>\n",
       "      <td> 16.918806</td>\n",
       "      <td> 16.720004</td>\n",
       "      <td> 17.767557</td>\n",
       "      <td> 16.464353</td>\n",
       "      <td> 16.534830</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td> 16.534830</td>\n",
       "      <td> 16.594450</td>\n",
       "      <td> 16.918806</td>\n",
       "      <td> 16.720004</td>\n",
       "      <td> 17.767557</td>\n",
       "      <td> 16.464353</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td> 16.464353</td>\n",
       "      <td> 16.534830</td>\n",
       "      <td> 16.594450</td>\n",
       "      <td> 16.918806</td>\n",
       "      <td> 16.720004</td>\n",
       "      <td> 17.767557</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>6 rows × 6 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "           0          1          2          3          4          5\n",
       "0  17.767557  16.464353  16.534830  16.594450  16.918806  16.720004\n",
       "1  16.720004  17.767557  16.464353  16.534830  16.594450  16.918806\n",
       "2  16.918806  16.720004  17.767557  16.464353  16.534830  16.594450\n",
       "3  16.594450  16.918806  16.720004  17.767557  16.464353  16.534830\n",
       "4  16.534830  16.594450  16.918806  16.720004  17.767557  16.464353\n",
       "5  16.464353  16.534830  16.594450  16.918806  16.720004  17.767557\n",
       "\n",
       "[6 rows x 6 columns]"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f = cal_f(n,P,K)\n",
    "from pandas import DataFrame\n",
    "DataFrame(f)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {},
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "検算"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false,
    "internals": {
     "slide_helper": "subslide_end"
    },
    "slide_helper": "subslide_end",
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div style=\"max-height:1000px;max-width:1500px;overflow:auto;\">\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>0</th>\n",
       "      <th>1</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td> 2.4465</td>\n",
       "      <td> 2.5535</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td> 1.0214</td>\n",
       "      <td> 3.9786</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>2 rows × 2 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "        0       1\n",
       "0  2.4465  2.5535\n",
       "1  1.0214  3.9786\n",
       "\n",
       "[2 rows x 2 columns]"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "DataFrame(cal_f(2, np.array([[0.5, 0.5], [0.2, 0.8]]), 4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false,
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "subslide_end",
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Counter({4: 23, 2: 21, 0: 20, 3: 15, 1: 13, 5: 9})\n",
      "[ 17.76755711  16.46435264  16.53482981  16.59445001  16.9188062\n",
      "  16.72000422]\n"
     ]
    }
   ],
   "source": [
    "import collections\n",
    "count = collections.Counter(ans)\n",
    "print count\n",
    "print f[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "subslide_end",
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "平均してみて、あっているかどうか。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false,
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "subslide_end",
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 167.74596838  166.67534215  166.64948632  166.72541605  166.70481025\n",
      "  166.49897685]\n",
      "1001.0\n",
      "[195.96666666666667, 137.98, 194.79666666666665, 137.38666666666666, 194.93333333333334, 139.93666666666667]\n",
      "1001.0\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAAERCAYAAACaUQc3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmYFNW5x/Fvwwgii4AY1gkDCMimgmxRhAbUoCRKNEZR\n75NLvGiiRuOOiZHRJMYQjRqvy40RVFQMUSQQDQ+KdAIiizIiiwgMoCyyKLsgDDN9/3ir7Z6ebqam\nZ6Z6qd/nefrpqtPV1aeo4a1Tb506BSIiIiIiIiIiIiIiIiIiIiIiIiIiUkPygbnASmAFcJNT3hx4\nC1gDzAaaOuUFwCGgyHk96WFdRUQkRa2AM5zpRsAnQDdgAnCnU34X8KAzXQAs97B+IiJSC6YD5wKr\ngZZOWStnHhTsRUSyXgHwKdAY2B1THoiZLwAOYCmcEDDIs9qJiEi1NQI+AEY587vjPt/lvNcDmjnT\nfYDPsIODiIikUZ6LZY4DXgMmY2kcgO1Y+mYb0BrY4ZQfcV4AS4FioLMz/Y1OnTqFi4uLq1VxEREf\nKgZOSeWLdSr5PAA8C6wCHo0pnwH82Jn+MdGDQAugrjPdEQv06yvUtriYcDic9HXkSPLPsuE1fvz4\ntNdB26bt0/bl3gvolEqgh8qD/dnA1cBQot0pR2C9b87Dul4OI9obZzCwzFnu78B1wJ6qVKi4GNq3\nh7lzbX7WLHj44aqsQURE4lWWxplP8gPCuQnKpjmvlI0eDY0awcKFMHQozJ8Pc+bAbbdVZ60iIv5W\nWcvec5s2wdix8OGHNr9hAyxbBqWl6a1XVQSDwXRXodbk8raBti/b5fr2VUcgTb8bdvJPFTRtCm+8\nAWPGwJo1cNZZ8N57sHIldO/ucS1FRDJIIBCAFON2xrXsDx6E3r1hyxbYv99a9t/5DhQVpbtmIiLZ\nK6OC/dGjlq5p0AB69IBFi2D3bhg5EpYurfz7IiKSWEYF+0OH4IQTIBCAgQNh0iT49rfhzDPVshcR\nqY6MCvYHD1qwB7jySnjlFejQwdI6RUWQJM0vIiKVyLhg36CBTQ8YAJ06WbBv2RIaNoSNG9NaPRGR\nrOVmuATPRNI4YKmcceOgSRObj7TuO3RIX/1ERLJVRgX72DQOwE9+Ep3u0wfef99y+H37el83EZFs\nlrFpnHi9e9uwCf37Ww8dERFxL6OCfWwaJ97gwXDzzfYeCnlaLRGRrJdRwT4+jROreXOYMAEuvBDe\necfbeomIZLuMC/bJ0jgRw4cr2IuIVFVGBftjpXEizjjDBkvbtevYy4mISFRGBftjpXEi6ta17pef\nfupNnUREckHWBXuwh5so2IuIuJdRwf7Qocpz9mB97T/7rPbrIyKSKzIq2KtlLyJSO7Iy2KtlLyJS\nNZUF+3xgLrASWAHc5JQ3B97CHjg+G2ga8527gbXAauD8qlTGbRpHLXsRkaqpLNiXALcAPYCBwA1A\nN2AcFuy7AHOceYDuwOXO+wjgSRe/8Q217EVEakdlgXgb4Dz6mwPAx0Bb4CLgeaf8eWCUM30xMAU7\nSGwE1gH93VbGbbBv1crGx/n6a7drFhHxt6rk7AuA3sAioCWw3Snf7swDtAE2x3xnM3ZwcMVtGqdO\nHcjPh7Vr3a5ZRMTf3Ab7RsBrwM3A/rjPws4rGdfPl3Lbsgd7ktWECW7XLCLib27Gsz8OC/STgelO\n2XagFZbmaQ3scMq3YBd1I9o5ZRUUFhZ+Mx0MBgkGg1UK9nfcAV27wrJlcPrp7r4jIpJNQqEQoRoa\n5jfg4vPngS+xC7URE5yyP2AXZ5s6792Bl7E8fVvgbeAUKrbuw+EED5Q99VSYPt3e3Rg3DurVg/vv\nd7e8iEg2CwQCUHncTqiyNM7ZwNXAUKDIeY0AHgTOw7peDnPmAVYBU533fwHXU0tpHIARI2DWLPfL\ni4j4VUpHiBqQsGXfogWsXm3vbhw5At/6ll2oPfnkGq6hiEiGqU7LPqOCfcOGsGOHvbt1ySXWO+fb\n34bu3eF//qcGaykikkFyItiXlkL9+tZar1OFDqHbtsFLL8Hnn8Orr8KGDRBI11aJiNSinAj2u3fb\nOPV79qS6QujYEWbMgF69aqCGIiIZpjYv0Hpm71448cTUvx8IwPe/DzNn1lydRERyRc4Ee4Af/ACe\nesr63ouISFROBfuhQ+GBB+DCCy2tA7BxI5SVVbt6IiJZLaeCPcDVV9vF3k2boLgYevZUX3wRkZwL\n9oEA9OsHS5bANddYn/3336/+ekVEslnOBXuwYP/Xv8L69fCHP8AHH9TMekVEslVGBfumTStfzo1+\n/Sx18+Mf27SCvYj4XUYF+5ps2QcCFuw7dLAxd7Zvh3377A5dERG/yclg36IFfPwxnHKKBf0zz7Qu\nmcOGwdln2w1cIiJ+kjHBfs+emgv2YGPdRzzyCPz73xboR46Ea6+tuHxJibv1/u53MHFitGuniEg2\ncPPwEk/UZMs+Xs+eMHeuTR86ZC3+pUuhTx8rC4dh0CD4zW/g/POTr+e99+CJJ2yEzV274Pbba6e+\nIiI1LWNa9rUZ7GM1aAB33QW/+lW0dV5UBIsXw5w5Nv/883YGcPRo+e/eeqs9CvGFF+xs4ciR2q+v\niEhNyJiB0Hr0gL/9zVrhte3wYRgwAL73PbsBa/lyu3gbDsO8eZbu2brVLu6ecYYF+I8+suEY1q+H\nunVh+HAYM8Zu4oooK4PzzoPLL4exYzX6pojUrJwY9bJdO0uT5Ocn+UYNW70a/uu/bIiF5cvhT3+C\nvn0tvTNoEKxaBW++acMnt21rB4ju3eGee+z7//wnjB9vN2xFgvqcOXDDDdCokaWLXn5Zz8cVyWX7\n9lks2LLFGos7d9rr4EFrLH71FXTpYg3Ixo0te9GyZeq/lxPBvnFj+wdr0iRNNcJa8ccdZ103n3zS\nyvbtg5tusjTP229DmzZWXlYG3brBM8/A4MFWdtVVMHAg3Hhj9Gauv/89PdsCdsBZuNAOaCKSuiNH\nLN27YkX0tXo1fPklBIPWEAwE7HpeixZw/PE2LtcJJ1jDccEC+PprW76kxB7QdOaZ0Ly5rbtlSzj3\nXFvHtm32OnDAnsT32WfWqBw4EB55JMuDfWmpPTi8pKRqDy6paf/7v/YQlLvvttZ5Zf7yF3j8cZg2\nzXZyQYGNx3PSSbB/P7Rvb+mfdu1qveoVlJXBFVdY3SZOtAfDnH++HcwaNrThJH7zG+u19NBDFb+/\nYYM99atuXSgstG1ZudK2Z/16+0P97nfTs20iNa2sDF5/HT75xFrhp54KmzfDP/4B8+dboG/f3s7U\ne/a0tHO3bvaEvLp1q/Y7Bw7YNcqiInuvV8/G8nrzTft/mp9vB4xGjewsoU0b6N3bMh/33lu7wX4i\nMBLYAUQeC3I68DTQENgIXAXsBwqAj4HVznLvYQ8dj1cu2Ff3wSXpEg5b75wJE+CWW6wV/be/RT+/\n/Xa7k/cXv4Arr6zaw9QTKS621nrPnpZWql8fvvjC7iG4+mqYOhWeftoOnuGwBfJ777Uzji5dLMAf\nPGgHgVAI7rwT7rvPUlft20d/Z948+OEP7UJ23bp2sNi61R4Kc9pptuyiRfDWW/aHeNZZNg5R//7V\n2z6R2rRjBzz4oKVdIjGnd29rVS9ZYv+XzjkH/vMfa5W3bm2t7QsusCBfnfRLTantNM45wAHgBaLB\nfglwKzAPGAN0AO7Fgv3MmOWSKRfsZ8+2oLRwYVWqnjkuucT+gGbOtNZuRFmZpX4eecQOZAsWpHbR\n9sABG+vnd7+z77dta63sHj2sBT5smP0bdu0Kzz1nB4HSUmt5xPriC8jLs/oOHGjDQd9xh52FPPWU\ntTKWL7fPX37ZLjYfS2mpnV7Omwd//KPdy3DqqVXfPpGIo0ft77RVq8qXDYet0fH555ZujaQ85s2z\nRsqoUZYteO45+7+4cCH8939bw6RlS0vDrFxpB4ETTrBG2XHH1fYWVo8XOfsCygfxPUBkJJt8YBbQ\nI8FyyZQL9mPHWpC47TaXtckwn3xiefpZsxKf0pWVWZ/+wkL7A3Tj6FGYPNlaHFOnWl6wsNAu8KxY\nYaeaq1bZPQMtWtg1hVNPdXfNIxyOHnR27LCDRV4erFtnf+wvvmhdT6ti0iT47W/tP9TJJ1ftuyIL\nFsCzz9pjRQ8fts4Sd91lKcM33rC/qyZNrPNEly7W0Fi40P5mu3e3a37bt9uZ5uDBdgb75pvWgLnx\nRgvmgwdbjjybpSPYvwtMAP6BtfALgSbOciuAtcBe4B5gfoL1fRPsS0rsdOmDD8qnEnLNzJnWe2fp\n0vLlixbZAHCRO363bLGgvmiRBfGRI63LZ8eOtVe3AwfszGDEiOqlmn75S7sgPW6ctZ5KSqy1dOON\nFc8yJPfs3m0Nn+JiO+Ps1MkCc7yDB+2MF6xlvWyZNW5+/nO49FKLBy+9BI8+aumTkSMtvbJ3r+XP\nN2ywRkmrVvZ3Vq+et9uZTukI9l2BPwMnATOAm4AWQD0sj78b6ANMx1r8++PWFx4/fjxgwW3OnCDr\n1wdTqX/WKCuzHOGMGdHumHv2WBAsLbW8eMeO1vIfNcpSKJGr89kiHLZUztix1gV16VL7z3nyydZy\nq8qFLMkMq1ZZqm/3bus10qmTnQ2uX2/pwi5dLPAuXGiNhq5doXNna7itWmX/v2OFw9YJY/hwC9Jd\nutjF/i5d7OxSyguFQoRCoW/m77vvPvA42MfqAkwGBiT4bC5wGxDXno227GfPtt4gs2e7rEkW+9Wv\n7FSzVy+7aHrnnfYHfuWVcN11lqt84AGbzmYrV9o2DRhg+/aHP7Q7l3/6U7u4nJ9vB7j162HNGrvo\n/L3vWfDIpoNbrvjsM+sSuGKF/Q0uWmRnZV27WmrlnnssIOflWcv9W9+y4Hz4sO2/gwftIn1+vvZf\nbUtHy/5kYCc23MJzwDvOewusVV8KdAT+A/TEcvyxvgn2r74Kr7xi77nuk08s2A0ebP+hTj/dLuw2\na2ansSUlFhRzzeHDls9fvNha+2edZT15Wra06wwHDtgFtl277Kzg3nvTXePcFw5byu2JJ6wF3rq1\n7YvWrS0HXqeOHYxHj7YeWJIZqhPs3Zw4TQGGYIF8EzAeaATc4Hz+GhboAQYD9wMlQBlwHRUDfTn7\n9qX3Riovde1qLacTT7T0Rt++dvMFWKspV09j69e3Pv1g+dklS6yHRLNmVhYOW5fVpk0tB9uund15\nOGeOzZ93XvTf6dNP7SARmZeqOXLE+nQ//DC8+6612keNyvxeKFJ9ab+p6tFH7YLLY4+lqSaSUd5/\n3wacC4ctFfTSS9bjYuFCC/6/+IVd23j9dXvfts1uHFuyBF57ze5tOHzYXmPGWPfUVMT2WKoN4bAd\nuFauhHfesTOdhg3t+kbknoyTTrIGQf36la8rkgtPZu9eu1nwxRdtvb162bRfGlq5orZb9rXKTy17\nqVzfvnZTS8TPfmb3KfTta7nixYvt+s7w4XYj2dNPw4UXWlAPhWyMow4dLB02dCjcf79dVGzQwLrt\nVWbBAru2sGGDpZ6uuebYd1OHQnYQuvLKaI+jyM12kyfbRemvvrLus5EeVwcP2nWMoiJLkfTvb7nx\nkhK7meeZZ+zi5c6dNn/xxdYza9kyS3m1amUHs379LIhPm2bLde1qPaKuuspSNNOm2dlTmza2zpEj\nraeMusb6U9pb9rffbqfld9yRpppIxguHrYU/alQ08P7f/1mwvPXWxN37wK6LPPywXRTfuDHayj98\n2ILz++/bPQbbtkW7/73xhq27Wze772P+fAueV19taaZdu6Kv2bOty/AFF1j9GjWyFvmuXXbj2x//\naOmRVassXdKtm10I3bTJutM++2zlqbuNG61OX3xhB7zGja2+eXmWhmnWzIL4GWfYzUQ33GC/X1Bg\nw10cOmS594suio7hJNkrqwdCu/Zau9Kf7T1QJLPt3WtBs0kTG5DqmWesRT1tml2UvOwya4UPGWJd\nByO++MLOLGbPthZ68+aWBmne3NJIt91mAf7AAWuJf/WVXZNp1658GmjNGhtyonlzOxCcdFLtbOdX\nX0VHblXPmNyT1cH+iivsNHX06DTVREQkS1Qn2Kf9SVXK2YuI1D4FexERH1CwFxHxAQV7EREfULAX\nEfGBtPbGCYetH/LBg/4aplREJBVZ2xvn0CEL9gr0IiK1K63BXikcERFvKNiLiPiAgr2IiA8o2IuI\n+EDag33jxumsgYiIP6Q12B8+rCcOiYh4Ia3BvrTUhpUVEZHapWAvIuIDboL9RGA7sDym7HTgPeAj\nYAYQm3m/G1gLrAbOP9aKFexFRLzhJthPAkbElf0VuBM4DXgdiDxUsDtwufM+AnjyWL+hYC8i4g03\nwX4esDuurLNTDvA2cKkzfTEwBSgBNgLrgP7JVqxgLyLijVRz9iuxwA5wGZDvTLcBNscstxlom2wl\nCvYiIt6o5Nn2Sf0E+DPwayxnf+QYy4YTFRYWFrJ4sT2kORQKEgwGU6yKiEhuCoVChEKhGlmX26Ey\nC4CZQK8En3UBJgMDgHFO2YPO+yxgPLAo7jvhcDjMY49BcTH8+c9VqrOIiC+lY4jjk2O+fw/wlDM/\nA7gCqAd0wHL7i5OtRGkcERFvuEnjTAGGAC2ATVhLvRFwg/P5a8BzzvQqYKrzfhS4niRpHFCwFxHx\niptgPzpJebLkywPOq1IK9iIi3tAdtCIiPqBgLyLiAwr2IiI+oGAvIuIDaQ32R48q2IuIeEEtexER\nH0h7sM9LdcAGERFxLe3BXi17EZHap2AvIuIDCvYiIj6gYC8i4gMK9iIiPqBgLyLiAwr2IiI+oGAv\nIuIDCvYiIj6gYC8i4gMK9iIiPqBgLyLiA26C/URgO7A8pqw/sBgoApYA/ZzyAuCQU14EPHmsFSvY\ni4h4w02wnwSMiCubAPwa6A3c68xHrHPKewPXH2vFCvYiIt5wE+znAbvjyj4HTnSmmwJbUvlxBXsR\nEW+kOpr8OGA+8BB2wPhOzGcdsBTOXuAeZ7mEFOxFRLyRarB/FrgJeB24DMvrnwdsBfKxM4E+wHSg\nB7A/fgWFhYUUF8PkyVC3bpBgMJhiVUREclMoFCIUCtXIugIulysAZgK9nPl9QJOYdewhmtaJNRe4\nDVgaVx4Oh8MMGgS//z2cc06V6iwi4kuBQADcx+1yUu16uQ4Y4kwPA9Y40y2ASGKmI9AZWJ9sJUrj\niIh4w00aZwoW2FsAm7DeN9cCTwD1sa6W1zrLDgbuB0qAMuA6rNWfkIK9iIg33AT70UnKByQom+a8\nXFGwFxHxhu6gFRHxAQV7EREfULAXEfEBBXsRER9QsBcR8YG0BvujRxXsRUS8oJa9iIgPpD3Y56U6\nOo+IiLiW9mCvlr2ISO1TsBcR8QEFexERH1CwFxHxAQV7EREfULAXEfEBBXsRER9QsBcR8YG0Bftw\n2F510nq4ERHxh7SF2tJSC/SBlB6dKyIiVZHWYK8UjoiIN9wE+4nAdmB5TFl/YDFQBCwB+sV8djew\nFlgNnJ9spQr2IiLecRPsJwEj4somAL8GegP3OvMA3YHLnfcRwJPJfkPBXkTEO26C/Txgd1zZ58CJ\nznRTYIszfTEwBSgBNgLrsLOAChTsRUS8k+oAw+OA+cBD2AHjO055G2BhzHKbgbaJVqBgLyLinVSD\n/bPATcDrwGVYXv+8JMuGExU++GAhX38NhYUQDAYJBoMpVkVEJDeFQiFCoVCNrMttx8cCYCbQy5nf\nBzSJWcceLK0zzil70HmfBYwHFsWtL7x1a5jevWHbthRqLSLiQwHrq55Sh/VUu16uA4Y408OANc70\nDOAKoB7QAeiM9dqpQGkcERHvuEnjTMECewtgE9b75lrgCaA+cMiZB1gFTHXejwLXkySNo2AvIuKd\ndN2/Gi4uDjN8OGzYkKYaiIhkmXSkcapNLXsREe8o2IuI+ICCvYiIDyjYi4j4gIK9iIgPKNiLiPiA\ngr2IiA+kLdgfPapgLyLilbS27PNSHYZNRESqRGkcEREfULAXEfEBBXsRER9QsBcR8QEFexERH1Cw\nFxHxAQV7EREfULAXEfEBBXsRER9QsBcR8QE3wX4isB1YHlP2ClDkvDY47wAF2APII589mWylCvYi\nIt5xMzrNJOBx4IWYsitiph8C9sTMrwN6V7ZSBXsREe+4CfbzsBZ7IgHgR8DQqv6wgr2IiHeqm7M/\nB0vxFMeUdcBSOCFgULIvKtiLiHinuoMMjwZejpnfCuQDu4E+wHSgB7A//oszZhSycycUFkIwGCQY\nDFazKiIiuSUUChEKhWpkXQGXyxUAM4FeMWV5wGYsqG9N8r25wG3A0rjy8GOPhVm7Fh5/3H1lRUT8\nLBAIgPu4XU510jjnAh9TPtC3ACLJmY5AZ2B9oi8rjSMi4h03wX4KsADoAmwCxjjllzufxRoMLMNy\n9n8HrqN8T51vHDkC9eunUGMREakyNzn70UnKxyQom+a8KvX113D88W6WFBGR6krbHbQK9iIi3lGw\nFxHxAQV7EREfSGuw1wVaERFvqGUvIuIDaQv2hw8r2IuIeEUtexERH1CwFxHxAQV7EREfULAXEfEB\nBXsRER9QP3sRER9Qy15ExAcU7EVEfEA3VYmI+IBy9iIiPpC2YF+nDuRV93HnIiLiStqCvVI4IiLe\ncRPsJwLbgeUxZa9gz5ktAjY47xF3A2uB1cD5yVaqYC8i4h03iZRJwOPACzFlV8RMP0T0oeLdsQeR\ndwfaAm9jDyovi1+pgr2IiHfctOznAbuTfBYAfgRMceYvdqZLgI3AOqB/oi/q4qyIiHeqm7M/B0vx\nFDvzbYDNMZ9vxlr4FahlLyLineoG+9HAy5UsE05UqGAvIuKd6nR+zAN+APSJKdsC5MfMt3PKKtix\no5DCQpsOBoMEg8FqVEVEJPeEQiFCoVCNrCvgcrkCYCbQK6ZsBHAXMDSmrDvW0u9P9ALtKVRs3YeH\nDw/z9tsp1FhExKcCgQC4j9vluEnjTAEWYL1qNgFjnPLLiV6YjVgFTHXe/wVcj9I4IiJpl9IRogaE\nL700zKuvpunXRUSyUG237GuFWvYiIt5RsBcR8YG0BXvdVCUi4h217EVEfEDBXkTEBxTsRUR8QMFe\nRMQH0hbs+yccC1NERGpD2h4MODdcSPC++yqUjx8ynsJgYYXywlAh9/1by2t5La/ltXwq0nYHbTic\ncBQFERFJIivvoBUREe8o2IuI+ICCvYiIDyjYi4j4gIK9iIgPKNiLiPiAgr2IiA8o2IuI+ICCvYiI\nD7gJ9hOB7cDyuPKfAx8DK4A/OGUFwCGgyHk9WSO1FBGRanET7CcBI+LKhgIXAacBPYGHYj5bB/R2\nXtfXQB2zTigUSncVak0ubxto+7Jdrm9fdbgJ9vOA3XFlPwN+D5Q48ztrslLZLpf/4HJ520Dbl+1y\nffuqI9WcfWdgMLAQCAF9Yz7rgKVwQsCgatRNRERqSKpDHOcBzYCBQD9gKtAR2ArkY2cCfYDpQA9g\nf7VrKiIita6A8hdo/wUMiZlfB5yU4HtzsaAfbx0Q1ksvvfTSq0qvdaQo1Zb9dGAY8G+gC1AP+BJo\ngbXqS7GWfmdgfYLvn5Li74qISC2ZgqVnDgObgDHAccBkrLX/ARB0lr0E64pZ5JSP9LiuIiIiIiLi\nhRHAamAtcFea61JTNgIfYWczi52y5sBbwBpgNtA0LTVLTaKb6I61PXdj+3M1cL5HdayORNtXCGwm\nejPgBTGfZdP25WPXyVZiZ9g3OeW5sv+SbV8hubH/jgcWAR8Cq7Du7ZCF+68udnGhAEsDfQh0S2eF\nasgGbGfEmgDc6UzfBTzoaY2q5xzshrjYYJhse7pj+/E4bL+uI/OH4Ei0feOBWxMsm23b1wo4w5lu\nBHyC/R/Llf2XbPtyZf8BnOC852Fd2wdRQ/vPyw3vj1VmI3Yz1ivAxR7+fm2KfwDwRcDzzvTzwChv\nq1MtiW6iS7Y9F2PXdEqw/boO28+ZLNH2QeKHOGfb9m3D/vMDHMCGM2lL7uy/ZNsHubH/AA467/Ww\nBvJuamj/eRns22IXeCM2E91R2SwMvA28D4x1ylpiqQKc95ZpqFdNSrY9bbD9GJHN+/TnwDLgWaKn\nydm8fQXYGcwicnP/FWDbt9CZz5X9Vwc7oG0nmrKqkf3nZbAPe/hbXjob+6O7ALgBSxPEivSPzRWV\nbU82butT2J3fZwCfAw8fY9ls2L5GwGvAzVS8oTEX9l8j4FVs+w6QW/uvDNuOdtgoBUPjPk95/3kZ\n7LdgF1gi8il/VMpWnzvvO4HXsdOo7Vh+EaA1sCMN9apJybYnfp+2c8qyzQ6i/4n+SvRUOBu37zgs\n0E/G7oeB3Np/ke17kej25dL+i9gLvAGcSRbuvzygGDv9qkduXKA9AWjsTDcE3sWuiE8g2ttoHNl1\ngRYq3jGdbHsiF4jqYS2rYhLnTjNNAeW3r3XM9C3Ay850tm1fAHgBeCSuPFf2X7Lty5X914JoCqoB\n8B9gOFm6/y7ArqCvw7oMZbsO2D/2h1hXsMg2Ncfy+NnY9TJyE90RojfRHWt7fontz9XAdz2taWri\nt+8nWAD5CMv5Tqf8NZZs2r5BWBrgQ6LdEEeQO/sv0fZdQO7sv17AUmz7PgLucMpzZf+JiIiIiIiI\niIiIiIiIiIiIiIiIiIiIiIiIiIhIVZ0I/KySZd6t5PMDKfzuLc73Irfsn409kvPqFNYlIiKVKKD8\neDcRAdyPFxI/QqQbQ7CxSm6JKftRCusREREXXsEe/FAETMXGYHoeG7MoMjJgpOU+HXv+wAqizyCA\naLBviI04+CF2ADlW8L4Ma9VHHknZGBuPRkREakF7oi379kApFZ/aEwnmzZz3Bs53msV9finwl5jv\nNTnG70YOBG8Bp2It/ZOrUnERL2T68xhF3ArETX9KtLUd72as1f4e1urvHPf5R8B5WHpmELDPxe+/\nBFyFjbi403WtRTyiYC+56qsk5UOwMcIHYk8EKgKOj1tmLdGHkv8W+HWSdbXChksGe6DGJWTQeOIi\nsfLSXQGRGrKf6INkjuVE7CHOX2Npl4EJlmntLPMS9sSga5Ksqx8wJ+b3V6AUjmQoBXvJFV9iXSuX\nAx+T+Fl4fPjUAAAATUlEQVScYWAW8FNgFXYR9724zwFOw54OVIY95CRRl85hQCFQH3seKtij8r6s\nxjaIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIhITfh/07B7xt4HqesAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f73b9461610>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAAERCAYAAACaUQc3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHQVJREFUeJzt3XmYHGWdwPHvkJCEKxBMgFw6CSRCMDAECCAIrXJ6AOoi\ncVV84rEILl6sByxuZkFWRGV1WfDYJSgCQbxYosIKWRoRhKBJyEUkGZKQEAiHOSYhJJNM7x+/Gqan\n6Zme6enMZLq+n+epp6vfqq5+366qX739VtVbIEmSJEmSJEmSJEmSJEmSpAqZDqwFFuSl3QHMTYbl\nyWuLy4ClwBLg9B7KoySpm94GHEXbYJ/v28AVyfgEYB6wO1ALLAN228n5kyRVSC3Fg30N8AxwcPL+\nMuAredPvBY7fqTmTJHVKd2rebyOaeBqS9yOA1XnTVwMju7F8SVKFdCfYfwi4vcQ8uW4sX5JUIf27\n8bn3AZPy0p4FRue9H5WktXHwwQfnGhoaCpMlSR1rAA4p98Pl1uxPBZ4E1uSl3Q1MAQYAY4BxwOzC\nDzY0NJDL5ap2mDZtWq/nwfJZvjSWr5rLlsvloPX8aFlKBfsZwCPAeGAVMDVJPz+Zlm8xcGfyeg9w\nMTbjSNIuoVQzzofaSZ/aTvq/JYMkaRfidfAVlslkejsLO5Xl69uquXzVXLZKqOmF78wl7U+SpE6q\nqamBbsRsa/aSlAIGe0lKAYO9JKWAwV6SUsBgL0kpYLCXpBQw2EtSChjsJSkFDPaSlAIGe0lKAYO9\nJKWAwV6SUqDcJ1V1W322nn998F9flz7tlGnUZ+qd3/md3/mdv535y2Gvl5LUB9jrpSSpJIO9JKWA\nwV6SUqBUsJ8OrAUWFKRfAjwJLAS+maTVAluAuclwY8VyKUnqllJX49wMXA/ckpf2duBs4AigCRiW\nN20ZcFQlMyhJ6r5SNfuHgHUFaRcB3yACPcCLlc6UJKmyymmzHwecDDwKZIFj8qaNIZpwssBJ3cyb\nJKlCyrmpqj8wBDgeOBa4ExgLrAFGE/8EJgF3AYcDjYULqK+vf208k8mQyWTKyIYkVa9sNks2m63Y\n8jpzgX4tMBOYmLy/B7gGeDB5vww4Dni54HMPAJcCcwrSvalKkrqoN26qugt4RzI+HhhABPqhQL8k\nfSzR3PN0uRmTJFVOqWacGcApwBuAVcC/EJdjTicux9wGXJDMezJwJXHithm4EFhf+SxLkrrKvnEk\nqQ+wbxxJUkkGe0lKAYO9JKWAwV6SUsBgL0kpYLCXpBQw2EtSChjsJSkFDPaSlAIGe0lKAYO9JKWA\nwV6SUsBgL0kpYLCXpBQw2EtSChjsJSkFDPaSlAIGe0lKAYO9JKVAqWA/HVhLPFw83yXAk8BC4Jt5\n6ZcBS4ElwOkVyqMkqZv6l5h+M3A9cEte2tuBs4EjgCZgWJI+ATg/eR0J3A+MB5ormF9JUhlK1ewf\nAtYVpF0EfIMI9AAvJq/nADOS9BXAMmByRXIpSeqWctrsxwEnA48CWeCYJH0EsDpvvtVEDV+S1MtK\nNeO095khwPHAscCdwNh25s0VS6yvr39tPJPJkMlkysiGJFWvbDZLNput2PJqOjFPLTATmJi8vwe4\nBngweb+MCPyfTN5fk7zeC0wDHitYXi6XK3oMkCS1o6amBjoXs4sqpxnnLuAdyfh4YADwEnA3MCV5\nP4Zo7pldbsYkSZVTqhlnBnAK8AZgFfAvxOWY04nLMbcBFyTzLiaadBYD24GLaacZR5LUs8r+S9AN\nNuNIUhf1RjOOJKmPMdhLUgoY7CUpBQz2kpQCBntJSgGDvSSlgMFeklLAYC9JKWCwl6QUMNhLUgoY\n7CUpBQz2kpQCBntJSgGDvSSlgMFeklLAYC9JKWCwl6QUMNhLUgoY7CUpBUoF++nAWuLh4i3qgdXA\n3GQ4M0mvBbbkpd9YwXxKkrqhf4npNwPXA7fkpeWA65Kh0DLgqMpkTZJUKaVq9g8B64qkl/2Ec0lS\nzyu3zf4S4AngJmC/vPQxRBNOFjipWzmTJFVMqWacYr4PXJmMXwV8B/gEsAYYTfwTmATcBRwONBYu\noL6+/rXxTCZDJpMpIxuSVL2y2SzZbLZiy+tMc0wtMBOY2MVpDwCXAnMK0nO5XK7zOZQkUVNTA91o\nQi+nGWd43vj7aL1SZyjQLxkfC4wDni43Y5KkyinVjDMDOIUI5KuAaUAGqCOuylkOXJjMezLRvNME\nNCfp6yueY0lSl/XGVTU240hSF/VGM44kqY8x2EtSChjsJSkFDPaSlAIGe0lKAYO9JKWAwV6SUsBg\nL0kpYLCXpBQw2EtSChjsJSkFDPaSlAIGe0lKAYO9JKWAwV6SUsBgL0kpYLCXpBQw2EtSChjsJSkF\nSgX76cBaYEFeWj2wGpibDGflTbsMWAosAU6vWC4lSd1S6uG1bwM2AbcAE5O0aUAjcF3BvBOA24Fj\ngZHA/cB4oLlgPh84LkldtLMfOP4QsK7Y9xZJOweYATQBK4BlwORyMyZJqpxy2+wvAZ4AbgL2S9JG\nEM07LVYTNXxJUi/rX8Znvg9cmYxfBXwH+EQ78xZtr6mvr39tPJPJkMlkysiGJFWvbDZLNput2PI6\n0/5TC8yktc2+vWlfTdKuSV7vJdr3Hyv4jG32ktRFO7vNvpjheePvo/VKnbuBKcAAYAwwDphdbsYk\nSZVTqhlnBnAKMBRYRdTUM0Ad0USzHLgwmXcxcGfyuh24mHaacSRJPavsvwTdYDOOJHVRbzTjSJL6\nGIO9JKWAwV6SUsBgL0kpYLCXpBQw2EtSChjsJSkFDPaSlAIGe0lKAYO9JKWAwV6SUsBgL0kpYLCX\npBQw2EtSChjsJSkFDPaSlAIGe0lKAYO9JKWAwV6SUqBUsJ8OrAUWFJl2KdAM7J+8rwW2AHOT4cbK\nZFGS1F39S0y/GbgeuKUgfTRwGrCyIH0ZcFRlsiZJqpRSNfuHgHVF0q8Dvlz57EiSdoZy2uzPAVYD\n84tMG0M04WSBk8rPliSpkko14xTaE7icaMJpUZO8riGad9YBk4C7gMOBxsKF1NfXvzaeyWTIZDJd\nzIYkVbdsNks2m63Y8mpKz0ItMBOYmAz3A68k00YBzwKTgRcKPvcAcRJ3TkF6LpfLlZldSUqnmpoa\n6FzMLqqrNfsFwIF575cDRwN/A4YStfodwFhgHPB0uRmTJFVOqTb7GcAjwHhgFTC1YHp+Ff1k4Ami\nzf7nwIXA+spkU5LUHWX/JegGm3EkqYu624zjHbSSlAIGe0lKAYO9JKWAwV6SUsBgL0kpYLCXpBQw\n2EtSChjsJSkFDPaSlAIGe0lKAYO9JKWAwV6SUsBgL0kpYLCXpBQw2EtSChjsJSkFDPaSlAIGe0lK\nAYO9JKVAqWA/HVgLLCgy7VKgGdg/L+0yYCmwBDi9EhmUJHVfqWB/M3BmkfTRwGnAyry0CcD5yeuZ\nwI2dWL4kqQeUCsYPAeuKpF8HfLkg7RxgBtAErACWAZO7mT9JUgWUU/M+B1gNzC9IH5Gkt1gNjCwz\nX5KkCurfxfn3BC4nmnBa1HQwf65YYn19/WvjmUyGTCbTxWxIUnXLZrNks9mKLa+jQN2iFpgJTEyG\n+4FXkmmjgGeB44CpSdo1yeu9wDTgsYLl5XK5oscASVI7ampqoHMxu6iuNuMsAA4ExiTDamASccXO\n3cAUYEAybRwwu9yMSZIqp1SwnwE8AowHVtFae2+RX0VfDNyZvN4DXEw7zTiSpJ5V9l+Cbqi6Zpxc\nDpqboV+/3s6JpGrV0804KvAf/wGjRsGgQfDxj/d2biSpOIN9N/zwh/C978Hvfw/r18Mf/wgzZ/Z2\nriTp9WzG6YI5c+Dv/g6+/GU46CC48EJ45BE4+OCYft99cNFFsGQJ9O/qRa2S1AGbcXrAnDnwnvdE\noP/kJ6P2/sUvwq23tgZ6gNNOiyadGTN6L6+SVIw1+074whfgueci4H/kIx3P+8ADcN55MVx/fdsa\n/vLlUFMDe+4Js2fH8iSpM7pbszfY53npJRg6tG1aLgfjxsHPfw5HHdW55SxZEv8ALrgA/uEf4JVX\n4LOfhf/5n5g+cCBs3AjZLEyaVNEi7LKeegr22w8OOKC3cyLtWnI5WLkStmyBbdugqQkGD4Zhw2Df\nfWG33eDll2Ho0O4Fe1uWEw0N8OY3w9SpcOONsPvukb5kCWzdCnV1nV/WoYfGidv3vAfe8Aa48kqY\nODFq9s89B8uWwdq18IlPxEndvfaKA8JPfxrnAt773ljBPampKX6DffaBkV3o0WjbNhgwoPX9mjVw\nxx2w//7wsY9F2q23wuc+Bzt2wNe+Bv/0T5XNe0d27Ih1uN9+XSuXKiuXg9Wr4ckn4cUX4dln4dFH\n4cQTYz/YY4/YjmqSULZtW2yP/frF5w49NLbNbDaaUQcOjM/W1UX688/D009HZe2YYyJt3brYj/fe\nu3fKu3lz7McbNsS5vQEDInivWgULFsS2OXAg/O//RtqQITFPv37Q2BiVz02bWtO6K1U1+/XrY6cv\n5pJLYkObPx/OPRc+//lIv+KKWFnXX9/177v11vjcBRfAxRe3bsgQG8PUqbB4cQT7P/8Z3vGO2LC3\nbYNvfxvOOKPr39kV69bBddfFdz/+ePw269bBlCnRXHXCCe1/NpeD7343TlYfcwysWBEb5ObN8fs9\n/njM89RTMHw43Htv1FZOPRWOOy6+Y/z4tuc8uiOXg3nzYOlSyGTiH8Ts2fDBD8YOv25dHICmTImT\n6MOHt35282aYPh1Gj4Yzz4zLaCvhqadg7lx4+GG4/3544YUoby4HZ58dv+/48XGeZ+lS+P73ozJw\nxRWx448cCYsWRa1v9Gg47LD2T/yvXx+/b34lYf16uP32+P6pU2HEiPiNjjkmKhX9+sW2l8tF4Mnl\nIp8jRsS2sH59VFZGjYrlPf10NFOuWgVjx8Jb3hL5GTu2NaBu3w5PPBG/+axZsT8tXhxBfq+9YMIE\nOPDACHonnhgXNfzyl1Gh2nffWF9HHhmXNG/eHPvMyJHxW27cGL/Z2WfH9zz8cCx706aoBR98cPzG\nc+e2/g47dsT3DRsW98KMGxfBc7fdovx1dVEJ69cvDiivvhrzDxkSAfeEE2J7WLo0DlJjxkSlMJeL\neTdtiorS8OGR17Vr42D0gx/EuoP4/Eknxfdv2BDlmTgxfrtXX4Wjj4Z3v7t4Ba+pKX6bPfaA/v1t\nxilpyxb41KfixOk998DpBY9VeeGFWNGLFsUGfvLJMT5wYGxAjz1WuaCUb/t2+PGPY2c67rjYwHI5\n+N3vohkok4kd/bLLorZfCY8+GgefurrYcc89N4aJE+GNb4zf4oYb4kA1YUJ87xFHRIDo3z92jF/8\nAv7yl9gBf/azqIGNHx8b+0EHxYa5cWMElsmT2wbPTZvin878+XHi+8gj4V3vit9i/PjYgfffP/J0\n6KHxmVdeiaGliW3LljiYLF8eNcRVq6IsW7fC4YfHv6XBg2O+//qvWFYuF993002xHbzznRHUHnww\nlnPaaXFAmDs38nHAAfEb1dXFdnDDDbGjDh0aO+BLL0XwPesseOaZ2G62bIka5cyZ8bd7+/bYyY88\nMnbmgw6K79q6FX796/gN/vrX+K322itquHvuGZf0btkSAbS5OcrUUtaWoDJ4cORzxYpYxtatsX5q\na2N9Adx9d2zrJ5wQ67OxMfL8yCOx/N12izytWRMVjMGDYztvbIxh330jeA0a1LpOTz01vqOhIfaR\n5uaooJx1VuTzzjtj2TU1UXk5+ujYjg47LNZrMRs2RFnXrImyr1kT28QHPtC2gtRZ27fHgWKffWK9\nP/10rK+amvhXvX17pG/dGhWC0aOjHKtWRVnXro1tYeDA+K0gDmgHHhhl3rAhyjhgQKw3iIPjoEHx\nW5xxBnz4w/HPvpL/0G2zL+HFFyO4DRkSG+NTT8FPftJ2no99LHbub30r3n/pS7EDH3AA/O1vcNtt\nPZbd1yxZAr/5TexYn/88/OhHsQOUo7k5aji33ho18e99L8p12mlRSylm61a4+ebYGR5/PNIuvxyu\nuirOXYweDdOmRRAo19at8KtfRXDu1y92pJNOisB5xx0R2PbeO4LibrtFQDnkEPjP/4xA11IrHjky\nDkbHHRc7dFNT7OCHHFL87+/f/ha/7Zo1EZAOPjhqsBA7+vLlMW3evBiefz7++Y0dG9vTgAEx/333\nxUHvkENi+xo0KILKGWfAm94UeevMJbgbN0ZgGTiwNW3Hjih3Sw0Q4nd54YUYf+mlCFwtNc1hwyJt\n1aoIUNu3R8AZNuz137dtW/yemzZFOUeNijI991wsL18uF4G9sTEOWsWC14svRgVg5cpYf+9+d3lB\neleUy7UtS3Nz/G577dW6beVy0SzT3BwHtfz1WEkG+w6sWAHHHhtNEtdeGzvDhAmxgbcEqQcfhI9+\nNP4OtvwV3bIlanRDhkQtrdgO05P++Meo5dx0UwT8hQvhoYeipjpoUASsefPictAhQ1o/19gYZb/n\nnvhbfeSR8feypebXWblcnFy+8cZYRstBcWdqaopAunVr5LemJpokFiyASy+NACulicG+A9/8ZtTS\nfvCD1rTTT4+a/kc+EjWcujr4+tfh/e9v+9mXX47gv7OO0l31f/8H//iPUcsaPjzKMWlSHNAefjhq\nZPfeG7XPV1+NJoo//zkOEtdeGwF78ODeLoWkchnsOzB5Mlx9dTRXtLjvvrgMcuHCOAn6hz/EX/q+\n8rdz8+Zo1y2W34aG1jbWlSuj2Wr8+J7Po6TKM9i345lnom35+edbL6OML492xcMPj/bi2bOjNixJ\nuzK7S2jHt74VzTX5gR6iRnzHHdFmf/XVBnpJ6VCVNfulS+Gtb41rewvviJWkvsiafYFcLm6a+cpX\nDPSS1KLPB/uVK9u+v+22uB655Q5YSVLpYD+deJj4gry0q4AngHnALGB0kl4LbAHmJsONlcxoMa++\nGleb5Af8H/4Q6uvtT16S8pUK9jcDZxakXQscCdQBdwHT8qYtA45KhosrlMd2/fWvca38Aw/E+2ee\niXb6wu4QJCntSgX7h4B1BWmNeeN7Ay9VKjOPPtq1PmAWLYqrbVqC/U03xc1R+b0wSpLK7+L4auCj\nwCvA8XnpY4gmnA3AFcAfu7LQmTOj171t26JrgxEjOp5/0aJ4SMisWdHR2axZ0YmYJKmtck/Q/jPw\nRuDHwL8naWuI9vujgC8CtwP7dGWhs2ZFB1C/+110ALZ+fcfzL1oUPRoOGBAdP82b19pToiSpVXdP\nY94OtNSltyUDwBygARiXjLdRX1//2ngmkyGTybBhQwTv88+PLn2bmqK3xfyuDgotXBh9aj/55K7T\nh40kVUI2myWbzVZseZ25QL8WmAlMTN6PA5Ym45cAk4kmnaFE+/4OYCzwB+AtQGH9vOhNVb/5TTwM\nY8qUaJI54oh4wPfXvlY8U6+8Et3Mbtz4+rtkJanadPemqlI1+xnAKUQgX0VcefMu4M1EUG8ALkrm\nPRm4EmgCmoELeX2gb9fixdF97lvfGg8CuPzy1/c7n2/OnKjVG+glqbRSwf5DRdKmtzPvr5KhLA0N\nEewnTIjuDpqa4lF+uVz0a37DDfD3f9/arPOnP3X82DxJUqtd5g7ahobWR/8NHdr6lKJly+LJUfPn\nx9U6LR55xGAvSZ21Swb7FkcfHbX6efPgM5+JGj9Ebd+avSR1Xq8H+3vuiZup1qyJ53bmO/roeDDz\nHnvEk+hbgv2KFdFVceH8kqTiej3Y/+xn8OlPx0OjC0+2TpoEd90Vjw4cMyae3N7UBL//fTx3ta88\nXUqSeluvB/vGRnjiieIPEZk0qfU5sQMGxAFh+XL47W/jCfaSpM7p9WC/cSMMG/b69nqAAw+MLhPq\n6uL9uHGwYAFks3BmYfdskqR29XpHwBs3wnXXRft8Md/9Lpx6aoyPGxePEqyrg/3377k8SlJf1+s1\n+8bGaK457LDi0887D4YMifETTogg39HNVpKk1+v1Z9COGhWXUY4e3cEnJCnl+vwzaBsbYfDg3s6F\nJFW3Xq3ZNzfH5ZbbtkG/fr2QE0nqI/p0zX7z5rhhykAvSTtXrwZ7m3AkqWf0arDfuBH26dKzrCRJ\n5ej1YG/NXpJ2PptxJCkFer1mbzOOJO18vR7srdlL0s5nM44kpUCpYD8dWAssyEu7CngCmAfMAvI7\nOrgMWAosAU4v9eU240hSzygV7G8GCjsTvhY4EqgD7gKmJekTgPOT1zOBG0stvxpr9tlstrezsFNZ\nvr6tmstXzWWrhFLB/iFgXUFaY9743sBLyfg5wAygCVgBLAMmd7Twamyzr/YNzvL1bdVcvmouWyWU\n25/91cBHgS20BvQRwKN586wGRna0EJtxJKlnlHuC9p+BNxLNPN/tYL5cscQpU+C974X774d99y0z\nB5Kkiqql7QnafG8EFibjX02GFvcCxxX5zDLiIODg4ODg0PlhGTtZLW2D/bi88UuAnybjE4grdAYA\nY4AGeqcLZUlSF80A1gDbgFXAx4FfEMF/HvBL4IC8+S8njj5LgDN6NKeSJEmSes6ZRK1/KfCVXs5L\npawA5gNzgdlJ2v7AfcBTwO+B/XolZ11X7Ca6jsrSpZvodgHFyldPXDk2NxnOypvW18o3GngAWESc\nS/tskl4t67C98tXT99fhIOAxosVkMfCNJL1Prrt+RBNPLbA7UajDejNDFbKcWCH5rgW+nIx/Bbim\nR3NUvrcBR9E2GLZXlpZzNLsT63QZu8AzjUsoVr5pwBeLzNsXy3cQcbMjxD0wfyX2sWpZh+2Vr1rW\n4Z7Ja3/iMvaTqOC668mCTyYytIK48eoO4kasalB4Ivps4CfJ+E+Ac3s2O2UrdhNde2Xp8k10u4Bi\n5YPiFxL0xfI9TwQAgE3Ak8S9LtWyDtsrH1THOnwleR1AVI7XUcF115PBfiRxkrdFyZuu+ogccD/w\nZ+BTSdqBRHMByeuBvZCvSmmvLCOIddiiL6/PS4j+nm6i9W9yXy9fLfEv5jGqcx3WEuVruZGzGtbh\nbsTBbC2tzVUVW3c9GexzPfhdPelEYqM7C/gM0VSQr+Ua2WpQqix9sZzfJy4VrgOeA77Twbx9pXx7\nE1fKfY623ZtAdazDvYmrAj9H1PCrZR02E2UYBZwMvL1gerfWXU8G+2dp20PmaNoemfqq55LXF4Ff\nE3+l1hLtiwDDgRd6IV+V0l5ZCtfnqCStr3mB1p3ov2n9K9xXy7c7Eeh/SnRUCNW1DlvKdyut5au2\ndbgB+C1wNH103fUnbrSqJdqkquEE7Z5AS+8+ewEPE2fFr6X1aqOv0ndO0MLrb6Jrryx99Sa6WtqW\nb3je+BeA25Pxvli+GuAW4N8L0qtlHbZXvmpYh0NpbX7aA/gD8E768Lo7iziDvoy4bKivG0P84POI\nS8FayrQ/0Y7f1y69LLyJbiodl6Wv3URX7CbBW4hLZ58gaor551f6WvlOIpoC5tF6GeKZVM86LFa+\ns6iOdTgRmEOUbT7wpSS9WtadJEmSJEmSJEmSJEmSJEmSJEmSJEnqqn2Bi0rM83CJ6ZvK+N4vJJ9r\nuV3/ROAvwEfKWJYkqYRa2vZ306KGzvcXUtg7ZGecQvRV8oW8tA+WsRxJUifcQTz4YS5wJ9H/0k+I\n/opaegZsqbnfRTx7YCGtzx+A1mC/F9Hj4DziANJR8D6PqNW3PI5yH6IvGknSTvAmWmv2bwJ28Pqn\n9rQE8yHJ6x7JZ4YUTP8A8KO8zw3u4HtbDgT3AYcSNf1hXcm41FN25ecxSp1VUzC+ktbadqHPEbX2\nPxG1/nEF0+cDpxHNMycBGzvx/bcBHyZ6W3yx07mWepDBXtVoczvppxB9hB9PPBFoLjCoYJ6ltD6U\n/OvA19pZ1kFEd8kQD9N4P7tYf+JSvv69nQGpAhppfYhMR/YlHuL8KtHscnyReYYn89xGPDHoE+0s\n61hgVt73L8QmHO3CDPaqBi8Tl1YuAJ6k+LM4c8C9wKeBxcRJ3D8VTAc4gng6UDPxkJNil3S+A6gH\nBhLPQoV4TN7L3SiDJEmSJEmSJEmSJEmSJEmSJEmSJEmqlP8HFGTWuOfvGV4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f73b7f47c10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAAERCAYAAACaUQc3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHkxJREFUeJzt3XmUFOW9//H3KIIIKCCIiiioEJfg+hNNXBiMGrzJ1Vxz\nNXpi4sHkZ4y5V6/xxuWayCQxiVFyf1nQLG4RoyQmLnGJYExshUTFRFyRCCMiuCCaYVOGtX9/fGvs\nnmYGemZ6pqe73q9z+nR1dXX3U13dn3rqqaeqQJIkSZIkSZIkSZIkSZIkSSUyDHgEeBF4ATg/GT8Q\n+CPwMvAQ0D/vNZcB84C5wAldVlJJUrvtDByUDPcF/gHsC1wNXJyMvwS4KhneD3gG2AYYDswHtuqi\nskqSSuQe4Dii1j4kGbdz8hiiVn9J3vTTgCO6rHSSpBa1pdY9HDgYeJII+iXJ+CXkgn9XYHHeaxYD\nQztWRElSRxUb9n2BO4ELgJUFz2WTW2s295wkqQv0KGKabYigv5VoxoGoze8MvAXsArydjH+d2Knb\nZLdkXDN77bVXtr6+vp1FlqTUqgf2bs8Lt1SzrwFuBOYAP8wbfy9wVjJ8FrmVwL3A6UBPYAQwEpi1\nSWnr68lms1V7mzhxYtnL4Lw5f85f9d2AvdoT9LDlmv2RwJnAc8DsZNxlRO+bO4AvAK8CpyXPzUnG\nzwHWA+dhM44kld2Wwn4mrdf+j2tl/HeTmySpm7APfCeora0tdxE6TTXPGzh/la7a568jasr0udmk\n/UmSVKSamhpoZ25bs5ekFDDsJSkFDHtJSgHDXpJSwLCXpBQw7CUpBQx7SUqBbhH2q1bBqaeWuxSS\nVL26xUFVr74Ke+8Na9bA1luXqUSS1M1V/EFVq1fDhg3wzjvlLokkVaduEfaNjXH/5pvlLYckVatu\nEfarV8e9YS9JnaNbhH1Tzf6tt8pbDkmqVt0i7K3ZS1Ln6hZhb5u9JHWubhH2q1dDz56GvSR1lm4R\n9o2NMHy4YS9JnaVbhP3q1TBihGEvSZ2lW4R9Y2Mu7Kv1aoXr1kFdHRx/PKxYUe7SSEqbbhP2gwbB\nsGHw5JPlLk3nmDYN7roLBg6EL38ZbrwRJk6E2bPLXTJJadAtwn71ath2WzjnHPjpT8tdms7xpz/B\n6afD9ddHzf7RR2O+x4+HSZNyWzQLF8JFF1XvFo6k8ijmhDo3AZ8A3gZGJ+MOBH4G9AFeBT4LrASG\nAy8Bc5PpHgfOa+E9m50I7atfhaFD4ayz4oRoCxfCDju0eV66tdGj4YYb4PDDm49fvBhOPBF23hk+\n+Un4zW/ghRdipfCZz5SnrJK6p84+EdrNwPiCcTcAFwMHAHcDX8t7bj5wcHJrKeg3sXo19O4dTTlD\nh8JrrxXzqvLLZmH9+paf27gxmmtOPRWmT4dFi+DQQzedbrfd4PHHY6vmhRfg4IPhnnvgv/8bbr4Z\nPvc5OO88eOWVLZdn7Vp45JFYafz7v0cZJAmKC/sZQEPBuJHJeICHgU93pBCNjdGMA1Gjr5QdmLfc\nAvvtB/PmxdZIvssvj/A+6CC49FI44wzo0aPl9+nbN1YK118P114Lxx4LP/wh/O53cMgh8Z2MHw/X\nXBO3lStzr81m4S9/gW9/O7qvXnxx7AR+8024+uqYZvnyTpn9LVq7tjyfK2lTrcTPFr0InAz8HjgV\nGJb33AhgNrAc+Dowc0tv1lSzB9h++/KFE0Rbeu/eMGbMlqe97TYYOTICvWfPqMl/5zvwj39Ek82c\nOTB4cAR/W33603Fr0r8/vPRSrBj32w+++EX44x9hyZJYiYwfD3/4Q5QF4JRTYnjHHWPL4K674F//\nte3laI+XXoLJk+EXv4D/+q9YETWtzCWVR3vD/mzgx8A3gHuBpjrcG0TwNwCHAPcA+xPt+c3U1dV9\nMLx4cS3bblsLlLdm/+CDsRN1993h2Wdhq7ztng0bml9YZelSmDUratC9e8fjE0+MWvqDD0ZtfvDg\n0pXtkktyw9Onwx13wBVXwE47wYEHQk1BK96wYXD++dE8dPXVMGFCrACuuabj+0PWrYPrroOjj44t\nj3z33Qdnnw1f+AK8+GKUe8894Sc/iZXX+vVxq6Twf++9uO/ZM34DW3Vxt4ZsNipE220Xj9esgalT\nozPDypXwkY/E/q4xY6JJ8K9/jX1D69bBuHGV9V2ruUwmQyaTKcl7FdvQPxy4j9wO2nyjgFuBw1t4\n7hHgIuDpgvHNdtB+/ONw4YVROz3nnGjb/tKXiixZCX3sY1E7//73I0BHj4YLLoBly6I2ff/9MGNG\n/Hn+/vfYAvnNb3Kvnzcvwm/CBPjf/2292aarrF4dWyrjx8cZRS+8EHr1gl/+svj3WLMmmmP69Yst\nhy9+Mba+dtop5nfQoGiC2muvaMq69toIomOOyb3Hk0/CaafFTuiFCyMsp06FsWNLPssl1dgIV14Z\nK7bVqyM8d989ektNmJALX4hAfvDBWOmOGxff9623RsXla1+L76utsll44IHoovv88/G7GzQIGhoi\nzC+4IPZxTZ8eTYovvwz77hvf69PJP27evFh2H/1ojO/XL8py0EGxtajNa9rvtdVWsTzefTcqd336\nNJ9u3TpYsCAuwNS7d+x7W7cOttlm8+/f0ABz58bW8Ny5cWtsjOW6997xO/rzn+O31r8/PPxw+3fQ\ntjfsBwNLiTb/XwJ/Tu4HEbX6DcCewGPAh4FlBe/XLOzHjoVvfhNqa2PH5E47RdvzSy9FYI4c2fYZ\na6v162HAgAijF1+MP+odd8QlE3/606iZ9u8fAX/oodE88t3vwi67NH+fVauidt8dvfde/Ag//OFo\nWmpph3G+bDZq4088EWF9++0wZUo8d/zxsbXzzDPw85/DP/8Jo0bBhz4UNftCK1bE9zp4MNTXR030\nK1+J+913L/28Qiyrxx+P5XHkkZtu/bTk/fejGa5vX/i3f4vw/MEPYksJ4ru45hqYOTN2nvfsCb/9\nbdT4t946wvjZZ+NPfuaZERY//3l8l7vvHt/ZWWfF99rYGOHbt2/c9+oVv78//AF+//toBhwyJP4b\nn/pU1OKXLs2tbAtls5vOY319rLAffDC+/1Wr4I03YiXwi1/AJz5RfKVkzZoo34gRWw6xclq+HB5+\nOJbRa6/FivHzn4/KRqHf/x5+9rNYLqNGxf/6xhtj2axfn9s/1qdPfL89euSOC9p553jdrrvCc8/F\nc0OGxNZ+r17RqWLMmKgMvftuPG5oiMrjihWxLN9/H/bZJ35n++wTtz594ntesCB+G+PHx+c0NMD4\n8Z0b9lOBsUSQLwEmAn2BryTP3wn8TzJ8CvAtYB2wEbgCeKCF92wW9mPGRJgefni0765ZEzWqc86J\nH/WVV7Zjztpo9mz47GfjD9bktNNi7XrzzVHDOvts+Na34KSTOr88naWhIWrVdXVw990Rgi2ZPDl6\nBS1fHtPOnBm1+r32Kk05Fi6MFc706fEDz2aj6euTn+z4e7/2WmxZTZkSNdglS+KPs9tusRO8aYsi\nk4Ef/zhqT2PGxJ/sZz+LZq7Fi+FHP4Jzz215JfHyy7FV19gYKwWILcFeveKPDREckGu6WrgwftuX\nXx5//N69I0yabmvWRHAcd1yE+yGHRGWimJVUWz32GPzHf8R8XnBBLNuhQ3PP19dH6L35ZoRTz55x\nYGC/frGyOeCA+A7eeivmef/9436ffeLAwe22i0rahz4U89BZW7mNjbECnjkzV9bp02N5HnMM7LFH\n9FD77W9ji2vp0gjpESPiP19fH1vyffvGSv611yIHdtst5nn77eNzVq2K32j//rECX7gQ3n47VvCv\nvRaf01R5WrcO/va3+D08+mjsLxsyJAJ94MB4jx12iOG2Lt+OdL3sFhccP+CAqEkfeGD8webPj/A/\n8sioiU6e3BkFiO6QTbXKyZOjRnb99blppk+Hf/mXOOjpwgtbrjlVqoceih/1uedGTW+XXeCyy+KH\nfNNNEYLf/35sbQ0a1HnlWLMmamHvvx81/YED48/yuc9F80VL3/eGDXGQ2o47brp1ctdd0QQ4YULs\ns9htt1huK1fGH++886LXUt++seV4+eURUA88kGtyGTGie2+hlVJ9PVx1Fdx5Z6xgRo+OrbQbbogt\nk333jcBbsSL+C0OGRMeApUujJjx0aHxXTZ0H5s6NCsJ778W4+vrYCpgwIaYfOTL+7z17xuc3NMRK\nZNttYyXXr9/my7t2bTSlPvporKyffjq2VI85Jv7L228fuVFYKWloiMrNsGER0osWxXI+6aRc55BK\nUPFhP3Jk/NlGjYr25EceifsBA6Km96tflb4A990Xm1fLlsUP7ZRTorfKhAn5hYyaz7Bhrb9PJZsx\nI/4A++wTNb2HHorQPfzwWPntu2/XlmfBggjlHj1iR/mxx0b7+KxZ0ay3/fYRQLfdFiuB11+PWvqu\nu8YK6Z13IrTuvrv1Jqq1a2PlsmxZ/NHTEOjFWLEiduSvWBEr0aOOin1YpfDEE7GVuGBBrABefRVO\nPjlW6tOnx47+99+PHcsjRsAJJ0RwDxwYFZE99oia++uvw1NPxeNx46Ii8tGPpmsZVnzYDxsWfcV3\n3z3+qLfcEjvFhg6NNsX77y/th2/YELWLRYuiLXPevGgqeuqpWMGkUTYbbbp77tl8x2O5LFkSPXnu\nvz9+F5MmxSbz5MmxEjj33Jhm2rSoib79dqworrgiap/qvt56K1bKffpED7am5bVhQzR/ZDIR8u+8\nE1uc9fXRHHfAAVEh3H//sha/rCo+7AcNijX+4MGxef6d70STwkknRQ1t5hZ76rfNtGnwjW9El7UB\nA6LJaObM6HEjSd1VR8K+zJ0DQ+ERtMuXRy1zzJjczq5Suv/+OJ3AiBHRbHPccQa9pOpW9rNeNh0w\nUhj2c+bEjpZlhZ02S/B5998f+wKOPDLaCs85p7SfIUndTdnDfv362NnW1G+3ac//iy/GzpdShf3k\nydG75Kmn4vP22y/2CUyZEv1YJamalb3NfuXK2AmzalU80dgYtfvttot+r0OHRu27IwdxbNwYTTbr\n18d7/ehHcZCFJFWSim6zzz8JGsSBGdlshPtOO+Waddra13vhwuiiBdGtsH//ODKxoSH65UpSmpS9\nGSd/5yxEE8sOO+S6V/Xv3/amnEmTogth0+uuvz76Zw8datBLSqduV7OHCPum3jFtDfsZM6KZZv/9\no89u795xtF21Xu5QkopRtpr9xo1xK6zZQ+ykza/ZNxRcOmXixJZPg5zNRv/5b387zqQ5a1Zc8nDS\npNw5LiQpjcoW9t/7XhwK/89/bhrEhx4aPXFg05r9e+/FychaOl/OrFlxeoMzz4TDDovTLLz+epwW\nQZLSrGxhf8stcdTq00/Hyc7yXX997opLhWFfXx87a3/4w+aX54M4+vakk+LcKmPGxFG5Z57Z/KIj\nkpRGZQv7tWvjwKkZMzZ/CcDCsJ8/Pw6GGj06etnke+yx3EUz9tgjdtKedVbpyy5JlaZsYf/5z8dJ\njR54oO1hv/fecWbGp57KjV+/Pi5UcdRR8bimJk5w1tVnbpSk7qhsvXEuvDDOfvfKK5u/EtWOO8aZ\nKTdujEuDzZ8fF3UYMiTOuf23v8XVZlavjnOX5/fH7+prhUpSd1W2OBwwAI44Imr1mwvl00+Pmv24\ncfC73+Vq9mPGxA7ZM86InbaNjXHFHUnSpsp6uoS1a6Nb5ZbOP75mTVyB6KKLcmfEHD48Llqx557R\n7l8tV5CSpNZU7OkSevYs7kITvXpFDb5Pn7hvunLUhAlxVRuDXpI2r+wnQmurd9/NXchZktKk4q9U\nJUnaso6Evf1VJCkFDHtJSoFiwv4mYAnwfN64A4HHgeeAe4F+ec9dBswD5gInlKaYkqSOKCbsbwYK\nL9x3A3AxcABwN/C1ZPx+wGeS+/HAdUV+hiSpExUTxDOAgpMMMzIZD/Aw8Olk+GRgKrAOeBWYD2zm\nZAiSpK7Q3lr3i0SwA5wKJD3f2RVYnDfdYmBoOz9DklQi7T2o6mzgx8A3iDb7tZuZtsU+lnV1dR8M\n19bWUltb286iSFJ1ymQyZDKZkrxXsf01hwP3AaNbeG4UcCtwOHBpMu6q5H4aMBF4suA19rOXpDYq\nRz/7wXmv/zrQdIXXe4HTgZ7ACKJtf1Y7P0OSVCLFNONMBcYCg4BFRE29L/CV5Pk7gV8mw3OAO5L7\n9cB5tNKMI0nqOp4uQZIqhKdLkCRtlmEvSSlg2EtSChj2kpQChr0kpYBhL0kpYNhLUgoY9pKUAoa9\nJKWAYS9JKWDYS1IKGPaSlAKGvSSlgGEvSSlg2EtSChj2kpQChr0kpYBhL0kpYNhLUgoY9pKUAoa9\nJKWAYS9JKWDYS1IKFBP2NwFLgOfzxo0BZgGzgaeAw5Lxw4HVyfjZwHWlKqgkqf1qipjmaGAVMAUY\nnYzLAN8DpgMnAhcD44iwvy9vutZks9ls20srSSlWU1MDxeX2Joqp2c8AGgrGvQnskAz3B15vz4dL\nkrpGj3a+7lJgJjCJWGF8JO+5EUQTznLg68l0kqQyam/Y3wicD9wNnEq06x8PvAEMI7YEDgHuAfYH\nVha+QV1d3QfDtbW11NbWtrMoklSdMpkMmUymJO9VbNvPcJq3xa8Ats97j2XkmnXyPQJcBDxdMN42\ne0lqo85us2/JfGBsMnws8HIyPAjYOhneExgJvNLOz5AklUgxzThTiWAfBCwCrgDOAa4FehFdLc9J\npj0G+BawDtgIfImo9UuSyqhdmwMlYDOOJLVROZpxJEkVxLCXpBQw7CUpBQx7SUoBw16SUsCwl6QU\nMOwlKQUMe0lKAcNeklLAsJekFDDsJSkFDHtJSgHDXpJSwLCXpBQw7CUpBQx7SUoBw16SUsCwl6QU\nMOwlKQUMe0lKAcNeklLAsJekFDDsJSkFign7m4AlwPN548YAs4DZwFPAYXnPXQbMA+YCJ5SmmJKk\njqgpYpqjgVXAFGB0Mi4DfA+YDpwIXAyMA/YDbifCfyjwMDAK2FjwntlsNtvBoktSutTU1EBxub2J\nYmr2M4CGgnFvAjskw/2B15Phk4GpwDrgVWA+sRUgSSqjHu183aXATGASscL4SDJ+V+CJvOkWEzV8\nSVIZtTfsbwTOB+4GTiXa9Y9vZdoW22vq6uo+GK6traW2tradRZGk6pTJZMhkMiV5r2LbfoYD95Fr\ns18BbJ/3HsuIZp1Lk3FXJffTgInAkwXvZ5u9JLVRZ7fZt2Q+MDYZPhZ4ORm+Fzgd6AmMAEYSvXYk\nSWVUTDPOVCLYBwGLgCuAc4BrgV7A6uQxwBzgjuR+PXAerTTjSJK6Trs2B0rAZhxJaqNyNONIkiqI\nYS9JKWDYS1IKGPaSlAKGvSSlgGEvSSlg2EtSChj2kpQChr0kpYBhL0kpYNhLUgoY9pKUAoa9JKWA\nYS9JKWDYS1IKGPaSlAKGvSSlgGEvSSlg2EtSChj2kpQChr0kpYBhL0kpYNhLUgoUE/Y3AUuA5/PG\n/RqYndwWJPcAw4HVec9dV6qCSpLar0cR09wM/ASYkjfu9LzhScCyvMfzgYM7XjRJUqkUE/YziBp7\nS2qA04BxpSqQJKn0OtpmfzTRxFOfN24E0YSTAY7q4PtLkkqgmJr95pwB3J73+A1gGNAAHALcA+wP\nrCx8YV1d3QfDtbW11NbWdrAoklRdMpkMmUymJO9VU+R0w4H7gNF543oAi4lQf6OV1z0CXAQ8XTA+\nm81miy+lJImamhooPreb6UgzznHASzQP+kHA1snwnsBI4JUOfIYkqQSKCfupwF+BUcAiYEIy/jPJ\nc/mOAZ4l2ux/C3yJ5j11JEll0K7NgRKwGUeS2qhczTiSpAph2EtSChj2kpQChr0kpYBhL0kpYNhL\nUgoY9pKUAoa9JKWAYS9JKWDYS1IKGPaSlAKGvSSlgGEvSSlg2EtSChj2kpQChr0kpYBhL0kpYNhL\nUgoY9pKUAoa9JKWAYS9JKWDYS1IKGPaSlALFhP1NwBLg+bxxvwZmJ7cFyX2Ty4B5wFzghNIUU5LU\nETVFTHM0sAqYAoxu4flJwDLgSmA/4HbgMGAo8DAwCthY8JpsNpttZ5ElKZ1qamqguNzeRDE1+xlA\nQ2ufDZwGTE0en5wMrwNeBeYDY9pTMElS6XS0zf5ooomnPnm8K7A47/nFRA1fklRGPTr4+jOIZpvN\nabG9pq6u7oPh2tpaamtrO1gUSaoumUyGTCZTkvcqtu1nOHAfzdvsexA190OAN5Jxlyb3VyX304CJ\nwJMF72ebvSS1UWe32bfmOOAlckEPcC9wOtATGAGMBGZ14DMkSSVQTNhPBf5K9KpZBExIxn+G3I7Z\nJnOAO5L7B4HzaKUZR5LUddq1OVACNuNIUhuVqxlHklQhDHtJSgHDXpJSwLCXpBQw7CUpBQx7SUoB\nw16SUsCwl6QUMOwlKQUMe0lKAcNeklLAsJekFDDsJSkFDHtJSgHDXpJSwLCXpBTo6AXH260uU8c3\nH/3mJuMnjp1IXW2d0zu90zu907cyfXt4pSpJqhBeqUqStFmGvSSlgGEvSSlg2EtSChQT9jcBS4Dn\nC8b/J/AS8ALw/WTccGA1MDu5XVeSUkqSOqSYsL8ZGF8wbhxwEnAA8GFgUt5z84GDk9t5JShjxclk\nMuUuQqep5nkD56/SVfv8dUQxYT8DaCgY92Xge8C65PHSUhaq0lXzD66a5w2cv0pX7fPXEe1tsx8J\nHAM8AWSA/5P33AiiCScDHNWBskmSSqS9R9D2AAYARwCHAXcAewJvAMOILYFDgHuA/YGVHS6pJKnT\nDaf5DtoHgbF5j+cDO7bwukeI0C80H8h68+bNm7c23ebTTu2t2d8DHAs8CowCegLvAoOIWv0GoqY/\nEnilhdfv3c7PlSR1kqlE88waYBEwAdgGuJWo7f8dqE2mPYXoijk7Gf+JLi6rJEmSpK4wHpgLzAMu\nKXNZSuVV4Dlia2ZWMm4g8EfgZeAhoH9ZStY+LR1Et7n5uYxYnnOBE7qojB3R0vzVAYvJHQx4Yt5z\nlTR/w4j9ZC8SW9jnJ+OrZfm1Nn91VMfy2xZ4EngGmEN0b4cKXH5bEzsXhhPNQM8A+5azQCWygFgY\n+a4GLk6GLwGu6tISdczRxAFx+WHY2vzsRyzHbYjlOp/ufwqOluZvIvDVFqattPnbGTgoGe4L/IP4\nj1XL8mtt/qpl+QFsl9z3ILq2H0WJll9XzvgYojCvEgdj/Ro4uQs/vzMVnl/6JOCWZPgW4FNdW5wO\naekgutbm52Rin846YrnOJ5Zzd9bS/EHL5wivtPl7i/jzA6wiTmcylOpZfq3NH1TH8gN4P7nvSVSQ\nGyjR8uvKsB9K7OBtspjcgqpkWeBh4G/A/03GDSGaCkjuh5ShXKXU2vzsSizHJpW8TP8TeBa4kdxm\nciXP33BiC+ZJqnP5DSfm74nkcbUsv62IFdoSck1WJVl+XRn22S78rK50JPGjOxH4CtFMkK+pf2y1\n2NL8VOK8/pQ48vsg4E3gB5uZthLmry9wJ3ABmx7QWA3Lry/wO2L+VlFdy28jMR+7EWcpGFfwfLuX\nX1eG/evEDpYmw2i+VqpUbyb3S4G7ic2oJUT7IsAuwNtlKFcptTY/hct0t2RcpXmb3J/oBnKbwpU4\nf9sQQX8rcTwMVNfya5q/X5Gbv2pafk2WAw8Ah1KBy68HUE9sfvWkOnbQbgf0S4b7AH8h9ohfTa63\n0aVU1g5a2PSI6dbmp2kHUU+iZlVP+a5r3BbDaT5/u+QNXwjcngxX2vzVAFOA/1cwvlqWX2vzVy3L\nbxC5JqjewGPAx6jQ5XcisQd9PtFlqNKNIL7sZ4iuYE3zNJBox6/ErpdNB9GtJXcQ3ebm53+I5TkX\n+HiXlrR9CufvbCJAniPafO+h+T6WSpq/o4hmgGfIdUMcT/Usv5bm70SqZ/mNBp4m5u854GvJ+GpZ\nfpIkSZIkSZIkSZIkSZIkSZIkSZIkqa12AL68hWn+soXnV7Xjcy9MXtd0yP6RxCU5z2zHe0mStmA4\nzc9306SG4s8XUniGyGKMJc5VcmHeuNPa8T6SpCL8mrjww2zgDuIcTLcQ5yxqOjNgU839HuL6Ay+Q\nuwYB5MK+D3HGwWeIFcjmwvtUolbfdEnKfsT5aCRJnWAPcjX7PYANbHrVnqYwH5Dc905eM6Dg+U8D\nv8h73fab+dymFcEfgX2Imv7gthRc6grd/XqMUrFqCoYXkqttF7qAqLU/TtT6RxY8/xxwPNE8cxSw\noojPvw34LHHGxaVFl1rqIoa9qtV7rYwfS5wj/AjiikCzgW0LpplH7qLkVwLfaOW9diZOlwxxQY1T\n6EbnE5fy9Sh3AaQSWUnuQjKbswNxEedGotnliBam2SWZ5jbiikFfaOW9DgP+lPf5L2ATjropw17V\n4l2ia+XzwEu0fC3OLDANOBeYQ+zEfbzgeYADiKsDbSQuctJSl85jgTqgF3E9VIhL5b3bgXmQJEmS\nJEmSJEmSJEmSJEmSJEmSJEml8P8BYxNORdDzZqkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f73c19f7790>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAAERCAYAAACaUQc3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHshJREFUeJzt3XmYFNW5x/HvsBoBFZgBQdABHFAUAirLDSiNC6IxcUsi\nJmpcrlfFi2i8LpiFMWokxGhMbtRHk1GJkVwjkYSo4BI7EjcSHQRkHbaAKMsAyjYLzLl/vNV2T9Mz\n3dPd0z3d9fs8Tz9dfWrpc7qq3jp96lQViIiIiIiIiIiIiIiIiIiIiIiIiIhImpQBm4HFEWl/AMq9\n11rvPWQKsApYDozLUB5FRCRFpwJDqR/sIz0A/MAbHggsBNoCxUAF0KqZ8yciImlSTOxgXwD8G+jn\nfZ4C3BExfi4wsllzJiIiCUml5n0q1sSz2vvcE9gYMX4jcFQKyxcRkTRJJdhfCjwbZxqXwvJFRCRN\n2qQw34XASRFpHwO9Iz738tLq6devn1u9enV0soiING41cGyyMydbsz8TWAZsikj7CzABaAf0AUqA\nBdEzrl69Gudc3r6mTp2a9TyofCqfH8uXz2VzzkH4/GhS4gX7mcDbQH9gA3CVl36JNy7SUuA57/1l\nYCJqxhERaRHiNeNc2kD6VQ2k/8R7iYhIC6J+8GkWCASynYVmpfLltnwuXz6XLR0KsvCdzmt/EhGR\nBBUUFEAKMVs1exERH1CwFxHxAQV7EREfULAXEfEBBXsRER9QsBcR8QEFexERH1CwFxHxAQV7EREf\nULAXEfEBBXsRER9QsBcR8YFkn1SVstJgKXf//e6D0qeOmUppoFTTa3pNr+k1fQPTJ0N3vRQRyQG6\n66WIiMSlYC8i4gMK9iIiPhAv2JcBm4HFUemTgGXAEuCnXloxsA8o916PpC2XIiKSkni9cZ4EfgXM\niEgbC3wdGAzUAkUR4yqAoenMoIiIpC5ezX4+sCMq7QbgfizQA2xNd6ZERCS9kmmzLwFOA94FgsAp\nEeP6YE04QWB0inkTEZE0SeaiqjZAZ2AkMAx4DugLbAJ6Y/8ETgJmAycAu6IXUFpa+sVwIBAgEAgk\nkQ0RkfwVDAYJBoNpW14iHfSLgTnAIO/zy8A04O/e5wpgBFAZNd8bwK3AB1HpuqhKRKSJsnFR1Wzg\ndG+4P9AOC/SFQGsvvS/W3LMm2YyJiEj6xGvGmQmMAboCG4AfYd0xy7DumDXAFd60pwE/xk7c1gHX\nATvTn2UREWkq3RtHRCQH6N44IiISl4K9iIgPKNiLiPiAgr2IiA8o2IuI+ICCvYiIDyjYi4j4gIK9\niIgPKNiLiPiAgr2IiA8o2IuI+ICCvYiIDyjYi4j4gIK9iIgPKNiLiPiAgr2IiA8o2IuI+ICCvYiI\nDyjYi4j4QLxgXwZsxh4uHmkSsAxYAvw0In0KsApYDoxLUx5FRCRFbeKMfxL4FTAjIm0s8HVgMFAL\nFHnpA4FLvPejgNeA/kBdGvMrIiJJiFeznw/siEq7AbgfC/QAW73384GZXvo6oAIYnpZciohISpJp\nsy8BTgPeBYLAKV56T2BjxHQbsRq+iIhkWbxmnIbm6QyMBIYBzwF9G5jWxUosLS39YjgQCBAIBJLI\nhohI/goGgwSDwbQtryCBaYqBOcAg7/PLwDTg797nCizw/6f3eZr3PheYCrwXtTznXMxjgIiINKCg\noAASi9kxJdOMMxs43RvuD7QDtgF/ASZ4n/tgzT0Lks2YiIikT7xmnJnAGKArsAH4EdYdswzrjlkD\nXOFNuxRr0lkK7Acm0kAzjoiIZFbSfwlSoGYcEZEmykYzjoiI5BgFexERH1CwFxHxAQV7EREfULAX\nEfEBBXsRER9QsBcR8QEFexERH1CwFxHxAQV7EREfULAXEfEBBXsRER9QsBcR8QEFexERH1CwFxHx\nAQV7EREfULAXEfEBBXsRER9QsBcR8YF4wb4M2Iw9XDykFNgIlHuv8V56MbAvIv2RNOZTRERS0CbO\n+CeBXwEzItIc8KD3ilYBDE1P1kREJF3i1eznAztipCf9hHMREcm8ZNvsJwEfAr8FjohI74M14QSB\n0SnlTERE0iZeM04sjwI/9obvAX4OXANsAnpj/wROAmYDJwC7ohdQWlr6xXAgECAQCCSRDRGR/BUM\nBgkGg2lbXiLNMcXAHGBQE8e9AdwKfBCV7pxziedQREQoKCiAFJrQk2nG6RExfCHhnjqFQGtvuC9Q\nAqxJNmMiIpI+8ZpxZgJjsEC+AZgKBIAhWK+ctcB13rSnYc07tUCdl74z7TkWEZEmy0avGjXjiIg0\nUTaacUREJMco2IuI+ICCvYiIDyjYi4j4gIK9iIgPKNiLiPiAgr2IiA8o2IuI+ICCvYiIDyjYi4j4\ngIK9iIgPKNiLiPiAgr2IiA8o2IuI+ICCvYiIDyjYi4j4gIK9iIgPKNiLiPiAgr2IiA/EC/ZlwGZg\ncURaKbARKPde50SMmwKsApYD4xpb8P79TcypiIgkLV6wfxIYH5XmgAeBod7rZS99IHCJ9z4eeKSx\n5Z96KqxalUSORUSkyeIF+/nAjhjpsZ5wfj4wE6gF1gEVwPCGFrx5M2zfnlgmRUQkNcm22U8CPgR+\nCxzhpfXEmndCNgJHNbSA3buhujrJbxcRkSZpk8Q8jwI/9obvAX4OXNPAtC5WYmlpKTt2wBNPQF1d\ngEAgkEQ2RETyVzAYJBgMpm15sZpjohUDc4BBccbd6aVN897nAlOB96LmcQcOOFq3hjlz4Lzzmppl\nERH/KSgogMRidkzJNOP0iBi+kHBPnb8AE4B2QB+gBFgQawF799q7mnFERDIjXjPOTGAMUAhswGrq\nAWAI1kSzFrjOm3Yp8Jz3vh+YSAPNOLt327uCvYhIZsQL9pfGSCtrZPqfeK9G7dlj71VV8aYUEZF0\nyMoVtKFgr5q9iEhmZCXYqxlHRCSzVLMXEfEB1exFRHxANXsRER9QsBcR8YGsNeO0aaNgLyKSKVmr\n2Xfpon72IiKZkrWafZcuqtmLiGRK1mr2Xbsq2IuIZEpWm3EU7EVEMiNrzTiq2YuIZI6acUREfEAn\naEVEfEA1exERH1A/exERH1AzjoiID6gZR0TEB+I9lrBZbN0KlZUK9iIimRKvZl8GbAYWxxh3K1AH\ndPE+FwP7gHLv9UhDC23bFtq3V7AXEcmUeDX7J4FfATOi0nsDZwHro9IrgKGJfLGCvYhI5sSr2c8H\ndsRIfxC4PZUvDgV751JZioiIJCKZE7TnAxuBRTHG9cGacILA6MYW0ro1tGoF+/cnkQMREWmSpp6g\nPRS4C2vCCSnw3jdhzTs7gJOA2cAJwK7ohZSWltqMBfDqqwHOPTfQxGyIiOS3YDBIMBhM2/IK4k9C\nMTAHGOS9XgP2euN6AR8Dw4EtUfO9gZ3E/SAq3Tmv7aZrV1ixAgoLk8m6iIh/FBQUQGIxO6am1uwX\nA90jPq8FTga2A4VYrf4A0BcoAdY0tjCdpBURyYx4bfYzgbeB/sAG4Kqo8ZGnV08DPsTa7P8IXAfs\nbGzhCvYiIpkRr2Z/aZzxfSOG/+S9EqZgLyKSGVm5XUKIgr2ISGYo2IuI+ICCvYiID2Q12B9yiO5p\nLyKSCarZi4j4QFaDfbt2UFubzRyIiPhD1oN9TU02cyAi4g8K9iIiPqBgLyLiAwr2IiI+oGAvIuID\nCvYiIj6gYC8i4gMK9iIiPqBgLyLiAwr2IiI+oGAvIuIDCvYiIj6gYC8i4gPxgn0ZsBlYHGPcrUAd\n0CUibQqwClgOjIv35Qr2IiKZES/YPwmMj5HeGzgLWB+RNhC4xHsfDzwSb/kK9iIimREv2M8HdsRI\nfxC4PSrtfGAmUAusAyqA4Y0tvG1bBXsRkUxIps3+fGAjsCgqvaeXHrIROKqxBalmLyKSGW2aOP2h\nwF1YE05IQSPTu1iJpaWlAKxZA5s3B4BAE7MhIpLfgsEgwWAwbctrLFCHFANzgEHe6zVgrzeuF/Ax\nMAK4ykub5r3PBaYC70Utzzlnx4D582HKFPjHP5LMvYiITxQUFEBiMTumpjbjLAa6A32810bgJKzH\nzl+ACUA7b1wJsKCxhakZR0QkM+IF+5nA20B/YAPh2ntIZDPNUuA57/1lYCINNOOEKNiLiGRG0n8J\nUvBFM85HH8E3vwlLl2YhFyIiOSTTzThp1ZSafW0tzJoF77xjnz//HB55pPnyJiKST3Im2H/3u3DL\nLXD//eAcXHkl3HgjbN3arFkUEckLWQ/2tbXxp5s/33rsvPoqvPsuvP46rFwJp5+unjwiIonIapt9\nZSX072/vDU8MX/kK/Pd/w7e/DUcdBcOHwymnQOvWsGULPPRQhnIuIpIled9m/9prsHMnTJgABQUw\nciT8+c/w1a/CaafBm29mJq8iIrmsxQf7Bx+0C69at7bPI0ZAjx4wZIjV7leuhO3bmz+vIiK5LKvB\nPnQjNNdIb/xFiyAQCH++6CK4916r5bdvD2ecAS+9ZE1BjS1HRMTPstpmDxbw9+6192j79kHnzrBn\nT7hmH+2pp+DZZ2HJEvj1r+HCC5sn09Ky1NbG3mZE8lVOt9lD4005a9fCMcc0HOgBzjvPeufs3Rvu\ng9/SxPv3Igerq4M77oDLLoMFC6C83NZ1ly4wfjwceij87Gfh6f/xD7j1Vpg2DV5+Gaqrs5d3kZao\nqXe9TLtQsO/Q4eBxa9ZA376Nz19YCI8/bvM/9ljz5DFZCxfCiSfaiWTn7PzDBx/AW2/BjBlWdrBr\nBd5/H0aNgk6dLK221pqwBgyAjh0PXvZTT9nB8M474UtfyliRMmLnTrjhBti4ES64AC65BHbtsua7\nhx+232/6dGvSe+EFOPJIS5s0yZrzfvITm2fgQPvtr70WSkqyXar0qKuz/eWQQxqfrqbG/hF37pyZ\nfLVUtbW2z/31r/DPf0JVlVUEqqttXL9+MGyY/UusqoJevayCefzxdm6woIn16EWL4I9/tOEBA+Dj\nj23/PuwwG373XdvHjzzS9uuOHaF3bxg92uYpKrK8lJZaJbZDB+jaNT3XE2W9Gad7d/uBunc/eMKH\nH4ZVq+B//zf+QnfsgKOPtkDR2D+BRFVXW8Dp1w++//2mz//pp7YShw2zjWryZLjrLttJ+/aFbt3g\n6afhmWfg5pttA2vb1jbKoiK47jp48UVo08ZOQocODGDz/fCHdrJ67lzrmVRWZrXdXFBXZztRrB2p\nvNyC+Pjx8MADtrEfOGDBK/qgVlVlQb6yEv7jP+z3DtmxAxYvtt/n8cftX0GPHjBvnm1ThYU2zxVX\nwNlnN32nTsX+/XaLkE8/tc4Fa9bY+7Jl8O9/w+DBdvBv3drK36GDbYfV1fC3v9mO36GDBY1vfAP6\n9LGyduli6cuWwSuv2PTDh1vX5REjrEPDYYdZHmpqrJm0Y0fb7l56yZbdty9861sW9JrL559bmY87\nrv5Bq6rKgtxrr9n+cOmlFhSPOw6OOML2hYZs3w5z5sDy5badLF1q28Ynn1jQ/drXLKB26GDn+tq3\nt+WtWGEHgQMHbP/5+GOrRC1davvtccfZttKpk223u3fba9cuW4+HHgqrV1tA7tHDbgFz9dWWp7Vr\nLSYVFlqZe/SweLB3L2zbFl7W8uWWh1atbN1WVsL111sM2LXLYlr37jBgQGrNOFkP9r17w9tv199R\nQyZPtpX+ve8ltuABA+D552HQoNQzOWGCrZR//QtmzoQxY5o2//TptgK3bIH77rMNrbragsr+/XZi\nubbWVuasWXDCCRbAH3sMvvxl2+nff9925osvtp2zVy/bYM47zy40O/54C2o33WQHhP/6L9vQMhm4\nmmL/ftshb7sNzjrLaurLl9uO+sQTFuj37IFHH7Waebrs3Gnnc6qqLLCfeKKtl7/9zSoSXbtaML35\nZgu0iXLO1stnn9m6HjLEgtMbb8Dhh9uOXlRk06xcaQHu+eftgH7kkbbNd+5swbqoyNbv0Ufb+aeR\nI23n37PHAsKKFRagRo2yfyk7d9rBYtYs2LTJ/sV89pltY3372vZVWGhlXLAA3nvPtqd9+yzf7dpZ\nUNy92w4qJ58M55xjB8Lnn7drWiZOtHzNm2cB6MorYd06W2737gcH36oqW97hh9u67NnTgtovf2lB\ntKrKttf1623c+vV2nc3gwfYvbuFCy/f//I/td/Pm2XpascJ+h5ISy+fQoZaPjz6yfFVW2u8xbpyt\ng717bbmjRtl6bZVkY/XWreHtc9cuK2+oNt6xo33etct+7+3bbZ8dMcJ+n+aQapt9NrhIffs6V1Hh\nYjrvPOdeeCH2uFi++13nfv1rG96zx7ny8sTmW7vWubFjnRszxrkVK5z761+dO/ZY5/bts+8//njn\namtt2rq6+MurrnaupMS5t99ueJpt25ybOtW5zz6rn75qlXMvvujc5s32+ZVXnAPnxo1zbvx457p0\nce6xx+rPs3+/c88+69yJJzr3y18mVuamWr/euZ07k5v3jTec+9rXnCssdG7ECOdmzbK8HnaYcyef\n7NywYc7NmOHcmjWJ/b7pVFtrv9199znXvbtz8+bFn2flSufKypwbPty51q2d69TJtp/CQhseMsS5\nwYPtc5cuznXu7NzAgc6deaZtn5980vzlStS+fc5t3Vo/bds25yZPdq5PH8v/GWfY+mvTxvaFHj2c\na9/euUsuse3tggtsP+7QwcrfpYuVt3t35/r1c+6ee5x7+WXbDhYudG7HjvB3v/eec7/5jf3ujf0u\n1dXOffCBTTtxonP33mv7xvvvO7dunXN79zbbT9RiEOcuwi1RvQIcd5xzS5ceXLDqatuAFi1K/Mf4\n059sp3POuSeesI01XvDYvt25r3zFuR/8wDbcI46wIPT66za+rs65QMC5226zjax9e+euv9652bMt\nyDrn3P/9n323c85VVdmOccEF6QlcdXXOBYOJLauiwrlu3SxwVVQ0/fu3b3fud79z7pln6qdfcYVz\nXbvass880wL0uHG2U5eWOnf33c5dc43tfHV1zn34oR2wNmxw7vHHnevVy4Ljhg3hZW7ebAeQluSV\nV2ybOfZY5wYMcO6cc2x9hlRX27ovKnLuO9+x36m21tKds7J/8kn9333jRnvlg5qa8PDOnXbguu46\nq3ysWuXc559bAF+zJnt5zGfkerAfPNiO9pEOHHBu9GgLmKEadSL27bNgvWmTcxddZDXijz5qePqy\nMuc6dnTuyivtO52zILlnT/3pliyxWsxtt1nNv7TUai4//akdVHr0sFrcW29Z7fsb36i/Y2TSmjXO\nTZhgAbakxP6lzJ/v3NVXN/wvacsW526/3WqgoVra009bjWnePAuA+/ZZEJ8927k337QDXHm5Bb8p\nU5x76CHnioudO+QQ5445xmqDRUX226xcmdGfICU1Nc4tW2br/MILnTv7bOduuMHKdPTRzp1/fvL/\ncERSQa4H+1NOcW7BgvqFWrrUAkcyNePLL7fmkcMPd+7iiy0gR7rqKueef97+qhYVHXygSdTatXZg\nOfZY55Yvt2DXs6f9/W3KAaq51NXZv5OiIjsQ3XOPBfFrr3XussvCTVy/+IUF+YkTwzXthQttnp49\nrZli1qzEvrO21rndu8Prra4u/O8nF+3a5dzDDzs3bZpzixc7989/5nZ5JLeRYrDP+gnaUaPsZOao\nUeEJnn7aelHMnNn0ha9ZYye3eve2E4BTp1p3p1atrC/2xRfbCcx+/aynQio3UVu40E4gt+Tube+/\nbyfJzj7bTuT97GfWE+jBB638y5ZZF68+fWLPX1NTvyeQiGRHqidosx7sx46FH/3I3kNuvBGOPdbu\nX5+M8nI7kz52rN0GeexYuPtu63N9zTXWHWv/futSmeyZ+ly3bRv8/vfWs6dfv2znRkTiSTXYx7uo\nqgz4KrAFCHVovAf4OvaXohK4Ens+bTGwDFjuTfcO9hzaRsW6gnbBAuv6layhQ8PDs2ZZ39Z166xL\n1uWXp6cffq4rLLSurSLiD/HqtU8C46PSpgNfBoYAs4GpEeMqgKHeK26gh4ODfVWVXdBw0kmJzB1f\nt252S+Q5c6xPtQK9iPhRvJr9fKzGHmlXxHBHYFsqGYgO9uvW2QUX6bwFwJAhdnGGbpwlIn6V7L1x\n7gMuB/YCIyPS+wDlwGfAD4C4Dw2MDvaVlXbVXrop0IuInyV7evL7wNHAU0CoP8smoDfWhPM94Fmg\nU7wFRQf7bdvs8nUREUmfVO96+Szwkjdc470APgBWAyXecD2lpaVfDG/bFqCmJvDF58pKBXsRkWAw\nSDAYTNvyEunGUwzMIdwbpwRY5Q1PAoZjTTqFwA7gANAXeBM4EdgZtbx6XS9vusluJHTzzfZ5+nRr\nX3/ggaYXRkQkXzV318uZwBgskG/Aet6cCwzAgvpq4AZv2tOAHwO1QB1wHQcH+oMUFVnTTYhq9iIi\n6Rcv2F8aI62sgWn/5L2apFs3u8ozpLLSLqgSEZH0yfr1o926WbNNiE7QioikX9aDfVFR/WCvZhwR\nkfTLerCPVbNvrie9iIj4VYsI9pEP01XNXkQk/bJ+10vn7Nman39uF1i1a2c3LNNtdUVEwpq762Wz\nKygI1+47drSnvyvQi4ikV9aDPYRP0lZXqwlHRKQ5tIhgHzpJW1OjYC8i0hxaTLDfutXuY3/yydnO\njYhI/sl6bxwI1+xfegnOPTfbuRERyT8tItj372+PD1ywwJ4ZKyIi6ZX1rpcABw7AxRfbIwnnzs1C\njkREWrhUu162iGAPUFsLe/bAEUdkIUciIi1c3gR7ERFpWKrBvkW02YuISPNSsBcR8QEFexERH1Cw\nFxHxAQV7EREfiBfsy4DNwOKItHuAD4GFwOtA74hxU4BVwHJgXPqyKSIiqYgX7J8ExkelTQe+DAwB\nZgNTvfSBwCXe+3jgkQSWn3eCwWC2s9CsVL7cls/ly+eypUO8YDwf2BGVtitiuCOwzRs+H5gJ1ALr\ngApgeOpZzC35vsGpfLktn8uXz2VLh2TvenkfcDmwj3BA7wm8GzHNRuCo5LMmIiLpkmwzy/eBo7Fm\nnl80Mp0ulRURyRHF1D9BG+loYIk3fKf3CpkLjIgxTwV2ENBLL7300ivxVwXNrJj6wb4kYngS8Dtv\neCDWQ6cd0AdYTXbuvSMiIk00E9gE1AAbgKuB57HgvxCYBXSLmP4u7OizHDg7ozkVEREREZHMGY/V\n+lcBd2Q5L+myDlgElAMLvLQuwKvASuAVIFfu0h/rIrrGypJrF9HFKl8p1nOs3HudEzEu18rXG3gD\n+Ag7l3aTl54v67Ch8pWS++vwEOA9rMVkKXC/l56T66411sRTDLTFCnV8NjOUJmuxFRJpOnC7N3wH\nMC2jOUreqcBQ6gfDhsoSOkfTFlunFbT8i+hilW8q8L0Y0+Zi+Y7ELnYEuwZmBbaP5cs6bKh8+bIO\nD/Xe22Dd2EeTxnWXyYIPxzK0Drvw6g/YhVj5IPpE9NeBp73hp4ELMpudpMW6iK6hsuTiRXSxygex\nOxLkYvk+xQIAwG5gGXatS76sw4bKB/mxDvd67+2wyvEO0rjuMhnsj8JO8obky0VXDngN+BdwrZfW\nHWsuwHvvnoV8pUtDZemJrcOQXF6fk7D7Pf2W8N/kXC9fMfYv5j3ycx0WY+ULXciZD+uwFXYw20y4\nuSpt6y6Twd5l8LsyaRS20Z0D3Ig1FUQK9ZHNB/HKkovlfBTrKjwE+AT4eSPT5kr5OmI95SZT//Ym\nkB/rsCPWK3AyVsPPl3VYh5WhF3AaMDZqfErrLpPB/mPq3yGzN/WPTLnqE+99K/AC9ldqM9a+CNAD\n2JKFfKVLQ2WJXp+9vLRcs4XwTvQbwn+Fc7V8bbFA/zvsRoWQX+swVL5nCJcv39bhZ8CLwMnk6Lpr\ng11oVYy1SeXDCdpDgU7ecAfgLeys+HTCvY3uJHdO0MLBF9E1VJZcvYiumPrl6xExfAvwrDeci+Ur\nAGYAD0Wl58s6bKh8+bAOCwk3P30JeBM4gxxed+dgZ9ArsG5Dua4P9oMvxLqChcrUBWvHz7Wul9EX\n0V1F42XJtYvoYl0kOAPrOvshVlOMPL+Sa+UbjTUFLCTcDXE8+bMOY5XvHPJjHQ4CPsDKtgi4zUvP\nl3UnIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiTXU4cEOcad6KM353Et97izdf6HL9UcD7wGVJ\nLEtEROIopv79bkIKSPx+IdF3h0zEGOxeJbdEpH0rieWIiEgC/oA9+KEceA67/9LT2P2KQncGDNXc\nZ2PPHlhC+PkDEA72HbA7Di7EDiCNBe9vYrX60OMoO2H3ohERkWZwDOGa/THAAQ5+ak8omHf23r/k\nzdM5avzFwOMR8x3WyPeGDgSvAsdhNf2ipmRcJFNa8vMYRRJVEDW8nnBtO9pkrNb+DlbrL4kavwg4\nC2ueGQ18nsD3/x74Dna3xa0J51okgxTsJR/taSB9DHaP8JHYE4HKgUOipllF+KHk9wI/bGBZR2K3\nSwZ7mMZFtLD7iYtEapPtDIikwS7CD5FpzOHYQ5yrsGaXkTGm6eFN83vsiUHXNLCsYcDrEd+/BDXh\nSAumYC/5oBLrWrkYWEbsZ3E6YC5wPbAUO4n7TtR4gMHY04HqsIecxOrSeTpQCrTHnoUK9pi8yhTK\nICIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIuny//Fb4bVz/e73AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f73b7f09850>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAAERCAYAAACaUQc3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAH8lJREFUeJzt3XmYFPWdx/H3wAwEGAeRQS5JBnCQSwSMCEGgFaOwZsUN\nCl67rLrRBFZcchBMXBmTPNEY88SND5p1FYMmoETUyLoYkdAKKocy3OdwKOfIMRwDE66p/eNbTfc0\n3TM9PT3d012f1/P009W/rq7+1dT0t371rV/VD0REREREREREREREREREREREREREJEE6AQuBdcBa\nYKJbfptbdhboH/aZh4EtwEbghuRUU0RE6qId0NedzgU2AT2A7kA3bEcQGux7AiuBHKAAKAEaJamu\nIiISRU2BeB8WvAHKgQ1AB6zVvjnC/KOAWcBpYAcW7AckoqIiIhK/2rS6C4B+wNJq5ukA7Ap5vQvo\nWPtqiYhIIsUa7HOB14GHsBZ+bTi1nF9ERBIsO4Z5coA5wB+Bt2qYdzd2UjfgEresiq5duzpbt26N\ntY4iImK2ApfG88GaWvZZwIvAeuDpauYJeBu4HWgCdAYKgWXhH9i6dSuO42TsY+rUqSmvg9ZN66f1\ny7wH0LWWMf6cmlr2g4G7gdVAsVv2E6Ap8AyQD7zjvjcS2ynMdp/PAONRGkdEJOVqCvaLid76j5bS\n+aX7EBGRBkJ94OuBz+dLdRXqTSavG2j90l2mr19dZNU8S71w3PyTiIjEKCsrC+KM22rZi4h4gIK9\niIgHKNiLiHiAgr2IiAco2IuIeICCvYiIByjYi4h4QMqC/ZkzqfpmERHvSVmw/93vUvXNIiLek7Jg\nf/Roqr5ZRMR7FOxFRDwgZcH+yJFUfbOIiPeoZS8i4gEK9iIiHlBTsO8ELATWAWuBiW75RcB8YDPw\nHnChW14AVGAjVxUDz0ZbsNI4IiLJU1OwPw1MAnoBA4EJQA9gChbsuwEL3NcBJUA/9zE+2oLVshcR\nSZ6agv0+YKU7XQ5sADoCNwMz3PIZwC21/WIFexGR5KlNzr4Aa60vBdoCpW55qfs6oDOWwvED10Rb\nmNI4IiLJU9OA4wG5wBzgIeBY2HuO+wDYg+X5y4D+2KDkvSJ8hvLyIh59FBo1snEjNXakiEhVfr8f\nv9+fkGXFMpZhDvC/wDzgabdsI+DD0jztsZO43SN8diHwA2BFWLmTl+fwxRfQsmUctRYR8aD6HIM2\nC3gRWE8w0AO8DYxzp8dhLXiAfKCxO90FKAS2RVpwXp5SOSIiyVJTGmcwcDewGsvDAzwMPAHMBu4D\ndgBj3PeGAj/DevFUAg8AhyMtuGVLnaQVEUmWmoL9YqK3/q+PUPaG+6hRXp6CvYhIsqTsCtqWLZXG\nERFJlpQFe7XsRUSSJ6XBXi17EZHkSGkaRy17EZHkUBpHRMQDFOxFRDwgZcG+eXM4cSJV3y4i4i0K\n9iIiHqBgLyLiASkL9i1aKNiLiCRLSlv2x4+n6ttFRLxFaRwREQ9QsBcR8QAFexERD1CwFxHxAAV7\nEREPqCnYd8LGkV0HrAUmuuUXAfOBzcB7wIUhn3kY2IKNU3tDtAU3awYVFVBZGV/FRUQkdjUNXNvO\nfawEcoHPgFuAe4ADwJPAj4FWwBSgJzATuAroCLwPdMOGKAzlOI5D8+Zw4IC18kVEpHr1OeD4PizQ\nA5QDG7AgfjMwwy2fge0AAEYBs7AxaHcAJcCAaAtXKkdEJDlqk7MvAPoBS4G2QKlbXuq+BugA7Ar5\nzC5s5xCRgr2ISHLUNOB4QC4wB3gIOBb2nuM+oon4XlFRESdOwJNPwq23+vD5fDFWRUTEG/x+P36/\nPyHLiiX3kwP8LzAPeNot2wj4sDRPe+wkbncsbw/whPv8LjAVOxoI5TiOQ//+8D//A1deGXf9RUQ8\noz5z9lnAi8B6goEe4G1gnDs9DngrpPx2oAnQGSgElkVbuNI4IiLJUVMaZzBwN7AaKHbLHsZa7rOB\n+7ATsWPc99a75euBM8B4qknx6M6XIiLJUVOwX0z01v/1Ucp/6T5qpJa9iEhypOwKWlCwFxFJlpQH\ne93TXkSk/qU82KtlLyJS/xTsRUQ8QMFeRMQDUhrs1fVSRCQ51LIXEfGAlAd79cYREal/KQ32eXlw\n+HAqayAi4g0pDfZdusD27TZa1d/+lsqaiIhktrjunpYAjuM4HD8O+fnwySfwjW9YSicrVTUSzzt1\nCtatg69+FVq3TnVtRM5Xl7texno/+3rRogW0agXz5tl4tAcOQJs2qayReMHRozBxIixfDjt2QHY2\nDB4M+/bZe/v3Q7t20L8/bNoE27bBt75lacfNm+HSS+Hqq6FHD7jqKmjcONVrJFKzlLbsAYYNsxb9\nZ5/BsmX24xGJxZ49FqB79YKmTWue/+hR+PWvYc4cuOYa+Pd/h86drUX/zjs2zz//s6UV166FFSug\ne3coKIDnnoNGjWDQINiwAT79FFavhi+/hHHj4JJL4JZbbEfx+edQWAhdu8Lp0/D66zB3Llx0ke0s\nTp+Gr38dRo2ysRwaNbJ16d7dpkWiqUvLPuXB/r77YPp0yM2Fl16CW29NUY2kwTl7FoqLobwcpk2D\nr30Njhyx4LxmjQXV9u2t9d2hgx0h9uwZ/PyJE5CTY489e2DkSOjbF26/HUaMSEzKcONGePllKC2F\nN96wo9V+/aCkxM5HZWfD9dfD6NG2HoWF0KwZLFwIf/2r7TDOnLEj2tOnYfhw24H17Qt9+kC3bnaE\nkZ3SY3BpKNI62D/+OPzkJ/Dtb1ve/gc/SFGNPOj0afj97+G3v4VDh6BlS2u5Tpxo2yLQSi0rs6Ov\nbt2stRpLK7o29uyBF16w72ne3NIjXbrAz39uwfQrX4EJE6xl3qqVve7a1eqYnQ2OAzNnwqRJ8B//\nYQ2In/0MXnnFdgb/8i82Itr3vgeTJ9ffeaEDByzYN2tmr8+ehb//3cqq4zhWp02bwO+3Oi9bZumj\ntWst1TRokB0FDx1qR79Nm8KxY/b3CqSRVq2ynUfr1vbdN94ITZrUz7pKaqRtzh4s/3nhhXZYvW1b\nqmuT/ioqLIjn5dnrsjJreQZ6Pd1zjwX1SZOsJTxwoKUZCgos4M+fD3feacGivNyCbsuWFlQ2b7Z5\n3nzTUhI9egTTDmfPwp/+ZC3xRx+1oBywdavN36qVBb+PP7YAduCABaznnoOxYy2AnzgBs2fD3r3W\nEp87t+ZWbVYW3HWXtYZ/9Svo1An+6Z9snT/+GBYvtiODm26qj794UH5+1deNG9cc6CG487nsMnsA\n3Hxz8P1Dh2DRIvjwQ9uZbdxoKZ/AjrBlS0sntWplO8DDh+3veN990LatbeN166x+o0fbtgj85o4d\ngwsusP8bv9+Ogho1siOTkydhwAA7Ukn0Dl6SL5Y9xHTgJuBL4HK37Arg90ALbKSqu7CByAuADdgY\ntQCfYKNVhTvXst+5036IAwdaOuftt+NaD885ftxavp99ZgF9zBgLEGPHwhdfBH/0mzYF0xcVFfBf\n/2WfHz/ejqgi/Ygdx37s7dpZIAj17ruWamvRwgLBrbfCa69ZWuKKK+Dyyy3ot2xpO/LLLrOdQ3a2\nHb0tWGBBqXt3S12cPGk7ni5dEve3OXbM0oKZ2rPr6FHbqfbrZ+taUWF/y7y8quu8Y4ftKPx+Sx/t\n22fThw/bjnDvXvv7Z2VZKmnQIPtcZaXNn5MDS5faTr5PH3t/wgTbKQccOGA7ku3bYf16S6cNG2b/\nCzr/kHj1ncYZApQDLxMM9suB7wOLgHuw8WYfxYL93JD5ojkX7AOKi+1E1+rVsVbdm86csXz1T39q\nAffee+35+ectEI8bZ2mYlSvtcL6gwE4eBpw8GWy9xauy0lryP/+5tdqHDbOWdaAVW14e3NEsW2Y7\noIoKu5bi0ktt5yOp5TjW0OrUyXYYLVpE71VUXm6Nivfes7Rfnz52pLdxowX7yy+37qq9esHu3bZD\nP3LETkJ3725pp8JC23F8+KF999VXg89nRxjvvGM7lOxsawh06mQ7ohMnbHnt2tk5kU8/tSPD3r3t\nkZ1t/3cFBdX/Px85Ykd4TZtaoycvz7p7V1QEz6Hk5dk6gf1/FxfbDrRJE/tfbtfO6ppqycjZF1A1\niB8GAqveCXgX6BVhvmjOC/ZlZbbRjhyJsUYes2WLdRWcNs0O2YcPh2eescAdEMj9itSXsjL7Pzx1\nyo7GovUg2rHDziFs2GApqM8/tyORb37T/kc/+sgC/7Fj1u114EBrQOzfbzuhSy6xINyhgwX+tm1t\n59Gihe14Nm2yWPHqq/bct6/9Lg4dsqDcvLkt6/Bh+40MGmTLX7HCdl5XXGFHn9u3WwPo2DHbyTRp\nArt2wcUXw8GDti6tW9t5pTZtrF7bttkO7847raHTsaP9PlessOkhQ2r+Hcb7W01FsP8IeBL4C9bC\nLwLy3PnWAluAI8Aj2Di24c4L9o5jG/fQoeAetrYOH7b87Le+Fd/nG6r33oO777aW0ODB8OCDOkQW\nCThxwlr9bdva49AhOwJo08YCeuPGwaOWykoLsuGB9tQpC9hnz1or/uKLLSaBzXv2rKU2d+607rqf\nfmrnwoqLLe7k5dnJ83XrbFlXX207wtJS2zlUVlrqde9eezRubD3HSkps53TddbYep05Z6rR9e9ux\nHD9uv/XBgy0u9u6d/GB/GfA7oDXwNjARyAeaYHn8MqA/8BbW4j8Wtjxn6tSp5174fD58Ph9t2lje\nr7oLq06csJzwuHHn9zR49ll47DH7Y86fb70R0lVJia3rM8/YidSZM+2fSUQalvJySxHl5NgOYtUq\na+Vv2mRBu0uX4Mn69u1tZ1JRYbGuRw9Lhc2fHzwiOXbMji727oWyMj/bt/vZudN2GAcOPAb13Iuy\nAFgT5b1uwNIo7y3Egn44J5KCAsfZti3iW47jOM6pU47Tu7fjXHyx4zz+eLD85El7b8QIx8nKcpzn\nn3cccJx9+6IvqyFbsMBx8vIcp0sXx5k82XHKylJdIxFpCAAnQjyNSbzJgEDbuxGWqnnOfZ0PBE7z\ndAEKgZg7VObm2l4ymt277ZBpyRJ46inL3QE88YTlyT76yK6A/OEPrXz58uq/r6LCTg45cf/5Euvk\nSfjgA/jud62P+Nat1pWwIZwYEpH0FkuwnwV8jKVudgL3AncAm7BulruAP7jzDgVWAcXAn4EHsJO5\nMcnNtUOYaPbutRM2nTvDiy9ar461ay1Pv3+/neS57TbrmnbXXdYTJJp9+yynduutwZ1GqpSWwgMP\nWD/oyZPh/vur9rMWEamrWC6quiNK+e8ilL3hPuJywQXVt+z37rWcF9h9RZYsgVmz7GTJqlW2s2jW\nzHL6ubmWwwcL/qFdyxzHug2OHm3ls2bZmf5kmznTuqn95S92DmLHDt1tUUTqR4Pq01FTGmffPju5\nETBypLXwc3Ot21TgUvo777S+vcuXw4wZdiQwYYLdu6R3bztD/+qrdlHRHXfYRUGVlfWzTv/3f7Yj\n+cd/tLQR2PMrr8CUKfbe0qXwm98o0ItI/Un57RJC1RTsQ1v2YH1nKyqsW1K49u2tL+xrr1kLeuJE\nC7x/+IN1eSors7RJfr4F/+XLrbtUIm3bBv/6r3bfl9des37GJ0/aehQW2pWlV16Z2O8UEYkk7YL9\ngAHB1zk51r0yWgpm4cKq0zk5wW6doUcIw4fbidFEBftVq2z5Y8bYla4332x3WVy82I5AOnfWXQxF\nJLnSKo2zd2/VIA12m4AHH6x52R06RO+/7/PZPUPCnTxZ83LDrVxpyysosMA+caKVN2liF04UFirQ\ni0jypUWwdxy7wGjfvqppHLBuiYFbysZr6FDrtnnmTLDsv//bLnAYPTr25bzzjuXmp02zXkJ//KNu\nXyAiDUODC/aRul6+8YYF5PCcfaLk59vAGIF++V9+CY88Yl03P/jALpGuzvHjlqqZNMny83feaXcG\njOX2tiIiydCggn20rpfz5llf+H377J4V9WHsWLvF8uHD1jq//347eXrjjfb91fnpT60X0Jo16X2L\nBhHJXA0q2EdL4yxYYINRtG5dfyPv/Nu/wZ//bEPIDRwIv/iFld90k/XiCRe4DeoPf2iDbTz3nAZ4\nEJGGq8EH+23b7ETpr39tLe760ratXXV73XXw9NPBXPuNN1oq54svglfsgvXbHz3a5luyRH3kRaRh\na1D9QkKDfWWl3dpzwQILwF27WpqlPk2bdn5Z69Z2r5qBA61v/ksvWet/3jwbaCUw/J+ISEPWIFv2\njmP90qdPt2A/fHhq6/XII3ZB1Btv2AVa3/mO1U2BXkTSRYNs2c+ZY6mRU6fsns9PPpnaerVoYd0q\nwUai+cY37GhDRCRdpKoXuHtr5qp27rRA2rEj/OhHdquBdu1sBJmG4uBBu72CLowSkWSry7CEDSpk\n5ebakGIHD1oa5x/+wW5u1pDoRKyIpKMG1bI/fdq6Vl5xhd12oLTUbkucn5+CGoqINDB1adnHcoJ2\nOlBK1WEJrwA+AVZjY9BeEPLew9iA4xuBG2pTmZwc66seuNlZ27YK9CIiiRBLsH8JGBFW9gIwGegD\nvAn8yC3vCYx1n0cAz8b4Hefk5la9s6WIiNRdLIF4EVAWVlbolgO8DwRuFzYKG8bwNLADKAFqFbpb\nt7Y+7SIikjjx9rNfhwV2gNuATu50B2xM2oBdQMfaLHjJEhtNSkREEifeYH8vMB74FMgFTlUz7/ln\nYqvR0HrfiIhkgni7Xm4CAvd37Abc5E7vJtjKB7jELTtPUVHRuWmfz4fP54uzKiIimcnv9+OPNLJS\nHGLtwlMAzAUud1+3AfZjRwZ/AP7mPvcEZmJ5+o5YPv9Szm/dR+x6KSIi0dX3RVWzgGFAPrATmIql\nbia478/BAj3AemC2+3wGS/UoqouIpFiDuqhKRESiq++LqkREJM0p2IuIeICCvYiIByjYi4h4gIK9\niIgHKNiLiHiAgr2IiAco2IuIeICCvYiIByjYi4h4gIK9iIgHKNiLiHiAgr2IiAco2IuIeICCvYiI\nByjYi4h4QCzBfjpQCqwJKRsALAOKgeXAVW55AVDhlhcDzyaqoiIiEr9YRjwZApQDLxMcg9YPPA78\nFRgJTAau5fyxaqPRSFUiIrVU3yNVLQLKwsr2Ai3d6QuB3fF8uYiIJEcsA45HMgVYDDyF7TAGhbzX\nGUvhHAEececTEZEUivcE7YvAROCrwCQsrw+wB+gE9AO+D8wELqhjHUVEpI7ibdkPAK53p18HXnCn\nT7kPgBXAVqDQna6iqKjo3LTP58Pn88VZFRGRzOT3+/H7/QlZVqyJ/gKqnnhdgbXoPwCGA09gPXLy\nsfz+WaAL8CHQGzgctjydoBURqaW6nKCNpWU/CxiGBfKdwKPA/cA0oCnW1fJ+d96hwM+A00Al8ADn\nB3oREUmyuPYQCaCWvYhILdV310sREUlzCvYiIh6gYC8i4gEK9iIiHqBgLyLiAQr2IiIeoGAvIuIB\nCvYiIh6gYC8i4gEK9iIiHqBgLyLiAQr2IiIeoGAvIuIBCvYiIh6gYC8i4gEK9iIiHhBLsJ8OlAJr\nQsoGAMuAYmA5NiRhwMPAFmAjcENiqikiInURy4gnQ4By4GWCY9D6gceBvwIjgcnAtUBPYCYW/DsC\n7wPdsCEKQ2mkKhGRWqrvkaoWYYOIh9oLtHSnLwR2u9OjsDFrTwM7gBLsKEBERFIolgHHI5kCLAae\nwnYYg9zyDsCSkPl2YS18ERFJoXiD/YvAROBN4DYsr//NKPNGzNcUFRWdm/b5fPh8vjirIiKSmfx+\nP36/PyHLijX3UwDMJZizPwrkhSzjMJbWmeKWPeE+vwtMBZaGLU85exGRWqrvnH0kJcAwd/o6YLM7\n/TZwO9AE6AwUYr12REQkhWJJ48zCAns+sBN4FLgfmAY0BSrc1wDrgdnu8xlgPFHSOCIikjxxHQ4k\ngNI4IiK1lIo0joiIpBEFexERD1CwFxHxAAV7EREPULAXEfEABXsREQ9QsBcR8QAFexERD1CwFxHx\nAAV7EREPULAXEfEABXsREQ9QsBcR8QAFexERD1CwFxHxAAV7EREPiCXYTwdKgTUhZa8Cxe5ju/sM\nNlZtRch7zyaqoiIiEr9YhiV8CXgGeDmk7PaQ6aewAccDSoB+da+aiIgkSizBfhHWYo8kCxgDXJuo\nComISOLVNWc/BEvxbA0p64ylcPzANXVcvoiIJEAsLfvq3AHMDHm9B+gElAH9gbeAXsCx8A8WFRWd\nm/b5fPh8vjpWRUQks/j9fvx+f0KWFeso5QXAXODykLJsYBcW1PdE+dxC4AfAirByx3Gc2GspIiJk\nZWVB7HG7irqkca4HNlA10OcDjd3pLkAhsK0O3yEiIgkQS7CfBXwMdAN2Ave45WPd90INBVZhOfs/\nAw9QtaeOiIikQFyHAwmgNI6ISC2lKo0jIiJpQsFeRMQDFOxFRDxAwV5ExAMU7EVEPEDBXkTEAxTs\nRUQ8QMFeRMQDFOxFRDxAwV5ExAMU7EVEPEDBXkTEAxTsRUQ8QMFeRMQDFOxFRDxAwV5ExANiCfbT\ngVJgTUjZq9hoVMXAdvc54GFgC7ARuCEx1RQRkbqIZcSTIUA58DJVBxwPeAobevAXQE9gJnAV0BF4\nHxvOsDLsMxqpSkSklup7pKpFQFm07wbGEByLdpQ7fRrYAZQAA+KpmIiIJE5dc/ZDsBTPVvd1B2BX\nyPu7sBa+iIikUHYdP38HlrapTsR8TVFR0blpn8+Hz+erY1VERDKL3+/H7/cnZFmx5n4KgLlUzdln\nYy33/sAet2yK+/yE+/wuMBVYGrY85exFRGqpvnP20VwPbCAY6AHeBm4HmgCdgUJgWR2+Q0REEiCW\nYD8L+BjrVbMTuMctH0vwxGzAemC2+zwPGE+UNI6IiCRPXIcDCaA0johILaUqjSMiImlCwV5ExAMU\n7EVEPEDBXkTEAxTsRUQ8QMFeRMQDFOxFRDxAwV5ExAMU7EVEPEDBXkTEAxTsRUQ8oK73s49bkb+I\nxz547LzyqcOmUuQr0vyaX/Nrfs0fZf546EZoIiJpQjdCExGRainYi4h4QCzBfjo2qPiasPIHsZGq\n1gK/cssKgAqg2H08m5BaiohIncQS7F8CRoSVXQvcDPQBegNPhbxXAvRzH+MTUMe0k6gBghuiTF43\n0Pqlu0xfv7qIJdgvAsrCyr4HPA6cdl/vT2Sl0l0m/8Nl8rqB1i/dZfr61UW8OftCYCiwBPADXw95\nrzOWwvED19ShbiIikiDx9rPPBloBA4GrsEHGuwB7gE7YkUB/4C2gF3CszjUVEZF6V0DVE7TzgGEh\nr0uA1hE+txAL+uFKAEcPPfTQQ49aPUqIU7wt+7eA64APgG5AE+AgkI+16s9iLf1CYFuEz18a5/eK\niEg9mYWlZ04CO4F7gBzgFay1/xngc+f9NtYVs9gtvynJdRURERERkWQYAWwEtgA/TnFdEmUHsBo7\nmlnmll0EzAc2A+8BF6akZvGJdBFddevzMLY9NwI3JKmOdRFp/YqAXQQvBhwZ8l46rV8n7DzZOuwI\ne6JbninbL9r6FZEZ2+8rwFJgJbAe694Oabj9GmMnFwqwNNBKoEcqK5Qg27GNEepJYLI7/WPgiaTW\nqG6GYBfEhQbDaOvTE9uOOdh2LaHh34Ij0vpNBb4fYd50W792QF93OhfYhP3GMmX7RVu/TNl+AM3d\n52ysa/s1JGj7JXPFB2CV2YFdjPUqMCqJ31+fwu9CdzMww52eAdyS3OrUSaSL6KKtzyjsnM5pbLuW\nYNu5IYu0fhD5ToLptn77sB8/QDl2O5OOZM72i7Z+kBnbD+CE+9wEayCXkaDtl8xg3xE7wRuwi+CG\nSmcO8D7wKfAdt6wtlirAfW6bgnolUrT16YBtx4B03qYPAquAFwkeJqfz+hVgRzBLycztV4Ct3xL3\ndaZsv0bYDq2UYMoqIdsvmcHeSeJ3JdNg7J9uJDABSxOECvSPzRQ1rU86rutz2JXffYG9wG+qmTcd\n1i8XmAM8xPkXNGbC9ssFXsfWr5zM2n6V2Hpcgt2l4Nqw9+PefskM9ruxEywBnai6V0pXe93n/cCb\n2GFUKZZfBGgPfJmCeiVStPUJ36aXuGXp5kuCP6IXCB4Kp+P65WCB/hXsehjIrO0XWL8/Ely/TNp+\nAUeAd4ArScPtlw1sxQ6/mpAZJ2ibAxe40y2Aj7Az4k8S7G00hfQ6QQvnXzEdbX0CJ4iaYC2rraRu\n9LPaKKDq+rUPmZ4EzHSn0239soCXgd+GlWfK9ou2fpmy/fIJpqCaAR8Cw0nT7TcSO4NegnUZSned\nsT/2SqwrWGCdLsLy+OnY9TJwEd0pghfRVbc+P8G250bgxqTWND7h63cvFkBWYznft6h6jiWd1u8a\nLA2wkmA3xBFkzvaLtH4jyZztdzmwAlu/1cCP3PJM2X4iIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIi\nIiJSWy2B79Uwz0c1vF8ex/dOcj8XuGR/MDYk591xLEtERGpQQNX73QRkEfv9QsLvEBmLYdi9SiaF\nlI2JYzkiIhKDV7GBH4qB2dg9mGZg9ywK3Bkw0HJ/Cxt/YC3BMQggGOxbYHccXIntQKoL3rdhrfrA\nkJQXYPejERGRevA1gi37rwFnOX/UnkAwb+U+N3M/0yrs/dHA8yGfy6vmewM7gvlAd6yl36Y2FRdJ\nhoY+HqNIrLLCpj8n2NoO9xDWav8Ea/UXhr2/Gvgmlp65Bjgaw/f/CbgLu+Pi/phrLZIkCvaSqY5H\nKR+G3SN8IDYiUDHwlbB5thAclPwXwH9GWVY77HbJYANqfJsGdD9xkVDZqa6ASIIcIziQTHVaYoM4\n/x1LuwyMME97d54/YSMG3RdlWVcBC0K+fy1K4UgDpWAvmeIg1rVyDbCByGNxOsC7wHeB9dhJ3E/C\n3gfog40OVIkNchKpS+d1QBHQFBsPFWyovIN1WAcRERERERERERERERERERERERERERERERERERFJ\nhP8HHXirIULKEScAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f73c19eb210>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAAERCAYAAACaUQc3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHqpJREFUeJzt3XmcFNW99/HPsAyIIiIgMjA4rHHDAKJyI0urqGiMJDHG\nLblujy8feVyiXjFoInPlIYqISdSrV6OgKOI1XmNEIxK8tmJkibI4iAjDogwou4DDMgNz7h+/aqdp\np+npfanv+/XqV3efqq4+Z6rnV6fOOXUKREREREREREREREREREREREREREQkRSYBG4CKsLQXgYXe\nY7X3HDIaWAEsA87JUB5FRCRJg4F+HBjswz0I/MZ7fTywCGgOlAGVQJM0509ERFKkjIaDfRHwBdDD\nez8auDNs+QxgYFpzJiIijZJMzXsw1sSz0ntfAlSFLa8COiexfRERSZFkgv1lwAsx1nFJbF9ERFKk\nWRKf+wnQPyxtHVAa9r6Ll3aAHj16uJUrV0Ymi4jIwa0Eeib64URr9sOAT4H1YWmvAZcCxUA3oBcw\nP/KDK1euxDlXsI8xY8ZkPQ8qn8rnx/IVctmcc1DfP5qQWMF+GvAB0BtYC1ztpV/iLQu3FHjJe34T\nGImacUREckKsZpzLoqRfHSX9d95DRERyiMbBp1ggEMh2FtJK5ctvhVy+Qi5bKhRl4Tud1/4kIiKN\nVFRUBEnEbNXsRUR8QMFeRMQHFOxFRHxAwV5ExAcU7EVEfEDBXkTEBxTsRUR8QMFeRMQHFOxFRHxA\nwV5ExAcU7EVEfEDBXkTEBxTsRUR8QMFeRMQHFOxFRHxAwV5ExAdi3ZYwbcqD5fz7u//+nfQxQ8dQ\nHijX+lpf62t9rR9l/UToTlUiInlAd6oSEZGYFOxFRHxAwV5ExAdiBftJwAagIiL9JuBTYAkw3ksr\nA3YDC73HYynLpYiIJCXWaJzJwCPAlLC0M4ALgZOAWqBD2LJKoF8qMygiIsmLVbOfDWyLSLsBuA8L\n9ACbUp0pERFJrUTa7HsBQ4C5QBAYELasG9aEEwQGJZk3ERFJkUQuqmoGtAUGAqcALwHdgfVAKXYm\n0B94FTgB2Bm5gfLy8m9fBwIBAoFAAtkQESlcwWCQYDCYsu01ZoB+GTAd6OO9fxO4H3jXe18JnAZs\nifjcO8DtwIKIdF1UJSISp2xcVPUqcKb3ujdQjAX69kBTL7071tyzKtGMiYhI6sRqxpkGDAXaAWuB\ne7DhmJOw4Zg1wL966w4B7sU6buuA64GvU59lERGJl+bGERHJA5obR0REYlKwFxHxAQV7EREfULAX\nEfEBBXsRER9QsBcR8QEFexERH1CwFxHxAQV7EREfULAXEfEBBXsRER9QsBcR8QEFexERH1CwFxHx\nAQV7EREfULAXEfEBBXsRER9QsBcR8QEFexERH1CwFxHxAQV7EREfULAXEfEBBXsRER+IFewnARuA\nioj0m4BPgSXA+LD00cAKYBlwToryKCIiSWoWY/lk4BFgSljaGcCFwElALdDBSz8euMR77gzMAnoD\ndSnMr4iIJCBWzX42sC0i7QbgPizQA2zynkcA07z0NUAlcGpKcikiIklJpM2+FzAEmAsEgQFeeglQ\nFbZeFVbDFxGRLIvVjBPtM22BgcApwEtA9yjruoYSy8vLv30dCAQIBAIJZENEpHAFg0GCwWDKtlfU\niHXKgOlAH+/9m8D9wLve+0os8P8f7/393vMMYAwwL2J7zrkGjwEiIhJFUVERNC5mNyiRZpxXgTO9\n172BYmAz8Bpwqfe+G9bcMz/RjImISOrEasaZBgwF2gFrgXuw4ZiTsOGYNcC/eusuxZp0lgL7gJFE\nacYREZHMSviUIAlqxhERiVM2mnFERCTPKNiLiPiAgr2IiA8o2IuI+ICCvYiIDyjYi4j4gIK9iIgP\nZCXYT50KNTXZ+GYREX/KykVVRx7pWL4c2rXLwreLiOShvLyoqrhYNXsRkUxSsBcR8QEFexERH1Cw\nFxHxAQV7EREfULAXEfEBBXsRER/ISrBv3lzBXkQkk1SzFxHxAQV7EREfULAXEfEBBXsRER9QsBcR\n8QEFexERH4gV7CcBG4CKsLRyoApY6D2Ge+llwO6w9MeibVTBXkQks5rFWD4ZeASYEpbmgIe8R6RK\noF+sL1WwFxHJrFg1+9nAtgbSk7rpSXEx1NYmswUREYlHom32NwGLgaeBI8LSu2FNOEFgULQPq2Yv\nIpJZsZpxGvI4cK/3eiwwEbgWWA+UYmcC/YFXgROAnZEbeP/9cvbtg7o6CAQCBAKBBLIhIlK4gsEg\nwWAwZdtrTHNMGTAd6BPnsneA24EFEelu/HjH5s3wwAONz6iIiJ9l4x60ncJe/4T6kTrtgabe6+5A\nL2BVQxtQM46ISGbFasaZBgzFAvlaYAwQAPpio3JWA9d76w7BmndqgTov/euGNqpgLyKSWbGC/WUN\npE2Ksu4r3iMmBXsRkczSFbQiIj6gYC8i4gMK9iIiPqBgLyLiAwr2IiI+oGAvIuIDCvYiIj6gYC8i\n4gMK9iIiPqBgLyLiAwr2IiI+oGAvIuIDCvYiIj6gYC8i4gMK9iIiPqBgLyLiA1kJ9k2b2s3G9+/P\nxreLiPhPVoJ9UZHV7mtrs/HtIiL+k5VgD2rKERHJJAV7EREfULAXEfEBBXsRER9QsBcR8YFYwX4S\nsAGoCEsrB6qAhd7jvLBlo4EVwDLgnINtuLgY9u6NM7ciIpKQWMF+MjA8Is0BDwH9vMebXvrxwCXe\n83DgsYNtv1Ur2LUrgRyLiEjcYgX72cC2BtKLGkgbAUwDaoE1QCVwarQNH3447NjRuEyKiEhyEm2z\nvwlYDDwNHOGllWDNOyFVQOdoG1CwFxHJnGYJfOZx4F7v9VhgInBtlHVdQ4nl5eWsWgVTpkC7dgEC\ngUAC2RARKVzBYJBgMJiy7TXUHBOpDJgO9Imx7Nde2v3e8wxgDDAv4jPOOcfNN0OPHnDLLfFmWUTE\nf4qKiqBxMbtBiTTjdAp7/RPqR+q8BlwKFAPdgF7A/GgbUTOOiEjmxGrGmQYMBdoDa7GaegDoizXR\nrAau99ZdCrzkPe8DRhKlGQegTRvYuDGJnIuISKMlfEqQBOec44kn4KOP4PLLoWNHOO64LORERCRP\nZKMZJyVCzThPPQVvvZWtXIiI+EMio3FSIhTsa2qgujpbuRAR8YesB/vt23UlrYhIumW9Geerr1Sz\nFxFJt6wG+61bYcsW1exFRNItq80469bZa9XsRUTSK2s1+9at61+rZi8ikl5ZvXlJy5bQtKlq9iIi\n6Za1Zhywq2g1r72ISPplrWYP1m7fvbtq9iIi6ZYTwV41exGR9MqJYK+avYhIemU12F99NZxzjmr2\nIiLplrVZL0P27rUa/t69WciJiEieyNtZL0OKi2HfPnuIiEh6ZD3YFxXBoYeqKUdEJJ2yHuxBY+1F\nRNItJ4L9oYdqRI6ISDrlRLBv1UrBXkQknXIi2KvNXkQkvXIi2KtmLyKSXjkR7FWzFxFJr5wI9qrZ\ni4ikV6xgPwnYAFQ0sOx2oA440ntfBuwGFnqPxxqbCdXsRUTSK9Z89pOBR4ApEemlwNnA5xHplUC/\neDOhmr2ISHrFqtnPBrY1kP4QMCpVmVDNXkQkvRJpsx8BVAEfN7CsG9aEEwQGNXaDqtmLiKRXvLcl\nbAXchTXhhIRmYVuPNe9sA/oDrwInADsjN1JeXv7t60AgQElJgDlz4syJiEgBCwaDBIPBlG2vMdNl\nlgHTgT7eYxYQanTpAqwDTgU2RnzuHawTd0FE+gFTHAPMmQM33wz//Gc8WRcR8Y9kpziOt2ZfAXQM\ne78aOBnYCrTHavX7ge5AL2BVYzZ64onw6aewfz80bRpnjkREJKZYbfbTgA+A3sBa4OqI5eFV9CHA\nYqzN/s/A9cDXjclE69bQoQOsatShQURE4pX1O1WFXHghXHUV/PSnmc+QiEiuy/s7VYWceCIsWZLa\nbUa7+9U//wkrV6b2u0REclnOBPsTToClS1O3vaefhvPOs36Au+6CZ56B0AnFv/0bjBuXuu8SEcl1\n8XbQpk3XrrB2bWq29fXXcPfdsHMnfPABPP+8BfoTT4Tu3eHDD+GTT6zm3yxn/gIiIumTM232a9bA\nkCHwxRfJbXzRIrj0UqvVv/cedOkCxxxjo3w6dICyMpg2Db78EiZMgDPOSO77REQyIdk2+5wJ9jU1\ncNhhsHt3csMvL7sMvv99+PWvbez+I4/AX/4CLVrA/ffb1Aw//jF88w1Mngx/+hMMHBh9e199BZs2\nQZ8+iedJRCRZBdNBW1wM7dpZcE1UXR3MmgVXXGHvTz8diopg6FA7a5g3D9avh1/+Em69FUaPttE/\nv/lN9G3+8Y8wbJgFfBGRfJUzwR6syaWqKr7POAd33gnLl8PChdZUU1pqy846y4J627ZWox83Dl5+\nGQ45xA4Cl18OFRXw2GOwbp1tI9Ls2dCrF9x+e/LlExHJlpxpxgGrZV9+OfzsZ43f2Pvvw/nnw9FH\nQ//+9vyHP8SXoRtvhI8+grlzrbY/bpwdDHbvhvbt4bPPbLTQ55/DEUfEt20RkVQomGYcSKxm/8AD\nMH58fQ3+uuvi/96RI61m//778NZbFvzr6mD+fBvB06ULnH02/PnP8W9bRCQX5FTNfsIEa7OfOLFx\nG9qxA0pKYONGmyY5uUxZbX77djtTGDgQtm61mv2ECfDXv9pB5YMPkvseEZFEFFzNfu5ceO65xq0/\nb5413SQb6MECPUCbNjB9Orz+ug0HvfNOS//hD2HDBqv9L1sGV15pI3327En+u0VE0i2ngv1JJ8G2\nbfCrXzVu/Tlz4Ac/SH0+jjzSOntnzbKaPdjFV3feCb/4BQwaBD17whtvwPDhdjYgIpLLcqoZB6yt\nvGVLu3NV8+YH39Dw4XDDDTBiRIpzGMXevdbEdM011hFcV2d9BHv2wNSpmcmDiPhTwVxUFa6kxCYr\n69w5+jp1dVYDX74cjjoqxTmMQ3W1deI+9phdtSsikg4F1WYf0rGjtY8fzKOP2lWt2Qz0YOP3n3nG\npmf+7LPs5kXSY8cO2Lw527kQSU5eBvuqKrj3XpvuIBcMHWojdYYNs7l5JDH791un95Ahdr3FkiXW\nSf7AA/Af/2HLM2nJEps7qaQEevSAH/0I/v53a7abPh3+67+saQ9sYr0JE2Ds2OSuAs+GurrEPrdm\njfVt7fzOXaYlF+XknI+hYL9zp93FCmz++RYtbMTOnDkweLB1kuaKq66yUUHDhsG558KUKdHn+Kmr\nsyt55861g9Zhh2U0qzlh504r/8qVFjzPOss65uvqbF6jZcvgzDNt3Ysvtumvp06F3/0OjjvOptYY\nP94O/D172nDZnj3td1NSAk3CqjE1NTbz6RFH2NXVc+dCt25wwQX166xcaZ3wXbvaHdOeftrmTRo3\nzq69qKuzysXo0fDxx3D88ZaHG26w3+s331jf0b590LevXZMxeHDq/241NVBba2eUkb75BsaMsQED\nt9wCxx5ro8oqKixfGzfa3+6VV2DxYst306Y22+w558CAAVZx2bPHrjupq7MydO9u29+92w5yjz4K\nW7bY2c7RR9vfq6TEvu/MMy1t2zY7QO/cafuifXs7gF95pZ0pVVfb0OYBA6yPTtIvJ9vsR42yf6Tf\n/97+gR591C58atvWanljxtgPcezYDOU4Drt2Wdv9iBFw220WpHr3rp9KedMm69Rdt85m4Fy92g4Q\ngwZBv372T1OUxr1SXd1woMi0H/3IasB9+lhAmDHD7jswalT9QXLOHAsE/frVB9tHHrGZUVu1ssA1\nYoQdGF56yQJhixYWXCZOtMD+3HPw6qu2jW3b7DFsGMycaWcQ27bBihUWvJo0sestamutE/7KK+3K\n6Ui1tbY/i4psrqWqKgtaoQPMzJk2P9N110GnTrZeba2dmdx4o302nqm19++HJ5+EZ5+1M8cmTWwK\njxNOsL9baan9xt57z357xx5r5e/UCQ4/3MoZmntqwQIL6BdeaH/Hmho7wL39tvWTzZ5tFaxQf9mM\nGbZ8yxZ79O9v94MoLbXRc82a2W/q889t23Pm2G+8rs7Keswx9nr9erjnHptevG1bm7KkdWu79/T3\nvmcXRZ5/vvXDge2zykqbxwrsAOJ3BdlBO3Gi1b7efdd28oMP2oyVLVrAm2/CRRfBz38Ol1ySoRzH\naeVKuyhr1CibZO3JJy1wTJ1qB62rrrJT/mbNLG3NGnjnHftHvu22g0/MFrr4KxEVFRb0Bgywf7DZ\ns+G00+Dkk22m0LPPTmy78dizx6aYvu8+y0+LFvFvY98+m96if//6EVvV1VZj7NTJynXPPVbDvegi\nq1F263bgNr780jrVTzrJaq59+9pBZtUqu9aiXbvkyvnxx1a7//prC9ZNmlgNOhi0CsEVV9hZyL59\nNqpsxAgLsJWVVsk57jh4+GELijt2WAC/+26rFDRvbkF10SIL2ps2WdAcONAOAmBNUAsWWLBMpvKw\nZYtVSDp0sDOBVNfCq6vtf33sWDvbqKmxg8Pxx9sB7PXXbb8MG2Z/txYt7Czh9NPtb9a2bWrzk8sK\nMtg//zzcdJMFoh49LDA99JD9ENavt/evvNJwrStXvP++1V5/9jOrOQaDFpzGj48eVNeutWD8hz9Y\nEDrpJEtfsMDOBNq1g2uvtYNBaGbPeFxzjdXSBg+2A9LJJ1uNbv16Oz3v2xd++1t7TsTWrVbjbtnS\ntrlmDZx6qp3+P/GE1RLnzbN/5IkTLXD5iXO2jw85BJ56ygJzq1ZWY58/34JZkya2bz//3Cbfa9PG\nzsRKS9N7xpcrduyw4B5+9rl9u1WKOne2A+fGjfC3v1kFqaTEfmNnn538gS3XJRvss8HFMnOmc+Dc\nqFHOPfywc9de61xxsXOHH+7cF18417Klc3v3xtxM1tXWWj47dHDuqaec69rVuf37D/6ZF15w7vzz\nnTv6aOfuvtu5f/zDuXbtnBs+3LmyMudGjnSuWzfnamoOvp2KCufuuMO5iROd27nTudWrnWvb1rlN\nmxpev7rauQkTnCstte3/4AfOffhh48u6Z49zxx7rXM+ezp17rnNHHOHciSc6d+ihzrVp49zllzv3\n2mvOff1147fpNxUVzv3tb9nORf7Yt8/+Zk8/7dwpp9j/zcKFztXVZTtn6QEcvJacg2IWavFiC/Yv\nvGCBv3Nn57p3d27wYOcefNCCSD6ZPNm5jh2du/fexn9myxbnjjnGgueTTx647KyznHv22eifXb7c\nuaOOsoPFxRc7V1JiB4qHH479vfv3O7d0qXNTpjh35JEW9OfPj53XkSOdGzHCuRkzbL9VV9uy7dud\n+/LL2N8rkoy9e50bP94qK126OHfFFc699VZhBX6SDPY52YyzYYO11S9dap04paV2mta7t53y3nij\ntfnmk7o6O8WM5zRz1iy44w5r+igurk9/+WXrB5g5097v3Wudl+3aWdPPBRdYJ9r119vyDz6wdtd4\nm36qq62j7Fe/suaojh2tCe2CC6wp4o03bNTRrFk2AumPf7Q2c5Fscc76Pd5+20YD7dljzZZdu9pI\noUGDDhyplU/S3WY/CfghsBGIvDHf7cAEoD2w1UsbDVwD7AduBmY2sM2YwX7/fgso06db+13r1jYn\nzcUX13e+5esOi1dDHbK7dllbZejq4f/8T3j8cUubNcv6BW67LXV5+OQT+J//saF2H31k+6BNG2tD\nvfZaOxBonn/JNc7Z6J/Fi62P6vXX7Tfco4fFl337rF+pbVs7IPTta/2Ce/dax/qKFdYvWFJilc/I\nodShW6iGV8TSKd3BfjDwDTCFA4N9KfAn4HvAyViwPx54ATgF6AzMAnoDkZdsxAz2kQYMsJE3d9wR\n18cK2hVX2CiM006z4YUvvmijMWprY88plKwtW6yjMdGOXJFscM6C+4oVNmX54Yfb8OfNm62StGKF\njfZp0cIqmD16WEUndM3PccfZUOEmTWyY6xdf2Hb79LHhrEOH2hluIiPMGiMTo3HKgOkcGOz/DIwF\n/kp9sB+NBfbx3jozgHJgbsT24g7248fbH/JgNwb3m8pK+4EuX27NLU88ke0ciRSuHTss8IcuUBs8\n2Gr9tbV2fUIwaFdXV1XZqKD27a3Wv3mzDYPt39+GAIdumbprlzXtHnaYHYRC14gcTDaC/QggANwK\nrKY+2D+CBfbQ/I9PAW8C/x2xvbiDvYhIPnjnHWvy3LrVmqNbt7aJEufOtbPvmhq7ZmLVKnvdubNV\n1jZtsvX69rX+tX/5F7v2pXVre7RqBQMGJBfs450uoRVwFxA+UvxgX95gVC8vL//2dSAQIBAIxJkN\nEZHcc8YZ9oh05ZV2Ed/OnXZ2UFZmAyrWrLE2/86d7XqaRYvsQsB334XHHw9SVRX8doqMZMVbs++D\ntcXv8pZ1AdYBpwFXe2n3e88zgDHAvIjtqWYvIhKnTE9xXAF0BLp5jyqgP7ABeA24FCj2lvUC5iea\nMRERSZ1YwX4a8AE2qmYt9bX3kPAq+lLgJe/5TWAkeXjFl4hIIcrJi6pERORABXmnKhERSS0FexER\nH1CwFxHxAQV7EREfULAXEfEBBXsRER9QsBcR8QEFexERH1CwFxHxAQV7EREfULAXEfEBBXsRER9Q\nsBcR8QEFexERH1CwFxHxAQV7EREfULAXEfEBBXsRER9QsBcR8QEFexERH1CwFxHxAQV7EREfiBXs\nJwEbgIqwtLHAYmAR8DZQ6qWXAbuBhd7jsVRmVEREEhcr2E8GhkekPQB8H+gLvAqMCVtWCfTzHiNT\nlMe8EgwGs52FtFL58lshl6+Qy5YKsYL9bGBbRNrOsNeHAZtTmqM8V+g/OJUvvxVy+Qq5bKnQLMHP\njQN+CewCBoald8OacLYDvwHeTyp3IiKSEol20N4NdAWeAX7vpa3H2u/7AbcBLwCtk8yfiIhkSBkH\ndtCG6wosibLsHaB/A+mVgNNDDz300COuRyVJSKQZpxewwns9Amu2AWiPte/vB7p7661q4PM9E/hO\nERFJo2lY80wNsBa4BngZq+kvAv4bOMpb96dYLX8h8BHww0xnVkREREREMmQ4sAxrBrozy3lJlTXA\nx9gZzXwv7Ujg78ByYCZwRFZyFr+GLqI7WFlGY/tyGXBOhvKYjIbKVw5UUX8x4Hlhy/KtfKVYX9kn\n2Fn2zV56oezDaOUrJ//3YUtgHtZishS4z0vPy33XFOtgKAOaY4U6LpsZSpHV2A4J9wAwynt9J3B/\nRnOUuMHYaKrwYBitLMdj+7A5tk8ryf3pNxoq3xhs9FikfCzf0djFjmDXwHyG/Y8Vyj6MVr5C2Yet\nvOdmwFxgECncd5ks+KlYhtYAtcCLWAdvISiKeH8h8Kz3+lngx5nNTsIauoguWllGYH06tdg+rcT2\ncS5rqHzw3f0H+Vm+r7AAAPAN8CnQmcLZh9HKB4WxD3d5z8VY5XgbKdx3mQz2nbFO3pAq6ndUPnPA\nLOBD4DovrSPWXID33DEL+UqVaGUpwfZhSD7vz5uw+Z6epv40Od/LV4adxcyjMPdhGVa+ud77QtiH\nTbCD2Qbqm6tStu8yGexdBr8rk07HfnTnAf8PayoIFxojWwhilSUfy/k4duV3X+BLYOJB1s2X8h2G\njZS7hQOnN4HC2IeHYaMCb8Fq+IWyD+uwMnQBhgBnRCxPat9lMtivo36GTLzXVVHWzSdfes+bgL9g\np1IbsPZFgE7AxizkK1WilSVyf3bx0vLNRur/iZ6i/lQ4X8vXHAv0z2ETFUJh7cNQ+Z6nvnyFtg+3\nA28AJ5On+64ZsBI7/SqmMDpoW1E/JcShwD+wXvEHqB9t9Gvyp4MWvnvFdLSyhDqIirFa1UoabjfN\nNWUcWL5OYa9vxab5gPwsXxEwhfopTEIKZR9GK18h7MP21Dc/HQK8B5xFHu+787Ae9Eps2FC+64b9\nwRdhQ8FCZToSa8fPt6GXkRfRXc3By3IXti+XAedmNKeJaegiwSnY0NnFWE0xvH8l38o3CGsKWET9\nMMThFM4+bKh851EY+7APsAAr28fAHV56oew7ERERERERERERERERERERERERERERERERERERiVcb\n4IYY6/wjxvJvEvjeW73PhS7XPx27JecvEtiWiIjEUMaB892EFNH4+UIiZ4dsjKHYXCW3hqX9PIHt\niIhII7yI3fhhIfASNv/Ss9h8RaGZAUM191exew8sof7+A1Af7A/FZhxchB1ADha8L8Zq9aHbUbbG\n5qIREZE0OIb6mv0xwH6+e9eeUDBv6z0f4n2mbcTyi4Anwz53+EG+N3Qg+DtwLFbT7xBPxkUyJZfv\nxyjSWEURrz+nvrYd6Ras1j4Hq/X3ilj+MXA21jwzCNjRiO+fClyBzba4qdG5FskgBXspRNVR0odi\nc4QPxO4ItBBoGbHOCupvSv7/gd9G2dbR2HTJYDfT+Ck5Np+4SLhm2c6ASArspP4mMgfTBruJ8x6s\n2WVgA+t08taZit0x6Noo2zoFeDvs+5egJhzJYQr2Ugi2YEMrK4BPafhenA6YAfxfYCnWiTsnYjnA\nSdjdgeqwm5w0NKTzTKAcaIHdCxXsNnlbkiiDiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIiIpMr/\nAlOuIOMuGhg0AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f73b7f509d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "n = 6 # a number of participants\n",
    "K = 1000 # max time\n",
    "\n",
    "\n",
    "N = 300 # trial for averaging\n",
    "\n",
    "p = np.random.rand(n) # random probability P_{X}: prob X\n",
    "s = np.sum(p)\n",
    "pp = p/s\n",
    "P = []\n",
    "for i in range(n):\n",
    "    P.append(np.roll(pp, i))\n",
    "P = np.reshape(P, (n,n))\n",
    "\n",
    "i = 0 # initial value\n",
    "\n",
    "f = cal_f(n, P, K)\n",
    "print f[i]\n",
    "print reduce(operator.add, f[i])\n",
    "\n",
    "_N = 0\n",
    "_ave = [[] for _i in range(n)]\n",
    "c = collections.Counter()\n",
    "while _N < N:\n",
    "    c = c + collections.Counter(circle_walk(P, n, i, K))\n",
    "    _N += 1\n",
    "    for k, v in dict(c).items():\n",
    "        _ave[k].append(float(v)/_N)\n",
    "last = [_ave[a][-1] for a in range(n)]\n",
    "print last\n",
    "print reduce(operator.add, last)\n",
    "        \n",
    "def plot_for_j(j):\n",
    "    plt.plot(range(1, N+1), _ave[j])\n",
    "    plt.plot([1, N], [f[i][j], f[i][j]], '--')\n",
    "    plt.xlabel(r'trials $N$')\n",
    "    plt.ylabel(r'')\n",
    "    plt.show()\n",
    "\n",
    "for j in range(n):\n",
    "    plot_for_j(j)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end",
     "slide_type": "subslide"
    },
    "slide_helper": "slide_end",
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "どう見ても別の値で収束している???"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_type": "subslide"
    },
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "ともかく(?)、このようにマルコフ連鎖としてこの状況をとらえたとき、距離の効果を入れるなら、遷移確率にその効果を埋め込むことが出来る。このとき確率$p(i,j)$と$p(j,i)$は一般には異なっていても良い。\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "internals": {
     "slide_helper": "subslide_end"
    },
    "slide_helper": "slide_end",
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "source": [
    "参考になる?[計算機実習 問題12.16 拡散が支配する1次元化学反応系のモンテカルロ・シミュレーション](http://nbviewer.ipython.org/github/ssh0/12-16_diffusion-dominant_chemical_reaction/blob/master/12-16_diffusion-dominant_chemical_reaction.ipynb)"
   ]
  }
 ],
 "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}