chap02.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from abc import ABC, abstractmethod\n",
    "import random\n",
    "import math\n",
    "\n",
    "\n",
    "# 離散確率分布の抽象クラス（現在は1次元のみ）\n",
    "class DiscreteProb(ABC):\n",
    "    @abstractmethod\n",
    "    def pmf(self, x: int):\n",
    "        # 確率量子関数\n",
    "        pass\n",
    "\n",
    "    @abstractmethod\n",
    "    def __call__(self):\n",
    "        # 確率変数\n",
    "        pass\n",
    "\n",
    "\n",
    "class Bern(DiscreteProb):\n",
    "    def __init__(self, mu: float):\n",
    "        self.mu = mu\n",
    "\n",
    "    def pmf(self, x: int):\n",
    "        if x == 0:\n",
    "            return 1 - self.mu\n",
    "        elif x == 1:\n",
    "            return self.mu\n",
    "        else:\n",
    "            return 0\n",
    "\n",
    "    def __call__(self):\n",
    "        return 1 if random.random() < self.mu else 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# サンプリングによる期待値計算とエントロピー計算\n",
    "\n",
    "\n",
    "class SimpleSampler(object):\n",
    "    def __init__(self):\n",
    "        pass\n",
    "\n",
    "    def __call__(self, p, L: int, f=None) -> float:\n",
    "        # 単純なサンプリングを用いて f(x) の期待値を計算する\n",
    "        # # p: 確率分布\n",
    "        # # f: 関数\n",
    "        # L: サンプル数\n",
    "        if f is None:\n",
    "            f = lambda x: x\n",
    "        return sum([f(p()) for _ in range(L)]) / L\n",
    "\n",
    "\n",
    "class Entropy(object):\n",
    "    def __init__(self, sampler=SimpleSampler()):\n",
    "        self.sampler = sampler\n",
    "\n",
    "    def __call__(self, p, n_samples: int = 10000) -> float:\n",
    "        return -1 * self.sampler(p, n_samples, lambda x: math.log(p.pmf(x)))\n",
    "\n",
    "\n",
    "class KLDiv(object):\n",
    "    def __init__(self, sampler=SimpleSampler()):\n",
    "        self.sampler = sampler\n",
    "\n",
    "    def __call__(self, *, p, q, L: int = 10000) -> float:\n",
    "        # p(X) と q(X) のKL距離\n",
    "        #    KL[q(x) || p(x)]]\n",
    "        # を計算する\n",
    "        return -1 * self.sampler(\n",
    "            q, f=lambda x: q.pmf(x) * log(p.pmf(x) / q.pmf(x)), L=n_samples\n",
    "        )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "近似値：0.63642175\n",
      "理論値：0.63651417\n"
     ]
    }
   ],
   "source": [
    "# ================================================\n",
    "# 2.1\n",
    "# ================================================\n",
    "\n",
    "\n",
    "# 2.1.5 サンプリングによる期待値の近似計算\n",
    "# 例題：単純な確率分布のエントロピー計算\n",
    "\n",
    "from math import log\n",
    "\n",
    "\n",
    "class MyDist(DiscreteProb):\n",
    "    def pmf(self, x: int) -> float:\n",
    "        assert x in [0, 1]\n",
    "        if x == 0:\n",
    "            return 2.0 / 3\n",
    "        else:\n",
    "            return 1.0 / 3\n",
    "\n",
    "    def __call__(self) -> int:\n",
    "        if random.random() < 2.0 / 3:\n",
    "            return 0\n",
    "        else:\n",
    "            return 1\n",
    "\n",
    "\n",
    "p = MyDist()\n",
    "H = Entropy()\n",
    "print(f\"近似値：{H(p):.8f}\")\n",
    "print(f\"理論値：{-(1/3 * log(1/3) + 2/3 * log(2/3)):.8f}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0]\n",
      "[1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]\n",
      "E[Bern(x|0.5)] = 0.4968\n",
      "E[Bern(x|0.9)] = 0.9048\n",
      "H[Bern(x|0.5)] = 0.6931471805600546\n",
      "H[Bern(x|0.9)] = 0.32662103059565345\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'H[Bern(x|μ)]')"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkIAAAG2CAYAAACTTOmSAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABavElEQVR4nO3deVxU5eIG8GdmYGbYd4ZFBEFxFxQU0VxKzLLNVtrUvK2WbXRvRd20rJtl1s9uebNMq1vdtMWyxcwkdzEUwRWVRQQVhk0Y9oGZ8/tjgCJBmWGGM8vz/XzmU4xnZh4YdR7P+573lQiCIICIiIjIAUnFDkBEREQkFhYhIiIiclgsQkREROSwWISIiIjIYbEIERERkcNiESIiIiKHxSJEREREDotFiIiIiBwWixARERE5LBYhIiIiclhWUYRWrFiBiIgIKJVKJCQkICMjo9tjp06dColEcsHtmmuu6cPEREREZA9EL0Lr1q1DSkoKFi1ahAMHDiAmJgYzZsxAWVlZl8evX78eJSUlHbcjR45AJpPh1ltv7ePkREREZOskYm+6mpCQgLFjx+Ldd98FAOj1eoSFheHRRx/Fs88+e8nHL1++HAsXLkRJSQnc3NwsHZeIiIjsiJOYL67VapGZmYnU1NSO+6RSKZKSkpCent6j51i9ejVuv/32i5ag5uZmNDc3d3yt1+tRVVUFPz8/SCQS078BIiIi6jOCIKC2thYhISGQSs0zqCVqEaqoqIBOp4NKpep0v0qlwvHjxy/5+IyMDBw5cgSrV6++6HFLlizBSy+91KusREREZB2Ki4vRr18/szyXqEWot1avXo2RI0di3LhxFz0uNTUVKSkpHV/X1NSgf//+KC4uhqenp6VjEhERkRloNBqEhYXBw8PDbM8pahHy9/eHTCaDWq3udL9arUZQUNBFH1tfX4+1a9di8eLFl3wdhUIBhUJxwf2enp4sQkQ2pqlFB01jC5pa9Ghq1aGpRYfmVj20rXoIAiDAMO1REACpRAJXhQxucie4dfzXCXIn0a8TIaJeMOe0FlGLkFwuR1xcHNLS0jBr1iwAhvk7aWlpWLBgwUUf+9VXX6G5uRl33313HyQlIktratEhr6wOZbVNKK1phlrThLLaJpRpmlFZr0VlfTOq6rSo1+p6/VoeSif4ucnh566An5sc/h4KhHq7oJ+PC0K9XRDq44JADyVkUs4hJLJ3og+NpaSkYO7cuYiPj8e4ceOwfPly1NfXY968eQCAOXPmIDQ0FEuWLOn0uNWrV2PWrFnw8/MTIzYR9ZKmqQWZheeRUViFjFNVOHSmGi26nl3EKpUASmeZ4eYkhcJZBrlMij//I1EikUCn16NBq0ODVof65lY0t+oBALVNrahtakVhZUO3ryGXSRHu54rIADdEBrgj0t8Ng1QeGBrsAYWTrFffOxFZD9GLUHJyMsrLy7Fw4UKUlpYiNjYWmzZt6phAXVRUdMHM8BMnTmDXrl3YvHmzGJGJyAT1za3YV1iF9PxK7MmvxNFzNdD/pff4uskR4q2EykMJlZfhv4GeirazN3L4uing5y6Hh8LJpFPjLTo96ppaDWeY6ppRVa9FRb0W5ZomnKluxNnzjThb3YiSmiZodXrkltUht6wOwB/D93KZFMNDPTE6zAdjwr0RF+6DYC+XXv50iEgsoq8jJAaNRgMvLy/U1NRwjhCRBdQ0tiCvrBa5akOROFhcjeziarT+pfmE+7liXIQvxg3wRcIAP4T5uljFkhatOj1KapqQX16HgvJ6nKqoR0FFHY6X1KKyXnvB8UOCPHDlMBWShqkwMtTLKr4HIntkic9vFiEWIaJe0ekF5JRosLegEnsLDENcZbXNXR7bz8cFE6L8MCHKH+Mj/RDkpezjtL0jCAKKqhqQVVSNA0XncaDoPI6d03Q6s6XyVGDaUBWmRAcgMcoPnkpn8QIT2RkWITNhESLqnaYWHX44eA6/HFUj41QlNE2tFxwT7KXEIJUHBgW6Y3CQBxIj/RDm6ypCWss6X6/F1hNl2JKjxvYT5Z0mc8ukEowO88akQQGYHO2PmH7ekHICNpHJWITMhEWIyDRnzjfgs71FWLevCOcbWjrud1c4YWyED8ZH+iE+wgfRKg94OOCZkKYWHfYWVOK342XYlVuBgor6Tr8e4qXEtTEhuG5UCEaEenIIjchILEJmwiJE1HONWh1251Xgq8xi/HpM3TEMFOrtgtvHhmFydACGh3jCSca1ef6quKoBu/IqsDO3HDtOVqCu+Y8zZxF+rrguJgQ3jg5FZIC7iCmJbAeLkJmwCBFdnFrThLScMqTlqLE7vwJNLfqOX5sQ5Ye5EyKQNFTFdXaM0NSiw7YT5fjh0Dmk5ag7/Uzjw31wa3w/zBwZ7JBn0oh6ikXITFiEiLqWVXQeb24+iV15FZ3uD/FS4srhQbgzoT+iVeZb2t5R1Te3YkuOGt9lncX2k+UdZ9lcnGW4emQQ/jZxAEaEeokbksgKsQiZCYsQUWc5JRq8ufkEtuSUAQAkEiCmnzeShgbiiiEqDA324HwWC1FrmvBt1ll8tb8Y+eV/zCmaEh2A+VOjkDDAlz97ojYsQmbCIkRkUFBeh//bkosfD51r25sLuHlMPzw2bZBdXuFlzQRBQFZxNT7ZU4gfDp7rOEs0pr835k8diGlDAnnFGTk8FiEzYREiR1dQXod3f8vDd9lnOz5wrx0VjCenRyOKE3dFV1TZgPd35OOrzDPQtm0LEhnghnsvG4CbRveDi5xbfJBjYhEyExYhclSnKurxzm+5+C7rjwKUNDQQT06PxvAQzkmxNmW1TVizqxCf7z2N2rYrznxcnTF7fDhmJ0YgwEMhckKivsUiZCYsQuRoquq1WLrpOL7cX9xRgKYNCcQTSdEY2Y8FyNrVNbfiy33FWLP7FM6cbwRg2PPslvh+eGhyFPr7cRiTHAOLkJmwCJGj0OsFrN1XjKW/HEd12wKI04YE4vGkQRjVz1vccGS0Vp0em4+psWpnAbKKqgEY5nVdFxOC+VOjMCSIf5+RfWMRMhMWIXIER87W4J/fHUF2cTUAw8agr8wagfgIX3GDkVlknKrCiq152H6yvOO+pKGBeHwaz/KR/WIRMhMWIbJnNY0teGvzCXy69zT0gmH7iyenR2NuYjhXf7ZDR87W4L1t+dh4pATtf5tfNTwIKVdGc80nsjssQmbCIkT2SK8X8PWBM3j95+OorNcCMAyZ/POaoVB52tYu72S8/D9dCSgIhrWgbogJwRNJ0YjwdxM7HpFZsAiZCYsQ2ZsjZ2vwwoYjHfNGBga646Xrh2PiQH9xg1GfO6muxf/9ehI/HykFADhJJXhyejTmT4niOkRk81iEzIRFiGxZU4sOeWV1OKmuxQl1LY6X1GJHbjkEAXCTy/B40iDcM2EA5E4cBnNkh8/UYNnmEx1ziC4fHIC3bouFj5tc5GREpmMRMhMWIbJFdc2t+PuXB/HLsVJ09af2+pgQPDdzKIK8OAxGBoIg4Kv9Z/DChiNobtUj1NsF7945GqP7+4gdjcgkLEJmwiJEtqZM04R7PtqHYyUaAIZF9aJVHoZbkAfG9PfmgojUrWPnNHj480wUVjbAWSbBczOH4p4JEdzDjGwOi5CZsAiRLckrq8XcNftwtroR/u5yfDAnHqPDvPkhRkapbWrBM98cwsbDhrlDSUMD8drNo+DvztWpyXZY4vObkwiIrFjGqSrc/F46zlY3YoC/G9bPn4gx/X1YgshoHkpnrLhzDF68bhjkMim25JThquU7sPV4mdjRiETFIkRkpX48dA53r/4dNY0tGN3fG9/Mn8CtFKhXJBIJ7pk4ABsWTMRglQcq6rSY9/E+vPDdETRqdWLHIxIFixCRlWnV6bHk5xws+F8WtK16XDlMhf/dNx6+vNqHzGRosCc2LJiIeRMjAACf7j2Na9/Zia3Hy+CAsyXIwXGOEOcIkRWprGvGo19kYU9+JQDgvssGIHXmUMi4/gtZyI6T5Xjqq4Mor20GAAwP8cSjVwzElcOCuO4QWR1OljYTFiGyRtnF1Zj/WSZKaprgKpdh6S2jcO2oELFjkQM4X6/Fe9vz8dne02hoGyIbFOiORy4fiOtjQliIyGqwCJkJixBZE0EQ8PnvRVj8wzFodXpE+rth5ew47hNFfa6qXouPdp/Cx3sKUdvUCgCYNMgfy5Nj4cery8gKsAiZCYsQWYsyTROe+eYQtp4wrP575TAVlt0WA0+ls8jJyJFpmlrw3z2FeHdrHppa9AjyVGLFXaMRF+4rdjRycCxCZsIiRNbgp0MleP67w6huaIHcSYqnZwzGvZcN4KXxZDVOlNZi/ueZKCivh5NUgmevHsLfoyQqFiEzYREiMdU0tGDR90fwXfY5AIbJqf+XHMuhMLJKdc2tSF1/GD8cNPx+vXKYCm/eFgMPnrUkEXBBRSIbV1Beh2vf3Ynvss9BJpXgsSsG4tuHJ7IEkdVyVzjh37fH4uUbhkMuk2LzMTXu+2Q/mlq47hDZBxYhoj6SXVyNW1amo7iqEf19XfH1Q4lIuXIwd4knqyeRSDA7MQLrHhwPD4UTfj9Vhce+yEKrTi92NKJe49/ARH1g24ky3PHBXlTVazEy1AvrH57AHcDJ5ozu74MP5sRD7mQ4M/T8t0e4ACPZPBYhIgtbf+AM7vtkPxpbdJg0yB9fPDCeG12SzUqM8sM7d4yGVAKs21+MN345IXYkol5hESKyoFU7CpDy5UG06gXcEBuC1XPHwl3hJHYsol6ZMTwIS24aCQD4z7Z8fLizQORERKZjESKykH+n5eJfG3MAGLbK+L/bYjkfiOxG8tj+ePqqwQCAV37KwWs/H0d9c6vIqYiMx7+VicxMEAQs++UE3vr1JADgHzMG45/XDuM2BWR35k+Jwv2TBgAAVm7PR9Jb2/HjoXOcN0Q2hUWIyIwEQcCSn4/j3a15AIB/XjMUj1w+UORURJYhkUjw3Myh+HBOPMJ8XVBS04QF/8vCnat+x0l1rdjxiHqERYjITARBwEs/HMMHOwzzJV66fjjumxQpcioiy5JIJEgapsKvT07BE0mDoHCSIr2gEle/vRPrD5wROx7RJbEIEZlBi06P1PWH8fGeQkgkwKs3jsTcCRFixyLqM0pnGZ5IisaWlClIGhoInV7AM98cQnp+pdjRiC6KRYiol87XazFndQbW7iuGVAK8cUsM7kzoL3YsIlGE+brig9nxuHZUMFp0Ah76LBP55XVixyLqFosQUS/kqmsx6z+7kV5QCTe5DO/Pjsctcf3EjkUkKqlUgmW3xmB0f2/UNLbgbx/vQ1W9VuxYRF1iESIy0W/H1bjxP3twurIB/XxcsP7hiZg+TCV2LCKroHSWYVXbJOrTlQ148NP9aG7l/mRkfayiCK1YsQIRERFQKpVISEhARkbGRY+vrq7GI488guDgYCgUCkRHR2Pjxo19lJYcnSAI+GBHPu79ZD/qmluRMMAX3y+4DIODuHEq0Z/5uyvw0T1j4aF0wr7C83j660O8tJ6sjuhFaN26dUhJScGiRYtw4MABxMTEYMaMGSgrK+vyeK1Wi+nTp6OwsBBff/01Tpw4gVWrViE0NLSPk5Mj0usFLP7xGF7deByCANyZ0B+f3psAXze52NGIrNLAQA+svDsOTlIJNmSfw+zVGTh8pkbsWEQdJILI9TwhIQFjx47Fu+++CwDQ6/UICwvDo48+imefffaC41euXIk33ngDx48fh7Ozs0mvqdFo4OXlhZqaGnh6evYqPzmO5lYd/v7VIfxw8BwAwxpB9142ABIJF0okupSvM88gdf0htOgMHznXjArG368cjAH+biInI1tiic9vUc8IabVaZGZmIikpqeM+qVSKpKQkpKend/mY77//HomJiXjkkUegUqkwYsQIvPrqq9Dpuh97bm5uhkaj6XQjMkZdcyvu/Xg/fjh4Ds4yCd6+PRb3TYpkCSLqoVvi+uG3p6bixtGhkEiAnw6VIOmt7Xju28PQNLWIHY8cmKhFqKKiAjqdDipV5wmmKpUKpaWlXT6moKAAX3/9NXQ6HTZu3IgXXngBb775Jl555ZVuX2fJkiXw8vLquIWFhZn1+yD7Vl7bjNs/SMeuvAq4ymVYPXcsbojlUCyRscJ8XfF/ybHY+NgkXD44ADq9gP/9XoSnvjzIuUMkGtHnCBlLr9cjMDAQH3zwAeLi4pCcnIznn38eK1eu7PYxqampqKmp6bgVFxf3YWKyZTUNLUj+IB1Hzmrg5ybH2gfGY3J0gNixiGza0GBPfDRvHD69dxycZRL8ekyNHw6ViB2LHJSoRcjf3x8ymQxqtbrT/Wq1GkFBQV0+Jjg4GNHR0ZDJZB33DR06FKWlpdBqu16nQqFQwNPTs9ON6FJ0egGPrs1CQXk9QryU+Hr+BIzq5y12LCK7MWlQQMdefIs2HEFFXbPIicgRiVqE5HI54uLikJaW1nGfXq9HWloaEhMTu3zMxIkTkZeXB71e33HfyZMnERwcDLmcV+6Q+SzddBw7TpZD6SzFqrnxnNRJZAEPTx2IIUEeON/QgkUbjoodhxyQ6ENjKSkpWLVqFT755BPk5ORg/vz5qK+vx7x58wAAc+bMQWpqasfx8+fPR1VVFR5//HGcPHkSP/30E1599VU88sgjYn0LZIe+yzqL99s2T112awyGh3iJnIjIPsmdpFh2awxkUgl+OlyCnw9ziIz6lpPYAZKTk1FeXo6FCxeitLQUsbGx2LRpU8cE6qKiIkilf/S1sLAw/PLLL3jyyScxatQohIaG4vHHH8czzzwj1rdAdubQmWo8880hAMAjl0fh2lEhIicism8jQr3w0JRIrNiajxc2HMH4SD/4cG0u6iOiryMkBq4jRN0pq23CDe/uRklNE6YNCcSqOfGQSnmJPJGlNbfqcO2/dyG3rA6zYkOw/PbRYkciK2R36wgRWZP65lY89GkmSmqaEBXghv+7PZYliKiPKJxkeOPWGEglwHfZ57Ah+6zYkchBsAgRwVCC5n20DweKquGpdMKqOfHwVJq2cjkRmSY2zBv3T44EADy5LhtrM4pETkSOgEWIHF57CcoorIKH0gn/vTcBkQHuYscickhPzxiC28eGQS8Az64/jJXb88WORHaORYgc2l9L0Kf3JiA2zFvsWEQOSyaVYMlNI/HQlCgAwGs/H8eSn3O48jRZDIsQOSyWICLrJJFI8OzVQ5B69RAAwPvbC5C6/jB0epYhMj8WIXJINQ0tuOejDJYgIiv24JQoLL15FKQSYO2+Yvzjq4PQswyRmbEIkcM5W92IW1buwb7C8yxBRFbutrFh+M9dYyCTSrA+6ywW/3iMw2RkVixC5FCOnqvBjSt2I7esDkGeSnz5YCJLEJGVu2pEMN68NQYA8PGeQizfkityIrInLELkMHbmliP5/b0oq21GtMod6x+egKHBXFCTyBbMGh2KxTcMBwC8nZaLNbtOiZyI7AWLEDmEb7POYN5H+1DX3Irxkb746qEJCPF2ETsWERlhTmIEnpoeDQBY/OMxfLW/WOREZA9YhMjuZZ4+j6e+PIhWvYDrY0Lwyd/GwcuFiyUS2aIFVwzEfZcNAAA8880hbD1eJnIisnUsQmTXapta8MS6LOgF4LqYECxPjoXCSSZ2LCIykUQiwfPXDMWtcf2gF4B/fH0QVfVasWORDWMRIrv24vfHUFzViFBvF/zrxhHcO4zIDkgkErxy4whEq9xRUafFwg1HxI5ENoxFiOzWDwfP4ZsDZyCVAMtvj+XeYUR2ROEkw7JbYyCTSvDjoRL8dKhE7Ehko1iEyC6dq27E898eBgA8cvlAjI3wFTkREZnbqH7eeHiqYSuOFzYcQUVds8iJyBaxCJHd0ekFPLkuG5qmVsSEeeOxaYPEjkREFvLoFYMwJMgDVfVavPDdES62SEZjESK788GOAvx+qgquchneTo6Fs4y/zYnsldxJijdvi4GTVIKfj5TiBw6RkZH4CUF25adDJXjjl+MAgBevG44IfzeRExGRpQ0P8cKjVxjO/C7ccARltU0iJyJbwiJEdiMtR43H1xoulU+OD8Ot8f3EjkREfeThy6MwPMQT1Q0tmLM6A2oNyxD1DIsQ2YVduRWY//mBjkUTX71pJCQSXipP5CicZVK8fftoBHgocLy0Fje/twcF5XVixyIbwCJENm9fYRXu/+9+aFv1uHKYCm/eZriklogcy8BAd6yfPwERfq44c74Rt65Mx+EzNWLHIivHIkQ27WBxNeZ9tA+NLTpMiQ7AO3eO5uRoIgcW5uuKr+dPwIhQT1TWa3H7B+nYnVchdiyyYvzEIJt1vl6LeR//sZHqyrvjuH0GEcHfXYEv7h+PCVF+qNfqMO+jfUjLUYsdi6wUixDZrKW/nEBVvRbRKnd8OHcsXOQsQURk4KF0xkfzxuKakcHQ6vR47Iss5KprxY5FVohFiGzSweJqrN1XBAB4ZdZIuCucRE5ERNZG4STD8ttjMT7SF/VaHR74NBM1jS1ixyIrwyJENkevF7BwwxEIAnDj6FCMG8DtM4ioa84yKVbcOQah3i44VVGPx9dmQafn6tP0BxYhsjnr9hfj4JkaeCickDpziNhxiMjK+bkr8P7sOCicpNh2ohxv/XpC7EhkRViEyKacr9di6SbDytFPTI9GoIdS5EREZAtGhHrh9ZtHAQBWbM3HxsPcioMMWITIpryx+QTON7RgsMoDcxPDxY5DRDZk1uhQ3D9pAADg718dxPFSjciJyBqwCJHNOHSmGl9kGCZIL75hOJy4XhARGemZq4Zg4kA/NGh1eOC/mahu0IodiUTGTxKyCTq9gBc2HIUgADfEhiAh0k/sSERkg5xkUrx7xxiE+bqgqKoBC/6XhVadXuxYJCIWIbIJ723Lw8HiargrnPDczKFixyEiG+bjJseqOfFwlcuwK68CS34+LnYkEhGLEFm97OJq/N+WXADAS9cPh8qTE6SJqHeGBHnizVtjAACrd53CN5lnRE5EYmERIqtW39yKJ9rW/bhmVDBuGhMqdiQishNXjwzGY1cMBACkfnsY2cXV4gYiUbAIkVV7+cdjKKxsQLCXEq/OGgmJhLvKE5H5PJEUjaShKmhb9Xjw0/0o0zSJHYn6GIsQWa1NR0qxdl8xJBLgzdti4OXqLHYkIrIzUqkE/5ccg4GB7lBrmvGPrw9BELjytCNhESKrpNY04dn1hwAAD0yOxIQof5ETEZG98lA6Y+XdYyCXSbH9ZDm+P3hO7EjUh1iEyOro9QL+/tVBVDe0YHiIJ56aPljsSERk5wYGemBB23yhxT8cw/l6ri/kKFiEyOq8tz0fO3MroHSW4u3bYyF34m9TIrK8h6ZEIVrljsp6LV75KUfsONRH+AlDVmVvQSXe3GzYEHHx9SMwMNBD5ERE5CjkTlIsuWkUJBLgmwNnsCu3QuxI1AdYhMhqlNc247EvsqAXgJvGhOLW+H5iRyIiBxMX7oPZ4w37GD737WE0anUiJyJLYxEiq6DTC3hiXRbKapsxKNAdr8wawUvliUgU/5gxGEGeShRVNWB52kmx45CFWUURWrFiBSIiIqBUKpGQkICMjIxuj/34448hkUg63ZRKrjRs6975LRe78yrh4izDf+4aA1e5k9iRiMhBeSid8fKsEQCAD3eewpGzNSInIksSvQitW7cOKSkpWLRoEQ4cOICYmBjMmDEDZWVl3T7G09MTJSUlHbfTp0/3YWIyt915FXg7zbCFxr9uHIFBKs4LIiJxTR+mwsyRQdDpBSzccIRrC9kx0YvQW2+9hfvvvx/z5s3DsGHDsHLlSri6umLNmjXdPkYikSAoKKjjplKp+jAxmVNpTRMeX5sFQQCS48Nw0xjOCyIi67DouuFwlctwoKgaG7K5tpC9ErUIabVaZGZmIikpqeM+qVSKpKQkpKend/u4uro6hIeHIywsDDfccAOOHj160ddpbm6GRqPpdCPxNbfq8NBnmaio02JIkAdeumG42JGIiDqoPJV45HLD2kJLfs5BfXOryInIEkQtQhUVFdDpdBec0VGpVCgtLe3yMYMHD8aaNWuwYcMGfPbZZ9Dr9ZgwYQLOnOl+5+AlS5bAy8ur4xYWFmbW74NM8+L3x5BdXA1PpRPenx0HpbNM7EhERJ3ce9kAhPm6QK1pxnvb8sWOQxYg+tCYsRITEzFnzhzExsZiypQpWL9+PQICAvD+++93+5jU1FTU1NR03IqLi/swMXVlbUYRvsgogkQC/PuO0Qj3cxM7EhHRBZTOMjw/cxgA4IOdBSiuahA5EZmbqEXI398fMpkMarW60/1qtRpBQUE9eg5nZ2eMHj0aeXl53R6jUCjg6enZ6UbiySo6j4UbDMOZT02PxtTBgSInIiLq3ozhKkwc6Adtqx6vbuSK0/ZG1CIkl8sRFxeHtLS0jvv0ej3S0tKQmJjYo+fQ6XQ4fPgwgoODLRWTzKi8thnzPzsArU6PK4ep8PDUgWJHIiK6KIlEgoXXDodMKsHPR0qxJ58rTtsT0YfGUlJSsGrVKnzyySfIycnB/PnzUV9fj3nz5gEA5syZg9TU1I7jFy9ejM2bN6OgoAAHDhzA3XffjdOnT+O+++4T61ugHtLrBTz6xQGUapoQFeCGN2+LgVTKRROJyPoNDvLA3Qn9ARg2ZW3V6UVOROYi+qp1ycnJKC8vx8KFC1FaWorY2Fhs2rSpYwJ1UVERpNI/+tr58+dx//33o7S0FD4+PoiLi8OePXswbNgwsb4F6qEfD5dgb0EV3OQyvD87Hh5KZ7EjERH12JPTo7Hh4DkcL63F578XYe6ECLEjkRlIBAdcJUqj0cDLyws1NTWcL9RHWnR6JL21HacrG/DU9Gg8Om2Q2JGIiIz26d7TeOG7I3BxluHnxychwp8XevQlS3x+iz40Ro5h3b5inK5sgL+7An+7bIDYcYiITHLXuP6YEOWHxhYdnvwym0NkdoBFiCyuQdvasYXGY9MGwk0h+ogsEZFJpFIJ3rg1Bh5KJ2QVVXNtITvAIkQW99HuQpTXNiPM1wW3j+0vdhwiol4J9XbByzcYNmV9Oy0Xh85UixuIeoVFiCyqukGLldsN/2J6avpgyJ34W46IbN8NsSG4ZlQwWvUCnlyXjUatTuxIZCJ+KpFFvbc9H7VNrRgS5IHrY0LEjkNEZBYSiQT/mjUCgR4K5JfX4/VNx8WORCZiESKLKa1pwse7CwEAT181mGsGEZFd8XaV441bYwAAH+8pxM7ccpETkSlYhMhi3k7LRXOrHmMjfHA5t9EgIjs0JToAcxLDAQCLvj/Kq8hsEIsQWcTxUg2+3G/Y3PaZq4ZAIuHZICKyT/+YMRi+bnIUlNdj3X5u6m1rWITI7ARBwAvfHYFOL+Cq4UGIj/AVOxIRkcV4KJ3x6BWGfROXb8lFg7ZV5ERkDBYhMrtvDpzFvsLzcHGWYeF13PqEiOzfXQnh6O/rivLaZny485TYccgILEJkVjUNLViyMQcA8HjSIIR4u4iciIjI8uROUvx9xmAAwPvb81FR1yxyIuopFiEyq2WbT6CyXouBge7420RupUFEjuPakcEYGeqFeq0O77Stpk/Wj0WIzObwmRp89vtpAMDLN4zg4olE5FCkUglSrx4CAPj89yIUVtSLnIh6gp9UZBY6vYB/fncYgmBYcTUxyk/sSEREfW7CQH9MiQ5Aq17AG5tPiB2HeoBFiMxi7b4iHDxTAw+FE56fOVTsOEREonn26iGQSICfDpXgYHG12HHoEliEqNfKapuwdJPhXz4pV0Yj0FMpciIiIvEMDfbEjaNDAQCp6w+juZX7kFkzFiHqFUEQ8Nz6w6hpbMHwEE/MHh8udiQiItE9c9UQ+Lg641iJBq//zCEya8YiRL3yzYGz2JJTBmeZBG/eFgMnGX9LERGpPJVY1rYP2Zrdp/DbcbXIiag7/NQik52rbsRL3x8FADw5PRpDgjxFTkREZD2mDVXhngkRAIC/f3UIZZomcQNRl1iEyCSCIOCZbw6htrkVsWHeeGBSpNiRiIiszrNXD8HQYE9U1Wvx5JfZ0OsFsSPRX7AIkUk+/70IO3MroHCSckiMiKgbSmcZ3rljNFycZdidV4n3dxSIHYn+gp9eZLSiyga82raNxtNXDUFUgLvIiYiIrNfAQHe8dP1wAMCbm08gq+i8yInoz1iEyCg6vYC/f30QDVodEgb4Yl7b+DcREXXv1vh+uHZUMFr1Qsc/JMk6sAiRUT7YUYCMU1Vwlcvwxi0xkEolYkciIrJ6EokEL1w7DE5SCfYVnseRszViR6I2LELUY0fO1uCtXw3rYSy6bhj6+7mKnIiIyHaoPJWYOTIYAPDxnkJxw1AHFiHqkUatDo+tzUKLTsBVw4NwW3yY2JGIiGzOvIkRAIDvs8+hoq5Z3DAEgEWIeuiVn46hoLweKk8Fltw0EhIJh8SIiIw1ur8PYsK8odXp8cXvRWLHIbAIUQ9sOabG521/YN+8NRY+bnKRExER2a72i0w+3XsaLTq9uGGIRYgurqy2CU9/cwgAcN9lA3DZIH+RExER2baZI4MR4KFAWW0zNh4uETuOw2MRom4JgoCnvz6EqnothgR54B9XDRY7EhGRzZM7SXF3gmGDak6aFh+LEHVr28lybDtRDrmTFP++YzQUTjKxIxER2YU7E/rDWSZBVlE1sourxY7j0FiEqEuCIOCdtFwAwJzx4YhWeYiciIjIfgR4KHDdqBAAwMe7T4mcxrE5GXOwRqMx+gU8PbkjuS1Kz6/EgaJqyJ2keGAyN1QlIjK3eyZGYH3WWfx0uATPzRyKQE+l2JEcklFFyNvb26jLpiUSCU6ePInISH6Q2pp//2Y4G3TH2DD+4SQisoBR/bwRF+6DzNPn8dGeQjxz1RCxIzkko4oQAHz99dfw9fW95HGCIGDmzJkmhSJx7Suswt6CKjjLJHhwSpTYcYiI7Nb9kyKReToTH+wowJXDVBjd30fsSA7HqCIUHh6OyZMnw8/Pr0fHR0ZGwtnZ2aRgJJ5/t80NuiUuDCHeLiKnISKyXzOGq3BdTAh+OHgOj6/NxsbHJ8FdYfQ5CuoFoyZLnzp1qsclCACOHDmCsDBuxWBLsoursTO3AjKpBA9P5dkgIiJLkkgkeGXWCIR6u6CoqgELNxwRO5LD4VVj1En7lWI3jg5FmC83VSUisjQvF2e8fXsspBJg/YGz2JB9VuxIDsXk82+nTp3Czp07cfr0aTQ0NCAgIACjR49GYmIilEpOrrVFR87WIO14GaQS4JHLB4odh4jIYcRH+OKxaYOwfEsu/vntEYzp78N/jPYRo4vQ559/jrfffhv79++HSqVCSEgIXFxcUFVVhfz8fCiVStx111145plnEB4ebonMZCHv/pYHALguJgQD/N1ETkNE5FgWXD4Qu3IrsP/0eTy+NgtfPpgIJxkHbizNqJ/w6NGj8e9//xv33HMPTp8+jZKSEmRmZmLXrl04duwYNBoNNmzYAL1ej/j4eHz11VeWyk1mlldWi01HSwHwbBARkRicZFL8X3IsPBROOFBUjZXb88WO5BCMKkKvvfYafv/9dzz88MNdToJWKBSYOnUqVq5ciePHj/d4/aAVK1YgIiICSqUSCQkJyMjI6NHj1q5dC4lEglmzZhnzbVAXPthRAAC4cpiKq0gTEYkkzNcVi2cNBwCs3F6A8/VakRPZP6OK0IwZMy55TGtrK+rq6uDn54e4uLhLHr9u3TqkpKRg0aJFOHDgAGJiYjBjxgyUlZVd9HGFhYX4+9//jkmTJvU4P3WtTNOE77LOAQDXDSIiEtms2FAMD/FEXXMr3m/7RypZjsmDj0VFRV3eNmzYgKFDh/b4ed566y3cf//9mDdvHoYNG4aVK1fC1dUVa9as6fYxOp0Od911F1566SWuWm0GH+0phFanR3y4D+LCuZgXEZGYJBIJUqZHAwA+2VOI8tpmkRPZN5OvGouIiOh2u43Ro0f36Dm0Wi0yMzORmpracZ9UKkVSUhLS09O7fdzixYsRGBiIe++9Fzt37jQuOHVS19yKz/aeBgDuKUZEZCWuGBKImDBvHCw2zBV64dphYkeyWyYXoaysrE5f63Q65OfnY/HixXjsscd69BwVFRXQ6XRQqVSd7lepVDh+/HiXj9m1axdWr16N7OzsHmdtbm5Gc/MfjdqUzWPt1dqMItQ2tSIywA1JQ1WXfgAREVlc+1mhuWsy8Nne03hgciRU3PfRIkwuQjExMRfcN2bMGPTr1w833XQT5syZ06tgXamtrcXs2bOxatUq+Pv79/hxS5YswUsvvWT2PLauRafH6l2nAAAPTIqEVNrzDXWJiMiyJg/yR3y4D/afPo8VW/Ow+IYRYkeyS2ZfoMDX17fT2ZeL8ff3h0wmg1qt7nS/Wq1GUFDQBcfn5+ejsLAQ1113HZycnODk5IT//ve/+P777+Hk5IT8/K4vNUxNTUVNTU3Hrbi42PhvzA79cPAcSmqa4O+uwKzRoWLHISKiP5FIJEi50jBXaG1GMc5WN4qcyD6ZXIQ0Gs0Ft9zcXDz33HN4/PHHO93fHblcjri4OKSlpXXcp9frkZaWhsTExAuOHzJkCA4fPozs7OyO2/XXX4/LL78c2dnZ3e5rplAo4Onp2enm6ARB6Lhkft7ECCidZSInIiKiv5oQ5Y/xkb7Q6vR497dcsePYJZOHxry9vbucLC0IAr799lssXrwYgiBAIpFAp9N1+zwpKSmYO3cu4uPjMW7cOCxfvhz19fWYN28eAGDOnDkIDQ3FkiVLoFQqMWJE51OD3t7eAHDB/XRx20+W43hpLVzlMtydwBXAiYis1VNXDsatK9Px1f4zmD9lIPr7cesNczK5CG3dutUsAZKTk1FeXo6FCxeitLQUsbGx2LRpU8cE6qKiIkilXGLc3NrPBt0+tj+8XJ1FTkNERN0ZG+GLSYP8sTO3Au/8los3br1wji6ZTiIIgiB2iL6m0Wjg5eWFmpoahxwmO3ymBte9uwsyqQTb/zEV/Xz4rwsiImuWefo8bn5vD5xlEux65gqHvYLMEp/fPNXigN7fYZhUfu2oYJYgIiIbEBfug/hwH7ToBHy0u1DsOHbFqKGxAQMGdLuI4sU88cQTPV5biCyrqLIBGw+XAAAenMztNIiIbMUDkyOx/9NMfP77aSy4YiDcFSbPbqE/Meqn+PHHH5v0IhERESY9jszvw10F0AvA5OgADAtxvGFBIiJblTRUhcgANxSU12NtRhHum8TdAMzBqCI0ZcoUS+WgPlBZ14wv9xvWUHqI22kQEdkUqVSC+ydFInX9YazZdQpzJ0TAWcYZLr1l8k+wvLy82187fPiwqU9LFvTf9NNoatFjZKgXEqP8xI5DRERGunF0KPzdFThX04QfD50TO45dMLkIjRw5Ej/99NMF9y9btgzjxo3rVSgyv0atDv9NLwQAPDgl0qS5XkREJC6lswz3TDCs/fb+9gI44IXfZmdyEUpJScHNN9+M+fPno7GxEWfPnsW0adOwdOlS/O9//zNnRjKDrzKLcb6hBf19XXHV8Au3LyEiIttw9/hwuMplOF5ai115FWLHsXkmF6Gnn34a6enp2LlzJ0aNGoVRo0ZBoVDg0KFDuPHGG82ZkXqpVafHqp2GBRTvnzQAThxTJiKyWd6uctwWb9hSqn1xXDJdrz4RBw4ciBEjRqCwsBAajQbJycldbpZK4vr5SCmKqxrh6ybHLXFd78dGRES2497LBkAmlWBnbgWOnqsRO45NM7kI7d69G6NGjUJubi4OHTqE9957D48++iiSk5Nx/vx5c2akXhAEoWMBxbmJEXCRc3NVIiJbF+bripkjgwHwrFBvmVyErrjiCiQnJ2Pv3r0YOnQo7rvvPmRlZaGoqAgjR440Z0bqha0nynDkrAYuzjLMSeTmqkRE9uLBtmVQfjxUguKqBpHT2C6Ti9DmzZvx2muvwdn5jw07o6KisHv3bjz44INmCUe9o9cLeHPzSQDAnMRw+LjJRU5ERETmMiLUC5MG+UOnFzrmgZLxTC5C3S2uKJVK8cILL5gciMxn09FSHD2ngbvCCQ9N4XYaRET2Zn7b3+3r9hWjoq5Z5DS2yagitHbt2h4fW1xcjN27dxsdiMxDpxfw1q+Gs0F/u2wAzwYREdmhxCg/xPTzQnOrHp/sKRQ7jk0yqgi99957GDp0KJYuXYqcnJwLfr2mpgYbN27EnXfeiTFjxqCystJsQck4G7LPIq+sDl4uzrhv0gCx4xARkQVIJBLMn2o4K/TJnkLUNbeKnMj2GFWEtm/fjtdffx2//vorRowYAU9PTwwaNAgjR45Ev3794Ofnh7/97W/o378/jhw5guuvv95SuekiWnR6LN+SCwB4aEoUPJXOl3gEERHZqiuHBSEywA2aplZ88XuR2HFsjlGbrgLA9ddfj+uvvx4VFRXYtWsXTp8+jcbGRvj7+2P06NEYPXo0pFIu2Cemr/afQVFVA/zdFZg7gVeKERHZM6lUgocmR+Hpbw7hw10FmDMhHAonLpXSU0YXoXb+/v6YNWuWGaOQOTS16PDOb4azQY9cHgVXuclvMRER2YgbRofgrV9PolTThO+yziJ5bH+xI9mMXp+60Wq1OHPmDIqKijrdSBz/+70IJTVNCPZS4o5x/INAROQIFE4y3HuZYT7o+9sLoNNzM9aeMrkI5ebmYtKkSXBxcUF4eDgGDBiAAQMGICIiAgMGcHKuGBq1OvxnWx4A4LFpg6B05qlRIiJHcUdCf3gqnVBQUY/NR0vFjmMzTB43ueeee+Dk5IQff/wRwcHBkEgk5sxFJvjh4DlU1GnRz8cFt8T1EzsOERH1IXeFE+ZOiMA7v+VhxbY8XDUiiJ/NPWByEcrOzkZmZiaGDBlizjzUC5/9fhoAcPf4cDhzh3kiIoczb+IArN51CkfOavDb8TJMG6oSO5LVM/nTctiwYaioqDBnFuqFQ2eqcehMDeQyKW7l2SAiIofk6ybH7PGGq4X/nZYLQeBcoUsxuQi9/vrrePrpp7Ft2zZUVlZCo9F0ulHf+nyvYYL61SOD4OeuEDkNERGJ5f7JkVA6S3HwTA22nywXO47VM3loLCkpCQAwbdq0TvcLggCJRAKdTte7ZNRjNY0t2HDwLADDsBgRETkuf3cF7k4Ix4e7TuHttFxMiQ7gXKGLMLkIbd261Zw5qBe+PXAGTS16DFZ5ID7cR+w4REQksgemROLTvaeRVVSNXXkVmDQoQOxIVsukItTS0oLFixdj5cqVGDRokLkzkREEQcBnbUuq3zW+P1s/EREh0EOJOxP646PdhXh7Sy4uG+jPz4dumDRHyNnZGYcOHTJ3FjLB76eqkFdWB1e5DDeODhU7DhERWYmHpkRB7iTF/tPnkZ7PTdC7Y/Jk6bvvvhurV682ZxYywedtZ4NuiA2BBzdXJSKiNipPJe4YGwYAeDstV+Q01svkOUKtra1Ys2YNtmzZgri4OLi5uXX69bfeeqvX4ejiymubselICQDgrgROkiYios4emhqFLzKK8fupKuwtqMT4SD+xI1kdk4vQkSNHMGbMGADAyZMnO/0axyH7xpf7i9GiExAb5o0RoV5ixyEiIisT7OWC28b2w2d7i/D+9nwWoS7wqjEbpdML+CKjbZJ0AjdXJSKirv1t4gB8trcIO3IrUFnXzLXm/qLX+zDk5eXhl19+QWNjIwBwFcs+8uuxUpw53wgvF2dcFxMidhwiIrJSkQHuGBnqBZ1ewMYj3Iz1r0wuQpWVlZg2bRqio6Mxc+ZMlJQY5qrce++9eOqpp8wWkC4kCAJWbM0HAMxJDOcu80REdFHXt/2D+fvssyInsT4mF6Enn3wSzs7OKCoqgqura8f9ycnJ2LRpk1nCUdd25VXg8NkaKJ2luGdChNhxiIjIyl0bEwyJBNhXeB5nqxvFjmNVTC5Cmzdvxuuvv45+/Tpv8Dlo0CCcPn2618Goe/9pOxt0x7j+HOslIqJLCvZyQcIAXwDADwfPiZzGuphchOrr6zudCWpXVVUFhYIfzpZyoOg80gsq4SyT4P5JkWLHISIiG3F9jGHR3Q3ZLEJ/ZnIRmjRpEv773/92fC2RSKDX67F06VJcfvnlZglHF2o/GzQrNhQh3i4ipyEiIltx9YggOMskyCnRIFddK3Ycq2Hy5fNLly7FtGnTsH//fmi1Wjz99NM4evQoqqqqsHv3bnNmpDYnSmuxJUcNicSwSBYREVFP+bjJMSU6AFtyyvD9wXN46srBYkeyCiafERoxYgROnjyJyy67DDfccAPq6+tx0003ISsrC1FR/JC2hPe25QEwtPqoAHeR0xARka1pX27l+4PnuNxNG5PPCAGAl5cXnn/+eXNloYsormrAD4cMSxQ8PHWgyGmIiMgWTR+mgouzDKcrG3DwTA1iw7zFjiS6Xi+oCBgmTq9ZswYrVqxAbi43drOE93fkQ6cXMDk6gNtpEBGRSVzlTpg+TAUA2MA1hQCYUISKioowZcoUeHh4YPr06SgqKsKYMWNw33334dFHH0VsbCx27Nhh1HOuWLECERERUCqVSEhIQEZGRrfHrl+/HvHx8fD29oabmxtiY2Px6aefGvtt2JSKumZ8uf8MAOBhzg0iIqJeuCHWMDz246ES6PQcHjO6CP3973+HVqvFypUr4erqihkzZmDQoEEoKSmBWq3G1VdfjRdffLHHz7du3TqkpKRg0aJFOHDgAGJiYjBjxgyUlZV1ebyvry+ef/55pKen49ChQ5g3bx7mzZuHX375xdhvxWZ8l3UW2lY9RvXz6lgHgoiIyBSTBgXA29UZ5bXN2FtQKXYc0UkEI2dLBQUF4fvvv8e4ceNQVVUFf39/7N69G4mJiQCAgwcPYtq0aaioqOjR8yUkJGDs2LF49913AQB6vR5hYWF49NFH8eyzz/boOcaMGYNrrrkGL7/8co+O12g08PLyQk1NDTw9PXv0GLEIgoCrlu/ECXUtXpk1AnePDxc7EhER2bjU9YfxRUYRbonrh2W3xogdp8cs8flt9BmhsrIyhIcbPox9fX3h6uoKlUrV8etBQUE4f/58j55Lq9UiMzMTSUlJfwSSSpGUlIT09PRLPl4QBKSlpeHEiROYPHlyt8c1NzdDo9F0utmKo+c0OKGuhdxJys1ViYjILG6JMyyu+P3BcyirbRI5jbhMmiwtkUi6/H9jVVRUQKfTdSpSAKBSqVBa2v0OuTU1NXB3d4dcLsc111yDd955B9OnT+/2+CVLlsDLy6vjFhYWZnLmvvZ1pmFu0IzhQfBycRY5DRER2YO4cF/EhftA26rHR7sLxY4jKpMun1+4cGHH9hparRb/+te/4OVluJKpoaHBfOm64eHhgezsbNTV1SEtLQ0pKSmIjIzE1KlTuzw+NTUVKSkpHV9rNBqbKEPNrTp81zar/5a4fpc4moiIqOcenByJBz7NxGd7T+PhqVHwUDrmP7aNLkKTJ0/GiRMnOr6eMGECCgoKLjimJ/z9/SGTyaBWqzvdr1arERQU1O3jpFIpBg40rKUTGxuLnJwcLFmypNsipFAobHL/s99yylDd0AKVpwKXDfQXOw4REdmRpKEqDAx0R15ZHb7IKMIDkx3zqmSji9C2bdvM9uJyuRxxcXFIS0vDrFmzABgmS6elpWHBggU9fh69Xo/m5maz5bIW7cNiN43pB5nU9CFIIiKiv5JKJXhgciSe/voQVu86hbkTIqBwkokdq8+ZZUHF3khJScGqVavwySefICcnB/Pnz0d9fT3mzZsHAJgzZw5SU1M7jl+yZAl+/fVXFBQUICcnB2+++SY+/fRT3H333WJ9CxZRVtuEbSfLAXBYjIiILGNWbChUngqoNc3YkOWYu9IbVYRSUlJQX1/f4+NTU1NRVVV10WOSk5OxbNkyLFy4ELGxscjOzsamTZs6JlAXFRWhpKSk4/j6+no8/PDDGD58OCZOnIhvvvkGn332Ge677z5jvhWr913WWej0Asb09+a+YkREZBFyJynuvWwAAMMOBnoHXGDRqHWEZDIZSktLERAQ0KPjPT09kZ2djcjISJMDWoK1ryMkCAJmLN+Bk+o6vHrjSNyZ0F/sSEREZKdqm1ow4bXfUNvUig9mx+HK4d3P0RWbJT6/jZojJAgCoqOje3zJvDFnj+gPh87U4KS6DgonKa6NCRY7DhER2TEPpTNmjw/Hf7blY+X2fEwfpurV0ji2xqgi9NFHHxn9An9dI4gurX2S9FUjguDpoJczEhFR37lnYgQ+3HUKB4qqsa/wPMY50HZORhWhuXPn9nhVZmsccrIFTS06fH/QMGGNk6SJiKgvBHoocfOYfvgiowgf7ChgEboYb2/vHp0y0+l0JgVydFty1KhpbEGwlxITorh2EBER9Y17LxuALzKK8NtxNc5WNyLU20XsSH3C6CK0devWjv8XBAEzZ87Ehx9+iNDQULMGc1Rf7jcMi90Sx7WDiIio7wwMdEdipB/SCyqxNqMIT105WOxIfcLoIjRlypROX8tkMowfP97qrgyzReeqG7Ezl2sHERGROO4eH24oQvuK8di0QXCWib7coMXZ/3doQ77JPANBAMZH+iLcz03sOERE5GCuHK5CgIcC5bXN+PWY+tIPsAMsQlZCrxfwVdvVYrfGWf+GsEREZH+cZVIkxxs+gz7be1rkNH3DLEXIkdYbsJTfT1WhqKoB7gonXD3SehezIiIi+3ZHQn9IJcCe/Erkl9eJHcfijJ4jdNNNN3X6uqmpCQ899BDc3DoP5axfv753yRzMV/uLAQDXxQTDVW7020JERGQWod4uuGJIILbklOF/vxfhhWuHiR3Jooz+xPXy8ur0tb1tdiqG2qYWbDxi2E/t1ngOixERkbjuGh+OLTll+DrzDP4xYzCUzva7K73RRciU1aXp4n48VIKmFj0GBrpjdJi32HGIiMjBTR4UgH4+LjhzvhE/HDxn1/9I52RpK/Bl27DYbfH9ON+KiIhEJ5NKOjb8/vz3IpHTWBaLkMhy1bXIKqqGTCrBjaO5dhAREVmH2+LD4CyTILu4GkfO1ogdx2JYhETWfsn85YMDEeChEDkNERGRgb+7AleNCAZg32eFWIRE1KLTY/0BQxG6LZ5ng4iIyLrcMc4wN2jj4RK06PQip7EMFiERbT6qRkWdFv7uclw+JFDsOERERJ0kDPCDv7sCNY0t2J1XIXYci2AREtEn6YUAgNvH9neI/VyIiMi2yKQSXDVCBcBwVsge8dNXJDklGmScqoJMKsFd4/uLHYeIiKhLM0ca5gltPqa2y+ExFiGRfLKnEABw1fAgBHu5iBuGiIioG4bhMTmqG1qwJ79S7DhmxyIkguoGLb7LPgsAmDshQtwwREREF2EYHjPsgbnxkP0Nj7EIiWDdvmI0tegxNNgTYyN8xI5DRER0Ue3DY78cK7W74TEWoT6m0wv4dO9pAMA9E8K5kjQREVm9Pw+PpdvZ8BiLUB9Ly1HjzPlGeLs644bYULHjEBERXZJMKsGM4W3DY3Z29RiLUB9rv2Q+eWyYXe/mS0RE9uWatuGxTUfta3iMRagP5aprsTuvElIJMHt8uNhxiIiIemzcAF/4udnf8BiLUB9qPxuUNFSFfj6u4oYhIiIygpNMihkj7G94jEWoj2iaWrD+gOGS+Xt4yTwREdmg9uGxX+xoeIxFqI9sO1GOBq0Okf5uSIzyEzsOERGR0RIG+MLXTY7zDS3YW2Afw2MsQn1k2/EyAMD0YSpeMk9ERDbJSSa1u6vHWIT6gF4vYNvJcgDA1MHcZZ6IiGxXx9VjR0rRagfDYyxCfeDQ2RpU1WvhrnBCPFeSJiIiGzY+8o/hsXQ7GB5jEeoDW9uGxSYN8oezjD9yIiKyXU4yacfeYz/Zwd5j/FTuA+3DYpdzWIyIiOzAtXa0uCKLkIVV1DXj0JlqAMCUwQHihiEiIjIDe1pckUXIwnacLIcgAMNDPKHyVIodh4iIqNfsaXiMRcjCtp7gsBgREdmfa0bZx/AYi5AFter02NFx2TyHxYiIyH4kDPCDv7scNY0t2J1XIXYck7EIWVB2cTVqGlvg5eKM2DBvseMQERGZjUwq6Rges+XFFVmELGjrCcNl85OjA+DEy+aJiMjOXDMyBADwy1E1tK22OTzGT2cL2nq8fX4Qh8WIiMj+jBvgC393hWF4LN82h8dYhCyktKYJx0o0kEgMZ4SIiIjsjUwqwcyRtn31mFUUoRUrViAiIgJKpRIJCQnIyMjo9thVq1Zh0qRJ8PHxgY+PD5KSki56vFi2nzQMi43q5w1/d4XIaYiIiCxjZtviipuPltrk8JjoRWjdunVISUnBokWLcODAAcTExGDGjBkoKyvr8vht27bhjjvuwNatW5Geno6wsDBceeWVOHv2bB8nvzgOixERkSMYG+GLAA8FNE2tNnn1mOhF6K233sL999+PefPmYdiwYVi5ciVcXV2xZs2aLo///PPP8fDDDyM2NhZDhgzBhx9+CL1ej7S0tD5O3j1tqx672n4zcP0gIiKyZzKpBDPbrh770QaHx0QtQlqtFpmZmUhKSuq4TyqVIikpCenp6T16joaGBrS0tMDX17fbY5qbm6HRaDrdLGl/YRXqmlvh5ybHyFAvi74WERGR2K4aYRge23qiDHq9IHIa44hahCoqKqDT6aBSqTrdr1KpUFpa2qPneOaZZxASEtKpTP3VkiVL4OXl1XELCwvrVe5L+TVHDQC4fEggpFKJRV+LiIhIbPERPvBQOKGqXotDZ2vEjmMU0YfGeuO1117D2rVr8e2330Kp7H4fr9TUVNTU1HTciouLLZZJEARsaStCSUNVlziaiIjI9jnLpJg40B8AsO1E13N8rZWoRcjf3x8ymQxqtbrT/Wq1GkFBQRd97LJly/Daa69h8+bNGDVq1EWPVSgU8PT07HSzlJPqOhRXNULuJMWkQf4Wex0iIiJr0r6V1La2PTZthahFSC6XIy4urtNE5/aJz4mJid0+bunSpXj55ZexadMmxMfH90XUHms/GzQxyg9uCieR0xAREfWNKW1F6OCZalTVa0VO03OiD42lpKRg1apV+OSTT5CTk4P58+ejvr4e8+bNAwDMmTMHqampHce//vrreOGFF7BmzRpERESgtLQUpaWlqKurE+tb6OTXY4YiNH3Yxc9oERER2ZNgLxcMCfKAIAA7c23nrJDopyySk5NRXl6OhQsXorS0FLGxsdi0aVPHBOqioiJIpX/0tffeew9arRa33HJLp+dZtGgRXnzxxb6MfoGy2iZkF1cDAKYN5WXzRETkWKYMDsDx0lpsO1GOG2JDxY7TI6IXIQBYsGABFixY0OWvbdu2rdPXhYWFlg9kot9yDBPEYvp5QeXZ/eRtIiIiezQ1OhDvby/AjpPl0OsFm7hyWvShMXvCq8WIiMiRxUf4wF3hhMp6LY6cs43L6FmEzKRRq8POXMNq0knDWISIiMjxGC6j9wNgO1ePsQiZya68CjS36hHqbZgsRkRE5Iimtm0tZSvrCbEImcmWjqvFVJBIrH9MlIiIyBLa1xPKLq5GdYP1X0bPImQGer2AtOOcH0RERBTs5YLBKg/oBWBHrvXvRs8iZAbZZ6pRUaeFh8IJ4wZ0v/krERGRI/hjlWnrHx5jETKD9mGxKYMDIHfij5SIiBxb+yrT7ZfRWzN+aptB+2Xz03m1GBEREeLDfeEml6GiTouj5zRix7koFqFeKqpswEl1HWRSCaZGczVpIiIiuZMUE2xkN3oWoV5qnyQ9NsIHXq7OIqchIiKyDh3zhE5a93pCLEK99NtxQ9OdNoTDYkRERO2mRBuKUFbRedQ0tIicpnssQr1Q39yK3wuqAACXD+GwGBERUbt+Pq4YGOgOvQDszrfey+hZhHphV14FtDo9+vu6IirATew4REREVmVqtPVfRs8i1Atb24bFrhgSyNWkiYiI/qL9MvrtJ8shCNZ5GT2LkIkEQcDWtobLYTEiIqILjY3whdJZCrWmGSfUtWLH6RKLkImOntNArWmGi7MMCVxNmoiI6AJKZxkSIw270W+30t3oWYRM1D4sdtkgfyidZSKnISIisk7tV49tt9LL6FmETJT2p/lBRERE1LUpgw2fk/sKq1Df3CpymguxCJmgsq4ZB89UAwAuH8wiRERE1J0IP1f093VFi05Aen6l2HEuwCJkgm0nyiEIwLBgTwR5KcWOQ0REZLUkEolVD4+xCJngtxMcFiMiIuqp9iK07WSZ1V1GzyJkpBadHjvaGu0VQ1mEiIiILiUxyg9ymRTFVY04VVEvdpxOWISMtL/wPGqbWuHrJkdMP2+x4xAREVk9N4UTxg7wAWB9w2MsQkZqX0RxanQAZFKuJk1ERNQT1jpPiEXISO27zXM1aSIiop6bEm343NxbUImmFp3Iaf7AImSE4qoG5JXVQSaVYHJbsyUiIqJLi1a5I8hTiaYWPTJOVYkdpwOLkBH25FcAAEaHecPLxVnkNERERLbDWi+jZxEywt4CQ4NNjPITOQkREZHtad+NfgeLkO0RBAF7CwwrYo6PZBEiIiIyVmKkHyQSILesDmWaJrHjAGAR6rHiqkaU1DTBWSbBmP4+YschIiKyOT5ucgwL9gQApBdYx3YbLEI91H42KKafN1zk3G2eiIjIFBMH+gMAdudViJzEgEWoh/aeMhShhEhfkZMQERHZrvZ5tnusZANWFqEeEAQBv7dNlOb8ICIiItONi/CFk1SCM+cbUVTZIHYcFqGeOHO+EWerG+EklSAunPODiIiITOWmcEJsmDeAP5alEROLUA+0zw8a1c8LrnInkdMQERHZtgnt84SsYHiMRagH9nJYjIiIyGwmtM0TSs+vgCAIomZhEeqB9jNCCSxCREREvTa6vzeUzlJU1GlxUl0nahYWoUsormrA2epGyKQSxHN+EBERUa8pnGQYG2G4ClvseUIsQpfwe9vGcKP6ecFNwflBRERE5jAhyjBPSOzL6FmELuH39mGxARwWIyIiMpf2eUJ7CyrRqtOLloNF6BLaF1Icz4UUiYiIzGZEqBc8lE6obWrF0XMa0XKwCF3E2epGFFe1zQ+KYBEiIiIyF5lU0nE19m4R5wlZRRFasWIFIiIioFQqkZCQgIyMjG6PPXr0KG6++WZERERAIpFg+fLlFsvVPiw2ItQL7pwfREREZFYTOy6jF2+ekOhFaN26dUhJScGiRYtw4MABxMTEYMaMGSgrK+vy+IaGBkRGRuK1115DUFCQRbO1XzY/fgDPBhEREZlb+8KK+wqr0NyqEyWD6EXorbfewv3334958+Zh2LBhWLlyJVxdXbFmzZoujx87dizeeOMN3H777VAoFBbN1n7FGBdSJCIiMr9Bge7wd1egqUWPrKJqUTKIWoS0Wi0yMzORlJTUcZ9UKkVSUhLS09PN9jrNzc3QaDSdbpdSUtOI05UNkEqA+AiuH0RERGRuEomk4+qxPXnizBMStQhVVFRAp9NBpVJ1ul+lUqG0tNRsr7NkyRJ4eXl13MLCwi75mIy2s0HDQ7zgoXQ2WxYiIiL6w8SBbUVIpHlCog+N9YXU1FTU1NR03IqLiy/5mP2F5wHwbBAREZElta/Td+hMDZpa+n6ekKiXQvn7+0Mmk0GtVne6X61Wm3UitEKhMHo+0f7ThiI0lpfNExERWUy4nyv83RWoqGvG4bM1ff65K+oZIblcjri4OKSlpXXcp9frkZaWhsTERNFyaZpacLzUMI+I+4sRERFZjkQiwdi20Zd9hVV9/vqiD42lpKRg1apV+OSTT5CTk4P58+ejvr4e8+bNAwDMmTMHqampHcdrtVpkZ2cjOzsbWq0WZ8+eRXZ2NvLy8syWKauoGoIA9Pd1RaCn0mzPS0RERBdqX7S4fVpKXxJ9lcDk5GSUl5dj4cKFKC0tRWxsLDZt2tQxgbqoqAhS6R997dy5cxg9enTH18uWLcOyZcswZcoUbNu2zSyZMtsaKc8GERERWV77GaH9hVXQ6wVIpZI+e23RixAALFiwAAsWLOjy1/5abiIiIiAIgkXz7OuYKM35QURERJY2LNgTrnIZNE2tyC2rw+Agjz57bdGHxqxNi06P7OJqALxijIiIqC84yaQY3d8bQN/PE2IR+otj5zRobNHBU+mEgQHuYschIiJyCPHh7fOEWIRE1X7ZfHyEb5+OURIRETmy9svm2z+H+wqL0F9knjY00ThOlCYiIuozsf29IZNKcOZ8I0pqGvvsdVmE/kQQhI6J0lxIkYiIqO+4K5wwLNgTQN9eRs8i9CfFVY0or22Gs0yCUf28xI5DRETkUNpHY/pynhCL0J+0z1QfGeoFpbNM5DRERESOpX00Zh/PCInjzxOliYiIqG+1L1tzvFQDTVNLn7wmi9CfcKI0ERGReFSeSvT3dYVeMGx31RdYhNpUN2hxUl0HgFtrEBERiSU+om/nCbEItTlQZBgWi/R3g5+7QuQ0REREjumPeUIsQn1qf8f+YjwbREREJJb2DVizi6vRotNb/PVYhNp0FKFwTpQmIiISS1SAO3xcndHUosfRcxqLvx6LEIDmVh0OnqkGwDNCREREYpJIJIjrw33HWIRg2Gi1uVUPXzc5Bvi7iR2HiIjIobUPj2WcYhHqEyfVtQCA4SGekEi40SoREZGYRvc3FCEOjfWR3LbL5gcGuouchIiIiAa1fR6frW5EfXOrRV+LRQhAXrmhCA0K9BA5CREREfm4yeHftpRNfttntKWwCIFnhIiIiKxN+1mh9s9oS3H4ItSgbcXZ6kYAf/zQiYiISFyDVG1FqIxFyKLyy+oBAH5ucvi4yUVOQ0RERMAfJyfyymot+joOX4Tyyg0/YA6LERERWY+BbfN2eUbIwvLKOD+IiIjI2rQPjRVVNaCpRWex13H4IsSJ0kRERNbHz00OH1dnCIJlrxxz+CLES+eJiIisj0Qi6fhszrPg8JhDFyFtqx6nKxsA8IwQERGRtRmosvwl9A5dhIqq6qHTC/BQOEHlqRA7DhEREf1Jx1pCFrxyzKGLUPul81GB7txjjIiIyMoM6oMrxxy6CBWUG4oQh8WIiIisT/uVY6crG9Dcapkrxxy7CFW0T5RmESIiIrI2gR4KeCidoNMLKKxosMhrOHQRar8cj2eEiIiIrI/hyjHLzhNy6CJ0qu2KMV46T0REZJ065glZ6Moxhy5CLa16KJykCPVxETsKERERdaF9npCl1hJy6CIEAFEB7pBJecUYERGRNRrIoTHL4vwgIiIi6zVIZRgaO1VRjxad3uzPzyLEIkRERGS1QryUcJPL0KITUFRl/ivHHL4I8dJ5IiIi6yWRSDpOWhRYYPNVhy9CPCNERERk3Qa2XTnWviOEOTl0EXKSShDu5yZ2DCIiIrqI9ivHeEbIzMJ8XSB3cugfARERkdVrn8aSX8EzQmYVGcCzQURERNaufVHFAhYh84ry54rSRERE1i7UxwVKZylaWu308vkVK1YgIiICSqUSCQkJyMjIuOjxX331FYYMGQKlUomRI0di48aNJr0uzwgRERFZP5lUgqgAy1zcJHoRWrduHVJSUrBo0SIcOHAAMTExmDFjBsrKyro8fs+ePbjjjjtw7733IisrC7NmzcKsWbNw5MgRo1+bRYiIiMg2WGq5G4kgCIJFnrmHEhISMHbsWLz77rsAAL1ej7CwMDz66KN49tlnLzg+OTkZ9fX1+PHHHzvuGz9+PGJjY7Fy5coevaZGo4GXlxdKy6ug8vcxzzdCREREFrNiax5e/z4LxctvQ01NDTw9Pc3yvE5meRYTabVaZGZmIjU1teM+qVSKpKQkpKend/mY9PR0pKSkdLpvxowZ+O6777p9nebmZjQ3N3d8XVNTAwBoaaqHRiPrxXdAREREfSHYVYC+2bCytDnP4YhahCoqKqDT6aBSqTrdr1KpcPz48S4fU1pa2uXxpaWl3b7OkiVL8NJLL11wf1hYmAmpiYiISEyVlZXw8vIyy3OJWoT6SmpqaqezSNXV1QgPD0dRUZHZfpBkGo1Gg7CwMBQXF5vtNCeZhu+F9eB7YV34fliPmpoa9O/fH76+vmZ7TlGLkL+/P2QyGdRqdaf71Wo1goKCunxMUFCQUccDgEKhgEKhuOB+Ly8v/qa2Ep6ennwvrATfC+vB98K68P2wHlKp+a71EvWqMblcjri4OKSlpXXcp9frkZaWhsTExC4fk5iY2Ol4APj111+7PZ6IiIioO6IPjaWkpGDu3LmIj4/HuHHjsHz5ctTX12PevHkAgDlz5iA0NBRLliwBADz++OOYMmUK3nzzTVxzzTVYu3Yt9u/fjw8++EDMb4OIiIhskOhFKDk5GeXl5Vi4cCFKS0sRGxuLTZs2dUyILioq6nQKbMKECfjf//6Hf/7zn3juuecwaNAgfPfddxgxYkSPX1OhUGDRokVdDpdR3+J7YT34XlgPvhfWhe+H9bDEeyH6OkJEREREYhF9ZWkiIiIisbAIERERkcNiESIiIiKHxSJEREREDstui9CKFSsQEREBpVKJhIQEZGRkXPT4r776CkOGDIFSqcTIkSOxcePGPkpq/4x5L1atWoVJkybBx8cHPj4+SEpKuuR7Rz1n7J+LdmvXroVEIsGsWbMsG9CBGPteVFdX45FHHkFwcDAUCgWio6P595SZGPteLF++HIMHD4aLiwvCwsLw5JNPoqmpqY/S2q8dO3bguuuuQ0hICCQSyUX3EG23bds2jBkzBgqFAgMHDsTHH39s/AsLdmjt2rWCXC4X1qxZIxw9elS4//77BW9vb0GtVnd5/O7duwWZTCYsXbpUOHbsmPDPf/5TcHZ2Fg4fPtzHye2Pse/FnXfeKaxYsULIysoScnJyhHvuuUfw8vISzpw508fJ7Y+x70W7U6dOCaGhocKkSZOEG264oW/C2jlj34vm5mYhPj5emDlzprBr1y7h1KlTwrZt24Ts7Ow+Tm5/jH0vPv/8c0GhUAiff/65cOrUKeGXX34RgoODhSeffLKPk9ufjRs3Cs8//7ywfv16AYDw7bffXvT4goICwdXVVUhJSRGOHTsmvPPOO4JMJhM2bdpk1OvaZREaN26c8Mgjj3R8rdPphJCQEGHJkiVdHn/bbbcJ11xzTaf7EhIShAcffNCiOR2Bse/FX7W2tgoeHh7CJ598YqmIDsOU96K1tVWYMGGC8OGHHwpz585lETITY9+L9957T4iMjBS0Wm1fRXQYxr4XjzzyiHDFFVd0ui8lJUWYOHGiRXM6mp4UoaeffloYPnx4p/uSk5OFGTNmGPVadjc0ptVqkZmZiaSkpI77pFIpkpKSkJ6e3uVj0tPTOx0PADNmzOj2eOoZU96Lv2poaEBLS4tZN9hzRKa+F4sXL0ZgYCDuvffevojpEEx5L77//nskJibikUcegUqlwogRI/Dqq69Cp9P1VWy7ZMp7MWHCBGRmZnYMnxUUFGDjxo2YOXNmn2SmP5jrs1v0laXNraKiAjqdrmNl6nYqlQrHjx/v8jGlpaVdHl9aWmqxnI7AlPfir5555hmEhIRc8JudjGPKe7Fr1y6sXr0a2dnZfZDQcZjyXhQUFOC3337DXXfdhY0bNyIvLw8PP/wwWlpasGjRor6IbZdMeS/uvPNOVFRU4LLLLoMgCGhtbcVDDz2E5557ri8i059099mt0WjQ2NgIFxeXHj2P3Z0RIvvx2muvYe3atfj222+hVCrFjuNQamtrMXv2bKxatQr+/v5ix3F4er0egYGB+OCDDxAXF4fk5GQ8//zzWLlypdjRHM62bdvw6quv4j//+Q8OHDiA9evX46effsLLL78sdjQykd2dEfL394dMJoNare50v1qtRlBQUJePCQoKMup46hlT3ot2y5Ytw2uvvYYtW7Zg1KhRlozpEIx9L/Lz81FYWIjrrruu4z69Xg8AcHJywokTJxAVFWXZ0HbKlD8XwcHBcHZ2hkwm67hv6NChKC0thVarhVwut2hme2XKe/HCCy9g9uzZuO+++wAAI0eORH19PR544AE8//zznfbGJMvq7rPb09Ozx2eDADs8IySXyxEXF4e0tLSO+/R6PdLS0pCYmNjlYxITEzsdDwC//vprt8dTz5jyXgDA0qVL8fLLL2PTpk2Ij4/vi6h2z9j3YsiQITh8+DCys7M7btdffz0uv/xyZGdnIywsrC/j2xVT/lxMnDgReXl5HWUUAE6ePIng4GCWoF4w5b1oaGi4oOy0F1SBW3f2KbN9dhs3j9s2rF27VlAoFMLHH38sHDt2THjggQcEb29vobS0VBAEQZg9e7bw7LPPdhy/e/duwcnJSVi2bJmQk5MjLFq0iJfPm4mx78Vrr70myOVy4euvvxZKSko6brW1tWJ9C3bD2Pfir3jVmPkY+14UFRUJHh4ewoIFC4QTJ04IP/74oxAYGCi88sorYn0LdsPY92LRokWCh4eH8MUXXwgFBQXC5s2bhaioKOG2224T61uwG7W1tUJWVpaQlZUlABDeeustISsrSzh9+rQgCILw7LPPCrNnz+44vv3y+X/84x9CTk6OsGLFCl4+/2fvvPOO0L9/f0Eulwvjxo0T9u7d2/FrU6ZMEebOndvp+C+//FKIjo4W5HK5MHz4cOGnn37q48T2y5j3Ijw8XABwwW3RokV9H9wOGfvn4s9YhMzL2Pdiz549QkJCgqBQKITIyEjhX//6l9Da2trHqe2TMe9FS0uL8OKLLwpRUVGCUqkUwsLChIcfflg4f/583we3M1u3bu3y7//2n//cuXOFKVOmXPCY2NhYQS6XC5GRkcJHH31k9OtKBIHn8oiIiMgx2d0cISIiIqKeYhEiIiIih8UiRERERA6LRYiIiIgcFosQEREROSwWISIiInJYLEJERETksFiEiIiIyGGxCBEREZHDYhEiIiIih8UiRER2ISIiAsuXL+90X2xsLF588UVR8hCRbWARIiIiIofFIkREREQOi0WIiIiIHBaLEBHZLZ1OJ3YEIrJyLEJEZDfUanXH/7e0tKC4uFjENERkC1iEiMhurFmzBlu2bEFubi6efPJJ1NTUID8/v1NBIiL6MxYhIrIb1113HR577DGMHDkSVVVVeOWVV7B+/Xps2bJF7GhEZKUkgiAIYocgIuqtiIgIPPHEE3jiiSfEjkJENoRnhIiIiMhhsQgRERGRw+LQGBERETksnhEiIiIih8UiRERERA6LRYiIiIgcFosQEREROSwWISIiInJYLEJERETksFiEiIiIyGGxCBEREZHDYhEiIiIih/X/NXGxAZxItDoAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# ================================================\n",
    "# 2.2 離散確率分布\n",
    "# ================================================\n",
    "\n",
    "b1 = Bern(0.5)\n",
    "print([b1() for _ in range(20)])\n",
    "\n",
    "b2 = Bern(0.9)\n",
    "print([b2() for _ in range(20)])\n",
    "\n",
    "sampler = SimpleSampler()\n",
    "print(f\"E[Bern(x|0.5)] = {sampler(b1, 10000)}\")\n",
    "print(f\"E[Bern(x|0.9)] = {sampler(b2, 10000)}\")\n",
    "\n",
    "H = Entropy()\n",
    "print(f\"H[Bern(x|0.5)] = {H(b1)}\")\n",
    "print(f\"H[Bern(x|0.9)] = {H(b2)}\")\n",
    "\n",
    "from matplotlib import pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "x = np.linspace(0, 1, 100)\n",
    "y = [H(Bern(mu), 100000) for mu in x]\n",
    "\n",
    "fig, ax = plt.subplots()\n",
    "ax.plot(x, y)\n",
    "ax.set_ylim(0, 0.7)\n",
    "ax.set_xlim(0, 1.0)\n",
    "ax.set_xlabel(\"μ\")\n",
    "ax.set_ylabel(\"H[Bern(x|μ)]\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}