14_マルコフ連鎖.ipynb

{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# マルコフ連載\n",
    "\n",
    "実数値確率変数の列$X=(X_n)_{(n=0,1,2,...)}$と集合$B \\subset R$に対して、\n",
    "$$\n",
    "P(X_{n+1} \\in B | X_n,X_{n-1},...,X_1,X_0) = P(X_{n+1} \\in B | X_n)\n",
    "$$\n",
    "となるとき、これを確率変数$X$の**マルコフ性**といい、このような性質を持つ$X$を**マルコフ連鎖**という。  \n",
    "要するに次の遷移の確率が完全に現在の状態にのみ依存し、それ以前の状態には依存しないことをいう。  \n",
    "$n$を時間の変数とみなすとき、$X_n$は時点$n$のおける何らかの\"状態\"を表すと考え、$X_n$がとる値の空間のことを**状態空間**ということがある。  \n",
    "また、ある時点$n$から$m$ステップ先の確率を$m$ステップ**推移確率**といい、特に$p_m^{(n)}$が時点$n$によらないとき、マルコフ連鎖$X$は**斉時的**であるという。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "天気予測結果：['sunny', 'sunny', 'sunny', 'sunny', 'sunny']\n"
     ]
    }
   ],
   "source": [
    "# 具体例(天気の予測)\n",
    "import numpy as np\n",
    "\n",
    "states = [\"sunny\", \"rainy\"]\n",
    "# S=sunny, R=rainnyでSSはsunny->sunnyとなる状態遷移\n",
    "transitions = [[\"SS\", \"SR\"],[\"RS\", \"RR\"]]\n",
    "# transitionsの各要素に対応する確率行列\n",
    "T = [[0.9, 0.1],[0.5, 0.5]]\n",
    "\n",
    "# 現在の天気\n",
    "today_weather = \"sunny\"\n",
    "\n",
    "# 予測したい日数\n",
    "forcast_days = 5\n",
    "\n",
    "# 現在の天気と上の遷移確率に基づいて推測\n",
    "weather_list = []\n",
    "n_days = 0\n",
    "prob = 1.0\n",
    "while n_days != forcast_days:\n",
    "    if today_weather == \"sunny\":\n",
    "        change = np.random.choice(transitions[0], p=T[0])\n",
    "        if change == \"SS\":\n",
    "            prob = prob * T[0][0]\n",
    "            weather_list.append(states[0])\n",
    "        else:\n",
    "            prob = prob * T[0][1]\n",
    "            weather_list.append(states[1])            \n",
    "    else:\n",
    "        change = np.random.choice(transitions[1], p=T[1])\n",
    "        if change == \"RS\":\n",
    "            prob = prob * T[1][0]\n",
    "            weather_list.append(states[0])\n",
    "        else:\n",
    "            prob = prob * T[1][1]\n",
    "            weather_list.append(states[1])        \n",
    "    n_days+=1\n",
    "\n",
    "print(f\"天気予測結果：{weather_list}\")"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 離散的で有限な状態空間\n",
    "\n",
    "応用上は、離散的な状態空間$S=\\{x_1,x_2,...\\}$をもつ斉時的マルコフ連鎖が重要である。  \n",
    "特に状態空間$S$が有限個の要素からなるとき、$X$を**有限マルコフ連鎖**という。  \n",
    "以下、一般性を失うことなく$S=\\{x_1,x_2,...\\}$とし、$X=(X_n)_{n=0,1,2,...}$を$S$上の斉時的有限マルコフ連鎖とする。  \n",
    "このとき、状態$i$から状態$j$への$m$ステップ推移確率を\n",
    "$$\n",
    "p_m(i,j) = P(X_{n+m}=j | X_n=i), \\quad i,j \\in S\n",
    "$$\n",
    "とおく。これを$i,j$成分において作った$N \\times N$行列\n",
    "$$\n",
    "Q(m) = \n",
    "\\begin{pmatrix}\n",
    "p_m(1,1) & p_m(1,2) & \\dots & p_m(1,N) \\\\\n",
    "p_m(2,1) & p_m(2,2) & \\dots & p_m(2,N) \\\\\n",
    "\\vdots & \\vdots & \\ddots & \\vdots \\\\\n",
    "p_m(N,1) & p_m(N,2) & \\dots & p_m(N,N)\n",
    "\\end{pmatrix}\n",
    "$$\n",
    "をマルコフ連鎖$X$の$m$ステップ**推移確率行列**という。  \n",
    "また、任意の行列成分の和が$1$で、各成分が$0$以上となるものを**確率行列**という。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "じゃんけんの推移確率行列：\n",
      "[[0.4 0.3 0.3]\n",
      " [0.6 0.1 0.3]\n",
      " [0.5 0.2 0.3]]\n"
     ]
    }
   ],
   "source": [
    "# 具体例(じゃんけんの推移確率行列(m=2))\n",
    "transition_matrix = np.zeros([3,3])\n",
    "\n",
    "# 0: rock, 1:scissors, 2:paper\n",
    "transition_matrix[0,0] = 0.4\n",
    "transition_matrix[0,1] = 0.3\n",
    "transition_matrix[0,2] = 0.3\n",
    "transition_matrix[1,0] = 0.6\n",
    "transition_matrix[1,1] = 0.1\n",
    "transition_matrix[1,2] = 0.3\n",
    "transition_matrix[2,0] = 0.5\n",
    "transition_matrix[2,1] = 0.2\n",
    "transition_matrix[2,2] = 0.3\n",
    "print(f\"じゃんけんの推移確率行列：\\n{transition_matrix}\")"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 推移確率の性質\n",
    "\n",
    "$X_n$に対して、$p_n(k):=P(X_n=k)(k \\in S)$と書くとき、\n",
    "$$\n",
    "π_n = (p_n(1), p_n(2),...,p_n(N))\n",
    "$$\n",
    "なる$N$次元ベクトルを時点$n$における$X$の**状態ベクトル**という。  \n",
    "特に$π_0$を$X$の**初期分布**という。また $\\sum_{k=1}^{N}p_n(k) = 1$となる。  \n",
    "このように、成分の総和が$1$となるベクトルを**確率ベクトル**という。  \n",
    "推移確率行列と状態確率ベクトルについては以下が成り立つ。  \n",
    "(1)任意の自然数$l,m$に対して、\n",
    "$$\n",
    "Q(m+l) = Q(m)Q(l) = Q^mQ^l\n",
    "$$\n",
    "(2)任意の$n=0,1,2,...$に対して、\n",
    "$$\n",
    "π_n = π_0Q^n\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Q(m+l):\n",
      "[[0.47488 0.22512 0.3    ]\n",
      " [0.4752  0.2248  0.3    ]\n",
      " [0.47504 0.22496 0.3    ]]\n",
      "\n",
      "Q(m)*Q(l) = Q^m*Q^l:\n",
      "[[0.47488 0.22512 0.3    ]\n",
      " [0.4752  0.2248  0.3    ]\n",
      " [0.47504 0.22496 0.3    ]]\n"
     ]
    }
   ],
   "source": [
    "# 具体例(じゃんけんの推移確率行列)\n",
    "m = 2\n",
    "l = 3\n",
    "\n",
    "m_transition_matrix = np.linalg.matrix_power(transition_matrix, m)\n",
    "l_transition_matrix = np.linalg.matrix_power(transition_matrix, l)\n",
    "\n",
    "# Q(m+l)\n",
    "ml_transition_matrix = np.linalg.matrix_power(transition_matrix, m+l)\n",
    "print(f\"Q(m+l):\\n{ml_transition_matrix}\\n\")\n",
    "\n",
    "# Q(m)*Q(l) = (Q**m)*(Q**l)\n",
    "ml_transition_matrix = np.dot(m_transition_matrix, l_transition_matrix)\n",
    "print(f\"Q(m)*Q(l) = Q^m*Q^l:\\n{ml_transition_matrix}\")"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 定常分布\n",
    "\n",
    "状態確率ベクトル$π_n$において$n \\rightarrow  \\infty$としたときの極値が存在する場合を考える。\n",
    "$$\n",
    "π:=\\lim_{n \\to  \\infty} π_n \n",
    "$$\n",
    "このとき両辺で$n \\rightarrow  \\infty$とすると、\n",
    "$$\n",
    "π = πQ\n",
    "$$\n",
    "が成り立つ。つまり、$π$は行列$Q$の固有値$1$に対する左固有ベクトルである。  \n",
    "この意味は、マルコフ連鎖$X$の初期分布が$π$であったとすると、$1$回の推移によって状態変化しないということであり、長い時間のあと状態確率が$π$に収束してしまうと、そのあと何回推移してもその状態確率は変化しない。  \n",
    "このような状態確率を$X$の**定常分布**といい、初期状態が定常分布なマルコフ連鎖のことを**定常マルコフ連鎖**という。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAACTgAAAIwCAYAAACoB5bfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABnKklEQVR4nOz9ebyWdYE//r8OHvbtAIqJoFhIuPXNBZLUMdxqNLdyxiwT3DNzXBhN+UwulSGZU35mmszBJZ1Smz6pETaaBu6OQhRaKO6hlYoKB9mX+/cHP24hOMB1OMvF4fl8PHg8rvvc7+v9ft+Hi3O/Or287ppKpVIJAAAAAAAAAABACbVr7Q0AAAAAAAAAAAA0RMEJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBSdgs7N8+fJMnz49N9xwQ84666zss88+6dChQ2pqalJTU5NPfOITjZ77gQceyEknnZTBgwena9eu6d27dz7ykY/kwgsvzLPPPtt0LwIAoAlsbrloxowZufDCC/ORj3wkvXv3TteuXTN48OCMHDkyDzzwQKP3CgCwMW6++eZqTtrYP6eddtpGzy8/AQCbE79XSmbPnp3vfOc7+fjHP57tttsunTp1yo477pjDDz88t956a5YuXdqoeYHmUVOpVCqtvQmAjXXXXXflC1/4QhYsWNDgmAMPPDCTJ08uNG99fX3OOOOM3HHHHQ2Oad++fa644opccsklheYGAGgOm1suuvLKK3PFFVes9xdDJ5xwQn74wx+me/fuhfYMALAxbr755px88smFzjn11FMzfvz49Y6RnwCAzY3fKyW//OUvc8opp+Stt95qcMxee+2V2267LYMHD97ovQLNp7a1NwBQxJw5c9Ybthpj6dKlOfbYY/Ob3/ym+rXdd989e+21VxYtWpSHH344f/nLX7J06dKMGTMmS5cuzaWXXtqkewAAKGpzykWXXnppvvGNb1Qfb7fddjnggAPSqVOnTJ06NX/4wx+SJLfddlvefvvtTJw4MbW1/ucqANB8hgwZkoMPPniD4z7+8Y+v93n5CQDYHG3pv1e67777cuyxx2bZsmVJki5duuTggw/ONttskxdffDEPPfRQKpVKfvvb3+bggw/O//7v/6Zfv36N+bYATakCsBm56aabKkkq2267beXTn/505Yorrqjcc889lXPPPbeSpJKkcuCBBxaa82tf+1r13E6dOlVuu+22NZ5fvHhx5cILL6yOqampqUyePLkJXxUAQHGbSy66//77q+OTVC688MLK4sWL1xjzk5/8pNKpU6fqmCuuuKLQvgEANsaq/JSkMnLkyCaZU34CADZHW/LvlWbPnl2pq6urjj/44IMrb7311hpjfve731V22GGH6piDDjqowHcCaC4+og7YrPz1r3/NkiVLssMOO6zx9csvvzxXXHFFkmK3zHzzzTfzwQ9+MPPnz0+SXHfddTnzzDPXOfZzn/tc9Zaaw4cPz2OPPdbIVwEAsOk2l1w0bNiwPPXUU9XzbrvttnWOu+6663LWWWclSbp3756XXnopW2+99UbtHQBgY6z+EXUjR47MzTffvEnzyU8AwOZqS/690kUXXZSrr746SfKhD30o06dPT5cuXdYaN3369Oyzzz7Vj8W79957c9hhhzW4V6D5tWvtDQAU8YEPfGCtsLUpfvSjH1XD1uDBg3PGGWc0OPbb3/522rVb+WPz8ccfz7Rp05psHwAARW0Oueipp56q/hKqXbt2+fa3v93gnGeeeWZ23nnnJMm8efNy6623Nup1AAC0FPkJANhcbam/V1q6dGn+8z//s/r461//+jrLTUnykY98JCNHjqw+/v73v9/g+kDLUHACtmh33XVX9XjUqFGpqalpcOwOO+yQgw46qPr4zjvvbM6tAQC0qObIRavPecghh2TAgAENzllTU7PGL41kLQCg7OQnAICVNpdcNHny5MyZMyfJyjs9ffazn21wzmTla1nlvvvuq5a4gNah4ARssRYtWpQnnnii+vgTn/jEBs8ZMWJE9fg3v/lNc2wLAKDFNVcumjRpUqPnfOyxx7J48eINngMA0BrkJwCAlTanXLT6nMOHD0/Hjh3XO+ewYcOqd3hatGhRHn/88Q3uA2g+ta29AYDW8txzz2XFihVJVja799xzzw2es9dee1WPZ8yY0Wx7AwBoSc2Vi1b/+urjG7L6usuXL8/MmTOzxx57bPA8AICi5syZk//+7//OH/7wh8ydOzc9evRIv379Mnz48Oyxxx7rvetAIj8BAKyyOeWionO2b98+e+yxR/73f/+3ev4hhxyywfOA5qHgBGyxnnvuuepx375906lTpw2es/rnEb/zzjt56623ss022zTL/gAAWkpz5KI333yzesvvJNlxxx03OGfnzp2zzTbb5K233kqSPPvss/4POgCgWdx99925++671/nczjvvnK9+9as55ZRTGiw6yU8AACttTrlo9b1uzJyr9rqq4PTss89u1DlA8/ARdcAW6+23364eb7vttht1zgc+8IE1Hr/zzjtNuicAgNbQHLlo9TkbO6+sBQC0hueffz6nnXZajjrqqMyfP3+dY+QnAICVNqdctKl7lbWgdSk4AVus9957r3rcuXPnjTrnb8etPgcAwOaqOXLR3z5uzLyyFgDQ1HbYYYeMHj0699xzT2bNmpVFixZl/vz5ee655/If//EfGTJkSHXsL3/5y3z+85+vfuTK6uQnAICVNqdctKl7lbWgdSk4AVusRYsWVY87dOiwUed07NhxjccLFy5s0j0BALSG5shFq8/Z2HllLQCgKR1zzDF5+eWX853vfCd///d/n/79+6djx47p0qVLBg8enLPOOiu///3vc/LJJ1fP+cUvfpGf/OQna80lPwEArLQ55aJN3ausBa1LwQnYYq3+GcBLlizZqHMWL168xuONbXcDAJRZc+Si1eds7LyyFgDQlOrq6tKu3fp/Jd6hQ4eMHz8+BxxwQPVr48aNW2uc/AQAsNLmlIs2da+yFrQuBSdgi9WtW7fq8cY2rv923OpzAABsrpojF/3t48bMK2sBAK2hXbt2ueyyy6qPn3nmmbz22mtrjJGfAABW2pxy0abuVdaC1qXgBGyx+vTpUz1+4403Nuqcv/71r2s87t27d5PuCQCgNTRHLlp9zsbOK2sBAK3l7/7u79K+ffvq4xkzZqzxvPwEALDS5pSLNnWvsha0LgUnYIv14Q9/uHr85ptvrvV5vuvypz/9qXrcu3fvbLPNNs2yNwCAltQcuahv376pq6urPn711Vc3OOeiRYvy1ltvVR8PGTJkg+cAADSH9u3bZ+utt64+nj179hrPy08AACttTrlo9b1uzJx/u1dZC1qXghOwxfrwhz+cdu1W/hisVCr53e9+t8Fzfvvb31aPd9lll+baGgBAi2quXLT616dNm1Zozq222iqDBw/e4DkAAM1l/vz51eOuXbuu8Zz8BACw0uaUi4rOuWzZsjz99NMb3CvQMhScgC1Wp06dsu+++1YfT548eYPnPPjgg9Xjgw46qDm2BQDQ4porF40YMaLRc3784x9Px44dN3gOAEBzeOmll1JfX1993K9fvzWel58AAFbanHLR6nM+/vjjWbJkyXrnfOqpp7JgwYIkK1/n8OHDN7gPoPkoOAFbtGOOOaZ6fPPNN6937KxZs/LAAw+s81wAgM1dc+Si1b9+//3357XXXlvvvKuvK2sBAK3pxhtvrB737NkzH/3oR9caIz8BAKy0ueSiT3ziE+nZs2eSpL6+Pj//+c83es5DDz10rbt6Ai1LwQnYoo0cObIaRp577rmMHz++wbFf/epXs3z58iTJ8OHDs9dee7XIHgEAWkJz5KKhQ4dm6NChSZLly5fn4osvbnDO66+/PjNnzkySdO/ePSeddFKjXgcAwLq89957Gz32scceyzXXXFN9/LnPfS61tbVrjZOfAABW2lxyUfv27XP66adXH1966aVZuHDhOsc+88wzaxSczj777AbXB1qGghOwRevbt28uuOCC6uN/+qd/yk9/+tM1xixdujQXX3xxbrvtturXxo4d22J7BABoCc2Vi1Z//sc//nEuvvjiLF26dI0xP/3pT3PeeedVH//zP/9ztt5668a8DACAdfrZz36WYcOG5ZZbbsncuXPXOWbRokX5v//3/+aQQw7JokWLkiR1dXW57LLL1jlefgIAWGlzykUXX3xx6urqkiTPP/98jj766Lz99ttrjJk+fXqOPPLI6kfYjRgxIp/85CfXu1eg+dVUKpVKa28CoIjDDz88f/7zn9f42l//+te88cYbSZKuXbtm0KBBa513zz33pF+/fmt9fenSpfnUpz6V3/zmN9Wv7bHHHtlrr72yaNGiPPTQQ/nLX/5Sfe6KK67IpZde2lQvBwCg0TaXXPS1r30t3/zmN6uP+/XrlwMOOCCdOnXK1KlT88wzz1SfO/TQQ3PPPfes8y4JAACNdfPNN+fkk09OktTW1mbIkCEZMmRIevXqleXLl+f111/P448/nvr6+uo5nTt3zv/8z//k7/7u7xqcV34CADZXW/Lvle699958+tOfzrJly5IkXbp0ySGHHJJtttkmL774Yh588MGsqlFsv/32efLJJ9f5moGW1aiC05tvvpknn3wyTz75ZJ566qk89dRT1VbjyJEjN/i5mo1x22235aabbsr06dMzZ86cbLvttjnggANy9tlnZ/jw4U2+HlBeAwcOzKuvvlr4vJdffjkDBw5c53Nz587NGWecsVabfHXt27fP5ZdfnjFjxhReGyCRoYCmt7nkokqlkiuvvDJf//rX1/qv7Fb3uc99Lj/84Q/To0ePjZoXaPvkJ6CprF5w2hjDhg3LzTffnF122WWDY+UnoGxkKGBjbOm/V5owYUJOOeWUzJ49u8Exe+65Z2677bZ8+MMf3qg5gebVqIJTTU1Ng881dTBauHBhjjvuuNxzzz3rfL5du3a59NJLG7xNMND2NEfgWuX+++/Pj370ozz++OP5y1/+kvbt22fAgAH55Cc/mVNPPXWjfqkF0BAZCmhqm1sumjFjRsaPH5/77rsvs2bNytKlS7Pddttl+PDhGTlyZA455JDCcwJtm/wENJXFixdn6tSpeeyxx/LEE0/kxRdfzOzZs/P2229nxYoV6dmzZ3baaafsu+++Oe6447L//vsXXkN+AspChgI2ht8rJW+99VZuuumm3HnnnXnppZeqBc3ddtstJ5xwQk444YS0b9++8LxA89jkgtMOO+yQIUOG5L777kvS9MHohBNOyO23355k5WdbnnvuuenXr1+efvrpfOtb38qLL76YJPnhD3+YM844o8nWBQBoajIUAEAx8hMAQHEyFADQFjWq4HTZZZdl6NChGTp0aLbddtu88sor2WmnnZI0bTD6zW9+k4MPPjhJcuSRR+bOO+/MVlttVX1+9uzZ2XvvvfOnP/0pdXV1eemll9KrV68mWRsAoKnJUAAAxchPAADFyVAAQFvUrjEnXXHFFfn0pz+dbbfdtqn3s4bvfOc7SZLa2tr8x3/8xxqhKEm23nrrjBs3LkkyZ86cjB8/vln3AwCwKWQoAIBi5CcAgOJkKACgLWpUwaklzJs3Lw888ECS5JBDDkn//v3XOe4zn/lMevTokSS58847W2x/AABlJEMBABQjPwEAFCdDAQAtrbQFp6eeeipLlixJkhx44IENjuvQoUP23Xff6jlLly5tkf0BAJSRDAUAUIz8BABQnAwFALS00hac/vjHP1aPhwwZst6xq55ftmxZnn/++WbdFwBAmclQAADFyE8AAMXJUABASyttwem1116rHjd0W8tVBgwYUD2eNWtWs+0JAKDsZCgAgGLkJwCA4mQoAKCl1bb2Bhoyb9686nG3bt3WO7Zr167V4/fee6/BcYsXL87ixYurj1esWJF33nknffr0SU1NzSbsFgAoq0qlknnz5qVfv35p16603e4mI0MBAE1hS8pQ8hMA0FRkqHWToQCAhhTJT6UtOC1atKh63KFDh/WO7dixY/V44cKFDY4bO3Zsrrjiik3fHACw2Zk1a9YG/2uytkCGAgCa0paQoeQnAKCpyVBrkqEAgA3ZmPxU2oJTp06dqsdLlixZ79jV29ydO3ducNwll1ySCy64oPp47ty52WGHHTJr1qz06NFjE3YLAJRVfX19BgwYkO7du7f2VlqEDAUANIUtKUPJTwBAU5Gh1k2GAgAaUiQ/lbbgtPrm13e7yiSZP39+9Xh9t8Hs2LHjGi3xVXr06CEYAUAbt6XcxlqGAgCa0paQoeQnAKCpyVBrkqEAgA3ZmPxU2g8AXv3WU6+99tp6x86aNat6PGDAgGbbEwBA2clQAADFyE8AAMXJUABASyttwWnXXXetHj/77LPrHbvq+dra2uy8887Nui8AgDKToQAAipGfAACKk6EAgJZW2oLT0KFD06FDhyTJgw8+2OC4JUuW5Iknnqie0759+xbZHwBAGclQAADFyE8AAMXJUABASyttwal79+45+OCDkyT3339/g7e3/PnPf576+vokybHHHtti+wMAKCMZCgCgGPkJAKA4GQoAaGmtVnC6+eabU1NTk5qamlx++eXrHPPP//zPSZJly5bl7LPPzvLly9d4fvbs2fnqV7+aJKmrq8tpp53WrHsGAGhtMhQAQDHyEwBAcTIUAFA2tY056ZFHHskLL7xQfTx79uzq8QsvvJCbb755jfGjRo1q1OYOOuigfO5zn8vtt9+eX/ziFzn00ENz3nnnpV+/fnn66adz5ZVX5k9/+lOSZNy4cenVq1ej1gEAaAkyFABAMfITAEBxMhQA0BY1quA0fvz4/OhHP1rnc48++mgeffTRNb7W2GCUJDfeeGPq6+tzzz33ZNKkSZk0adIaz7dr1y5f+9rXcsYZZzR6DQCAliBDAQAUIz8BABQnQwEAbVGrfUTdxurcuXMmTpyYH//4xzn00EPTt2/fdOjQIQMGDMjnP//5PPLIIw3eGhMAYEslQwEAFCM/AQAUJ0MBAC2lplKpVFp7E62lvr4+PXv2zNy5c9OjR4/W3g4A0Ay83zc931MAaPu83zct308A2DJ4z29avp8A0PYVeb8v/R2cAAAAAAAAAACALZeCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApbXJBadXX301o0ePzpAhQ9K1a9f07t07Q4cOzdVXX50FCxY0xR7zyiuv5Ktf/Wr23nvv1NXVpX379undu3c+/vGP5+tf/3refPPNJlkHAKClyFAAAMXITwAAxclQAEBbUVOpVCqNPXnChAk58cQTU19fv87nBw8enIkTJ2bQoEGN3uCtt96aM888MwsXLmxwTO/evXP77bfn0EMPLTR3fX19evbsmblz56ZHjx6N3iMAUF5lfL+XoQCAsivb+738BABsDsr2ni9DAQBlV+T9vtF3cJo2bVqOP/741NfXp1u3brnyyivz2GOP5YEHHsjpp5+eJJk5c2aOOOKIzJs3r1FrPProoxk1alQWLlyYdu3a5eSTT85dd92VJ598Mj/72c9y5JFHJkneeeedHH300XnppZca+3IAAFqEDAUAUIz8BABQnAwFALQ5lUY64IADKkkqtbW1lccee2yt57/97W9XklSSVC677LJGrXHEEUdU5/j+97+/zjEXXHBBdczZZ59daP65c+dWklTmzp3bqP0BAOVXtvd7GQoA2ByU6f1efgIANhdles+XoQCAzUGR9/tGfUTdk08+mY997GNJkjPPPDPXXXfdWmNWrFiR3XffPTNmzEhdXV3efPPNtG/fvtA6vXv3zrvvvps+ffpk9uzZ6xwzd+7c1NXVJUn22muvTJ06daPnd2tLAGj7yvR+L0MBAJuLsrzfy08AwOakLO/5MhQAsLlo9o+ou+uuu6rHJ5988ronbtcuJ510UpJkzpw5mTRpUuF1lixZkiTZaaedGhzTs2fPbL311muMBwAoIxkKAKAY+QkAoDgZCgBoixpVcHrkkUeSJF27ds3ee+/d4LgDDzywevzoo48WXufDH/5wkuTll19ucEx9fX21Fb5qPABAGclQAADFyE8AAMXJUABAW9SogtOMGTOSJIMGDUptbW2D44YMGbLWOUV86UtfSpK8/fbb67x9ZpJ84xvfWGs8AEAZyVAAAMXITwAAxclQAEBb1HCqacCiRYuqTev+/fuvd2yvXr3StWvXzJ8/P7NmzSq8uVNOOSWPPPJIbrnllpx99tmZOnVqjjrqqGy33Xb505/+lFtvvbV6m83/83/+Tw455JDCawAAtAQZCgCgGPkJAKA4GQoAaKsKF5zmzZtXPe7WrdsGx68KRu+9917RpbLVVlvlRz/6UY488sh861vfyvjx4zN+/Pg1xowYMSJjxozZqFC0ePHiLF68uPq4vr6+8J4AABpDhgIAKEZ+AgAoToYCANqqRt3BaZUOHTpscHzHjh2TJAsXLiy6VJKVt8S85ZZb8vTTT6/z+ccffzw33HBDdtlll2y//fbrnWvs2LG54oorGrWPxhh48cQWW4vye+WqI1p7CwC0IhkKAKAY+QkAoDgZCgBoq9oVPaFTp07V4yVLlmxw/KqmdefOnYsulYcffjjDhw/PhAkTsv322+fWW2/NX//61yxZsiSzZs3K97///XTp0iW33357hg0blj/84Q/rne+SSy7J3Llzq38ac7tNAIDGkKEAAIqRnwAAipOhAIC2qvAdnLp371493pjbVc6fPz/Jxt0Gc3WLFy/OCSeckLlz5+YDH/hAnnjiiXzgAx+oPt+/f/98+ctfzoEHHph99tknf/7znzNy5MhMmTKlwTk7duxYbaIDALQkGQoAoBj5CQCgOBkKAGirGnUHpz59+iRJXnvttfWOfffdd6vBaMCAAYXW+Z//+Z+8/vrrSZJzzjlnjVC0ut122y0nnnhikmTq1Kn5/e9/X2gdAICWIEMBABQjPwEAFCdDAQBtVeGCU5LsuuuuSZIXXnghy5Yta3Dcs88+Wz3eZZddCq0xY8aM6vFee+213rF77733OtcEACgTGQoAoBj5CQCgOBkKAGiLGlVw2n///ZOsvG3l1KlTGxz34IMPVo/322+/QmvU1r7/6XnrC19JsnTp0nWeBwBQJjIUAEAx8hMAQHEyFADQFjWq4HTMMcdUj2+66aZ1jlmxYkVuueWWJEldXV1GjBhRaI2ddtqpevzwww+vd+zqAWz18wAAykSGAgAoRn4CAChOhgIA2qJGFZyGDRuWAw44IElyww035PHHH19rzDXXXFO9PeW5556b9u3br/H85MmTU1NTk5qamowaNWqt8w8++OB06dIlSfKDH/wgTz/99Dr38qtf/Sp33nlnkmT77bfPRz/60ca8JACAZidDAQAUIz8BABQnQwEAbVGjCk5Jcu2116Zz585ZtmxZDjvssIwdOzZPPPFEJk2alDPPPDMXXXRRkmTw4MEZPXp04fnr6upy8cUXJ0nmzZuXj3/84xkzZkwmTZqU3/3ud7n33nvz5S9/OUcddVRWrFiRJLnqqqvSrl2jXxIAQLOToQAAipGfAACKk6EAgLam0R90u+eee+aOO+7IiSeemPr6+owZM2atMYMHD87EiRPTvXv3Rq3xL//yL3nnnXdy7bXX5r333svYsWMzduzYtca1b98+3/rWt3LiiSc2ah0AgJYiQwEAFCM/AQAUJ0MBAG3NJtWkjzzyyEyfPj3nn39+Bg8enC5duqSuri777LNPxo0bl2nTpmXQoEGNnr+mpibf/e5389RTT+VLX/pSdt9993Tv3j1bbbVVevbsmb333jsXXHBBnnnmmfzzP//zprwUAIAWI0MBABQjPwEAFCdDAQBtSU2lUqm09iZaS319fXr27Jm5c+emR48eTT7/wIsnNvmcbL5eueqI1t4CwBapud/vt0S+pwDQ9nm/b1q+nwCwZfCe37R8PwGg7Svyfu+DbgEAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQ2ueD06quvZvTo0RkyZEi6du2a3r17Z+jQobn66quzYMGCpthj1f33359Ro0Zl0KBB6dq1a3r27JnBgwfnuOOOyw9+8IO89957TboeAEBzkaEAAIqRnwAAipOhAIC2oqZSqVQae/KECRNy4oknpr6+fp3PDx48OBMnTsygQYMavcEkeffdd3PyySfn7rvvXu+4adOm5aMf/ehGz1tfX5+ePXtm7ty56dGjxybtcV0GXjyxyedk8/XKVUe09hYAtkjN/X7fGDIUAFB2ZXu/l58AgM1B2d7zZSgAoOyKvN/XNnaRadOm5fjjj8/ChQvTrVu3XHLJJRkxYkQWLlyY22+/Pf/5n/+ZmTNn5ogjjsiUKVPSvXv3Rq0zd+7cHHrooZk6dWqS5Nhjj81xxx2XD33oQ9lqq60ya9asPPjgg/l//+//NfalAAC0GBkKAKAY+QkAoDgZCgBoaxp9B6e/+7u/y8MPP5za2to89NBDGT58+BrPX3311bnooouSJJdddlkuv/zyRm3wpJNOyq233pqOHTvmpz/9aY466qh1jqtUKlm+fHlqaze+s+UOTrQkd3ACaB1l+y+9ZCgAYHNQpvd7+QkA2FyU6T1fhgIANgdF3u/bNWaBJ598Mg8//HCS5NRTT10rFCXJ6NGjs8suuyRJrr322ixdurTwOo888khuvfXWJMk3v/nNBkNRktTU1BQKRQAALU2GAgAoRn4CAChOhgIA2qJGFZzuuuuu6vHJJ5+87onbtctJJ52UJJkzZ04mTZpUeJ1///d/T5L07NkzX/nKV4pvFACgRGQoAIBi5CcAgOJkKACgLWpUwemRRx5JknTt2jV77713g+MOPPDA6vGjjz5aaI0lS5bk7rvvTpIceuih6dSpU5Jk+fLlmTVrVl555ZUsWrSo6NYBAFqNDAUAUIz8BABQnAwFALRFjSo4zZgxI0kyaNCg9d5OcsiQIWuds7F+//vfV4PPHnvskfr6+px33nnZeuuts8MOO2SnnXZKz549c+ihh2by5MnFXwQAQAuToQAAipGfAACKk6EAgLaocMFp0aJFmT17dpKkf//+6x3bq1evdO3aNUkya9asQuv88Y9/rB6vWLEi++yzT6699trMmTOn+vUlS5bk/vvvz0EHHZRx48YVmh8AoCXJUAAAxchPAADFyVAAQFtVuOA0b9686nG3bt02OH5VMHrvvfcKrfPOO+9Uj8eNG5fnn38+n/rUp/Lkk09m0aJFefPNN/ODH/wgPXv2TKVSycUXX1y9FWZDFi9enPr6+jX+AAC0BBkKAKAY+QkAoDgZCgBoqxp1B6dVOnTosMHxHTt2TJIsXLiw0Drz589fY81DDz00v/zlLzN06NB07Ngx22yzTb70pS/ll7/8Zdq1W/kyLrnkklQqlQbnHDt2bHr27Fn9M2DAgEJ7AgBoLBkKAKAY+QkAoDgZCgBoqwoXnDp16lQ9XrJkyQbHL168OEnSuXPnRq+TrGx/b7XVVmuN23///fOZz3wmycrPB3766acbnPOSSy7J3Llzq3+K3m4TAKCxZCgAgGLkJwCA4mQoAKCtKlxw6t69e/V4Y25XuarBvTG3wWxonW222SZ77rlng2M/+clPVo+feuqpBsd17NgxPXr0WOMPAEBLkKEAAIqRnwAAipOhAIC2qlF3cOrTp0+S5LXXXlvv2HfffbcajIreRnL18f3799/osW+99VahdQAAWoIMBQBQjPwEAFCcDAUAtFWFC05JsuuuuyZJXnjhhSxbtqzBcc8++2z1eJdddim0xm677VY9Xr58+XrHrv58bW1toXUAAFqKDAUAUIz8BABQnAwFALRFjSo47b///klW3rZy6tSpDY578MEHq8f77bdfoTV23HHH7LDDDkmSV155JZVKpcGxL774YvV4++23L7QOAEBLkaEAAIqRnwAAipOhAIC2qFEFp2OOOaZ6fNNNN61zzIoVK3LLLbckSerq6jJixIjC63z2s59NktTX1+eBBx5ocNzPf/7z6vGq0AYAUDYyFABAMfITAEBxMhQA0BY1quA0bNiwHHDAAUmSG264IY8//vhaY6655prMmDEjSXLuueemffv2azw/efLk1NTUpKamJqNGjVrnOuedd146deqUJLngggtSX1+/1pj/+q//yuTJk5MkRxxxROHPCAYAaCkyFABAMfITAEBxMhQA0BY1quCUJNdee206d+6cZcuW5bDDDsvYsWPzxBNPZNKkSTnzzDNz0UUXJUkGDx6c0aNHN2qNHXbYIV//+teTJE8//XSGDRuWm266KVOnTs2kSZNyzjnnVENVjx498t3vfrexLwcAoEXIUAAAxchPAADFyVAAQFtT29gT99xzz9xxxx058cQTU19fnzFjxqw1ZvDgwZk4cWK6d+/e6A1eeOGFeeeddzJu3Lg899xzOeWUU9Ya07dv39x1113ZeeedG70OAEBLkKEAAIqRnwAAipOhAIC2ptF3cEqSI488MtOnT8/555+fwYMHp0uXLqmrq8s+++yTcePGZdq0aRk0aNAmb3Ls2LF59NFH88UvfjEDBw5Mx44d07NnzwwdOjTf+MY3MnPmzAwfPnyT1wEAaAkyFABAMfITAEBxMhQA0JbUVCqVSmtvorXU19enZ8+emTt3bnr06NHk8w+8eGKTz8nm65WrjmjtLQBskZr7/X5L5HsKAG2f9/um5fsJAFsG7/lNy/cTANq+Iu/3m3QHJwAAAAAAAAAAgOak4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAadW29gYAAGB1Ay+e2NpboEReueqI1t4CAAAAAADQytzBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEqrtrU3AAAAAAAAANCSBl48sbW3QIm8ctURrb0FADbAHZwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASqu2tTcAAAAAADTewIsntvYWKJFXrjqitbcAAAAATc4dnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEprkwtOr776akaPHp0hQ4aka9eu6d27d4YOHZqrr746CxYsaIo9rmXBggX54Ac/mJqamtTU1GTgwIHNsg4AQHORoQAAipGfAACKk6EAgLaidlNOnjBhQk488cTU19dXv7ZgwYJMmTIlU6ZMyfjx4zNx4sQMGjRokze6uksvvTQvv/xyk84JANBSZCgAgGLkJwCA4mQoAKAtafQdnKZNm5bjjz8+9fX16datW6688so89thjeeCBB3L66acnSWbOnJkjjjgi8+bNa7INT5s2Ld/73vfSqVOndO/evcnmBQBoCTIUAEAx8hMAQHEyFADQ1jS64HTuuedm4cKFqa2tzX333ZcxY8Zk+PDhOeigg3L99dfn29/+dpKV4eiaa65pks0uX748p59+epYvX54xY8akd+/eTTIvAEBLkaEAAIqRnwAAipOhAIC2plEFpyeffDIPP/xwkuTUU0/N8OHD1xozevTo7LLLLkmSa6+9NkuXLt2EbaY6z9SpU/PhD384X/3qVzd5PgCAliRDAQAUIz8BABQnQwEAbVGjCk533XVX9fjkk09e98Tt2uWkk05KksyZMyeTJk1qzFJVr776ai699NIkyXXXXZcOHTps0nwAAC1NhgIAKEZ+AgAoToYCANqiRhWcHnnkkSRJ165ds/feezc47sADD6weP/roo41ZqurLX/5y5s+fny9+8Yv5xCc+sUlzAQC0BhkKAKAY+QkAoDgZCgBoixpVcJoxY0aSZNCgQamtrW1w3JAhQ9Y6pzFuv/323HPPPenVq1eTfQ4wAEBLk6EAAIqRnwAAipOhAIC2qHDBadGiRZk9e3aSpH///usd26tXr3Tt2jVJMmvWrEZsL3n33Xdz3nnnJUmuuuqqbLPNNo2aBwCgNclQAADFyE8AAMXJUABAW9VwbbsB8+bNqx5369Ztg+O7du2a+fPn57333iu6VJLkwgsvzBtvvJHhw4fn9NNPb9QcqyxevDiLFy+uPq6vr9+k+QAANpYMBQBQjPwEAFCcDAUAtFWNuoPTKh06dNjg+I4dOyZJFi5cWHSpPPTQQ7nxxhtTW1ub6667LjU1NYXnWN3YsWPTs2fP6p8BAwZs0nwAABtLhgIAKEZ+AgAoToYCANqqwgWnTp06VY+XLFmywfGrmtadO3cutM7ixYtzxhlnpFKp5Nxzz81HPvKRYhtdh0suuSRz586t/mns7TYBAIqSoQAAipGfAACKk6EAgLaq8EfUde/evXq8MbernD9/fpKNuw3m6q688so899xzGTBgQK644opim2xAx44dq010AICWJEMBABQjPwEAFCdDAQBtVeGCU6dOndKnT5+8/fbbee2119Y79t13360Go6K3kRw3blyS5JBDDsmECRPWOWbV3PPnz8/tt9+eJOnbt28OOuigQmsBADQ3GQoAoBj5CQCgOBkKAGirCheckmTXXXfNww8/nBdeeCHLli1Lbe26p3n22Werx7vsskuhNVbdNvOmm27KTTfdtN6xs2fPzgknnJAkOfDAAwUjAKCUZCgAgGLkJwCA4mQoAKAtateYk/bff/8kKxvXU6dObXDcgw8+WD3eb7/9GrMUAECbIUMBABQjPwEAFCdDAQBtUaMKTsccc0z1uKFW9ooVK3LLLbckSerq6jJixIhCa1QqlQ3+2XHHHZMkO+64Y/VrkydPbsxLAgBodjIUAEAx8hMAQHEyFADQFjWq4DRs2LAccMABSZIbbrghjz/++FpjrrnmmsyYMSNJcu6556Z9+/ZrPD958uTU1NSkpqYmo0aNasw2AAA2KzIUAEAx8hMAQHEyFADQFq37Q3c3wrXXXpv99tsvCxcuzGGHHZYxY8ZkxIgRWbhwYW6//fZcf/31SZLBgwdn9OjRTbZhAIDNmQwFAFCM/AQAUJwMBQC0NY0uOO2555654447cuKJJ6a+vj5jxoxZa8zgwYMzceLEdO/efZM2CQDQVshQAADFyE8AAMXJUABAW9Ooj6hb5cgjj8z06dNz/vnnZ/DgwenSpUvq6uqyzz77ZNy4cZk2bVoGDRrUVHsFAGgTZCgAgGLkJwCA4mQoAKAtqalUKpXW3kRrqa+vT8+ePTN37tz06NGjyecfePHEJp+TzdcrVx3R2lsA2CI19/v9lkiGoiXJUACtQ4ZqWvITLUl+Amg9MlTTkqFoSTIUQOso8n6/SXdwAgAAAAAAAAAAaE4KTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJRWbWtvAAAAAAAAoLkMvHhia2+BEnnlqiNaewsAADSCOzgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGnVtvYGAAAAAAAAAAAoj4EXT2ztLVAir1x1RGtvwR2cAAAAAAAAAACA8lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLRqW3sDAAAAZTbw4omtvQVK5JWrjmjtLbgmWUMZrkkAAAAAaG7u4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFq1rb0BAAAAAAAAANiSDbx4YmtvgRJ55aojWnsLUDru4AQAAAAAAAAAAJSWOzgBAAAAANBk3H2A1bn7AAAA0BTcwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKa5MLTq+++mpGjx6dIUOGpGvXrundu3eGDh2aq6++OgsWLNikuRcsWJCf//znOeusszJ06ND06tUr7du3T58+fTJ8+PBcfvnl+etf/7qpLwEAoMXJUAAAxchPAADFyVAAQFtRuyknT5gwISeeeGLq6+urX1uwYEGmTJmSKVOmZPz48Zk4cWIGDRpUeO7p06dnv/32y3vvvbfWc++8806eeOKJPPHEE/nud7+b66+/Pscff/ymvBQAgBYjQwEAFCM/AQAUJ0MBAG1Jo+/gNG3atBx//PGpr69Pt27dcuWVV+axxx7LAw88kNNPPz1JMnPmzBxxxBGZN29e4fnr6+uroWi//fbL2LFj8+tf/zq//e1vc++99+bMM89Mu3btUl9fny984Qv51a9+1diXAgDQYmQoAIBi5CcAgOJkKACgrWn0HZzOPffcLFy4MLW1tbnvvvsyfPjw6nMHHXRQdt5551x00UWZOXNmrrnmmlx++eWF5m/Xrl3+8R//MZdddll23XXXtZ4/7LDD8vd///c59thjs3z58pxzzjl5/vnnU1NT09iXBADQ7GQoAIBi5CcAgOJkKACgrWnUHZyefPLJPPzww0mSU089dY1QtMro0aOzyy67JEmuvfbaLF26tNAaH//4x3PHHXesMxStcvTRR+czn/lMkuTFF1/MtGnTCq0BANCSZCgAgGLkJwCA4mQoAKAtalTB6a677qoen3zyyeueuF27nHTSSUmSOXPmZNKkSY1ZaoNGjBhRPX7xxRebZQ0AgKYgQwEAFCM/AQAUJ0MBAG1RowpOjzzySJKka9eu2XvvvRscd+CBB1aPH3300cYstUGLFy+uHm+11VbNsgYAQFOQoQAAipGfAACKk6EAgLaoUQWnGTNmJEkGDRqU2traBscNGTJkrXOa2oMPPlg9XnUrTQCAMpKhAACKkZ8AAIqToQCAtqhwwWnRokWZPXt2kqR///7rHdurV6907do1STJr1qxGbG/9fv/732fixIlJkj322EMwAgBKS4YCAChGfgIAKE6GAgDaqoZr2w2YN29e9bhbt24bHN+1a9fMnz8/7733XtGl1mvx4sU57bTTsnz58iTJlVdeuVHnrH4rzPr6+ibdEwBAQ2QoAIBi5CcAgOJkKACgrWrUHZxW6dChwwbHd+zYMUmycOHCokut11e+8pVMmTIlSTJy5MgceeSRGzxn7Nix6dmzZ/XPgAEDmnRPAAANkaEAAIqRnwAAipOhAIC2qnDBqVOnTtXjJUuWbHD8qqZ1586diy7VoLFjx2b8+PFJkqFDh+b73//+Rp13ySWXZO7cudU/zXG7TQCAdZGhAACKkZ8AAIqToQCAtqrwR9R17969erwxt6ucP39+ko27DebG+OEPf5gxY8YkSYYMGZJ77rmn+vnAG9KxY8dqEx0AoCXJUAAAxchPAADFyVAAQFvVqDs49enTJ0ny2muvrXfsu+++Ww1GTXEbydtuuy1f/vKXkyQ77rhjfv3rX2frrbfe5HkBAJqbDAUAUIz8BABQnAwFALRVhQtOSbLrrrsmSV544YUsW7aswXHPPvts9XiXXXZpzFJVv/jFL3LSSSdlxYoV2W677fLAAw+kf//+mzQnAEBLkqEAAIqRnwAAipOhAIC2qFEFp/333z/JyttWTp06tcFxDz74YPV4v/32a8xSSZIHHngg//iP/5hly5alT58++fWvf50PfehDjZ4PAKA1yFAAAMXITwAAxclQAEBb1KiC0zHHHFM9vummm9Y5ZsWKFbnllluSJHV1dRkxYkRjlspjjz2Wo48+OosXL07Pnj1z7733ZrfddmvUXAAArUmGAgAoRn4CAChOhgIA2qJGFZyGDRuWAw44IElyww035PHHH19rzDXXXJMZM2YkSc4999y0b99+jecnT56cmpqa1NTUZNSoUetc53e/+12OOOKIzJ8/P127ds3EiROz9957N2bLAACtToYCAChGfgIAKE6GAgDaotrGnnjttddmv/32y8KFC3PYYYdlzJgxGTFiRBYuXJjbb789119/fZJk8ODBGT16dOH5X3zxxXzyk5/MnDlzkiTf/OY307NnzzzzzDMNntO3b9/07du3Ua8HAKAlyFAAAMXITwAAxclQAEBb0+iC05577pk77rgjJ554Yurr6zNmzJi1xgwePDgTJ05M9+7dC8//8MMP580336w+Pv/88zd4zmWXXZbLL7+88FoAAC1FhgIAKEZ+AgAoToYCANqaRn1E3SpHHnlkpk+fnvPPPz+DBw9Oly5dUldXl3322Sfjxo3LtGnTMmjQoKbaKwBAmyBDAQAUIz8BABQnQwEAbUlNpVKptPYmWkt9fX169uyZuXPnpkePHk0+/8CLJzb5nGy+XrnqiNbeAsAWqbnf77dEMhQtqQwZyjXJ6lyTlE1zXZMyVNOSn2hJ3qsoG9ckZdOc16QM1bRkKFqS9yvKxjVJ2ZThd1CbdAcnAAAAAAAAAACA5qTgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUloITAAAAAAAAAABQWgpOAAAAAAAAAABAaSk4AQAAAAAAAAAApaXgBAAAAAAAAAAAlJaCEwAAAAAAAAAAUFoKTgAAAAAAAAAAQGkpOAEAAAAAAAAAAKWl4AQAAAAAAAAAAJSWghMAAAAAAAAAAFBaCk4AAAAAAAAAAEBpKTgBAAAAAAAAAAClpeAEAAAAAAAAAACUVm1rbwBoOQMvntjaW6BEXrnqiNbegmuSNZThmgQAAAAAAADKxx2cAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtBScAAAAAAAAAACA0lJwAgAAAAAAAAAASkvBCQAAAAAAAAAAKC0FJwAAAAAAAAAAoLQUnAAAAAAAAAAAgNJScAIAAAAAAAAAAEpLwQkAAAAAAAAAACgtBScAAAAAAAAAAKC0FJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8EJAAAAAAAAAAAoLQUnAAAAAAAAAACgtDa54PTqq69m9OjRGTJkSLp27ZrevXtn6NChufrqq7NgwYKm2GOS5Fe/+lWOPfbY9O/fPx07dkz//v1z7LHH5le/+lWTrQEA0BLkJwCA4mQoAIBi5CcAoC2p3ZSTJ0yYkBNPPDH19fXVry1YsCBTpkzJlClTMn78+EycODGDBg1q9BorVqzIGWeckRtuuGGNr7/++ut5/fXXc9ddd+W0007LD3/4w7Rr54ZUAEC5yU8AAMXJUAAAxchPAEBb0+g0MW3atBx//PGpr69Pt27dcuWVV+axxx7LAw88kNNPPz1JMnPmzBxxxBGZN29eozf4f/7P/6kGoz333DO33XZbnnzyydx2223Zc889kyTjx4/Pv/zLvzR6DQCAliA/AQAUJ0MBABQjPwEAbVGj7+B07rnnZuHChamtrc19992X4cOHV5876KCDsvPOO+eiiy7KzJkzc8011+Tyyy8vvMbMmTPzne98J0myzz775KGHHkrnzp2TJEOHDs1RRx2VAw88MFOmTMnVV1+dU045ZZOa5gAAzUl+AgAoToYCAChGfgIA2qJG3cHpySefzMMPP5wkOfXUU9cIRquMHj06u+yyS5Lk2muvzdKlSwuv873vfS/Lli1Lkvzbv/1bNRit0qVLl/zbv/1bkmTZsmX57ne/W3gNAICWID8BABQnQwEAFCM/AQBtVaMKTnfddVf1+OSTT173xO3a5aSTTkqSzJkzJ5MmTSq0RqVSyd13350kGTJkSPbdd991jtt3333z4Q9/OEly9913p1KpFFoHAKAlyE8AAMXJUAAAxchPAEBb1aiC0yOPPJIk6dq1a/bee+8Gxx144IHV40cffbTQGi+//HL+/Oc/rzXP+tZ5/fXX88orrxRaBwCgJchPAADFyVAAAMXITwBAW9WogtOMGTOSJIMGDUptbW2D44YMGbLWORvrj3/84zrnaep1AABagvwEAFCcDAUAUIz8BAC0VYULTosWLcrs2bOTJP3791/v2F69eqVr165JklmzZhVa57XXXqseb2idAQMGVI+LrgMA0NzkJwCA4mQoAIBi5CcAoC1ruLrdgHnz5lWPu3XrtsHxXbt2zfz58/Pee+812zqrAliS9a6zePHiLF68uPp47ty5SZL6+vpCe9tYKxYvaJZ52Tw113VWhGuS1bkmKZvmuiZXzVupVJpl/o2xOeenRIaidXm/omxck5SNDPW+MmUo+YnW5L2KsnFNUjbNeU22dobanPNTIkPRurxfUTauScqmDL+DKlxwWrRoUfW4Q4cOGxzfsWPHJMnChQubbZ1Va2xonbFjx+aKK65Y6+urt8ehufT8XmvvANbkmqRsmvuanDdvXnr27Nm8izRgc85PiQxF6/J+Rdm4JikbGep9ZcpQ8hOtyXsVZeOapGxa4ppsrQy1OeenRIaidXm/omxck5RNGX4HVbjg1KlTp+rxkiVLNjh+VdO6c+fOzbbO6m3u9a1zySWX5IILLqg+XrFiRd5555306dMnNTU1hfbHxqmvr8+AAQMya9as9OjRo7W3A65JSsc12fwqlUrmzZuXfv36tdoeNuf8lMhQrcHPBsrGNUnZuCabnwy1/jU2tI781PL8XKBsXJOUjWuyZbR2htqc81MiQ7UGPxsoG9ckZeOabH5F8lPhglP37t2rxxtzy8r58+cn2bhbYTZ2nVVrbGidjh07rtEUT5K6urpC+6JxevTo4R88peKapGxck82rte46sMrmnJ8SGao1+dlA2bgmKRvXZPOSoRpeY0PryE+tx88FysY1Sdm4Jptfa2aozTk/JTJUa/KzgbJxTVI2rsnmtbH5qV3RiTt16pQ+ffokSV577bX1jn333XerwaXo7SP79+9fPd7QOrNmzaoeu00lAFA28hMAQHEyFABAMfITANCWFS44Jcmuu+6aJHnhhReybNmyBsc9++yz1eNddtmlUWv87TxNvQ4AQEuQnwAAipOhAACKkZ8AgLaqUQWn/fffP8nK20pOnTq1wXEPPvhg9Xi//fYrtMZOO+1U/Yy91edZl4ceeihJsv3222fgwIGF1qF5dezYMZdddtlatxSF1uKapGxck1sO+Yki/GygbFyTlI1rcsshQ7Gx/FygbFyTlI1rcsshP1GEnw2UjWuSsnFNlktNpVKpFD3pySefzMc+9rEkyZlnnpnrrrturTErVqzI7rvvnhkzZqSuri5vvvlm2rdvX2idL3/5y/nBD36QJHn88cez7777rjXmiSeeyPDhw6vjv//97xd9OQAAzU5+AgAoToYCAChGfgIA2qpG3cFp2LBhOeCAA5IkN9xwQx5//PG1xlxzzTWZMWNGkuTcc89dKxhNnjw5NTU1qampyahRo9a5znnnnZetttoqSXLOOedk4cKFazy/cOHCnHPOOUmS2tranHfeeY15OQAAzU5+AgAoToYCAChGfgIA2qpGFZyS5Nprr03nzp2zbNmyHHbYYRk7dmyeeOKJTJo0KWeeeWYuuuiiJMngwYMzevToRq0xePDgXHjhhUmSKVOmZL/99ssdd9yRKVOm5I477sh+++2XKVOmJEkuvPDC7Lzzzo19OQAAzU5+AgAoToYCAChGfgIA2qJGfUTdKhMmTMiJJ56Y+vr6dT4/ePDgTJw4MYMGDVrrucmTJ2fEiBFJkpEjR+bmm29e5xwrVqzI6aefnhtvvLHBfZx66qm5/vrr065do/taAAAtQn4CAChOhgIAKEZ+AgDamk1KE0ceeWSmT5+e888/P4MHD06XLl1SV1eXffbZJ+PGjcu0adPWGYwKbbBdu9xwww2ZOHFijj766PTr1y8dOnRIv379cvTRR+eee+7J+PHjBaM26vLLL6/eBpUty6hRo1JTU5OBAwe29lYAmpT8REuQobZcMhTQVslQNDf5acslPwFtlfxES5ChtlwyFNAaajd1gh133DH/+q//mn/9138tdN4nPvGJFLl51OGHH57DDz+86PYAAEpHfgIAKE6GAgAoRn4CANoSlWkAAAAAAAAAAKC0FJyAUrr55ptTqVTyyiuvtPZWAAA2GzIUAEAx8hMAQHEyFNAaFJwAAAAAAAAAAIDSUnACAAAAAAAAAABKS8FpC3X55ZenpqYmNTU1SZK5c+fmG9/4Rvbcc8/U1dWlpqYmN998c3X8e++9l6uuuirDhw9P796907Fjx/Tv3z/HHXdcfvnLX27UmitWrMhtt92Wz372s9lhhx3SuXPndO7cOYMHD84XvvCF/OxnP8vSpUsLv5bf/e532XbbbVNTU5Ptttsu06dPLzwHG/bnP/85F198cfbaa6/07Nkz7du3z7bbbps99tgjJ5xwQm6++ebU19ev89x58+blmmuuyUEHHZQPfOAD6dChQ3r06JE999wz55xzTh599NG1zhk1alRqamoycODABvd055135phjjkn//v3TsWPHdO/ePR/84AdzwAEH5Gtf+1qefPLJJn8tSfL000/njDPOyM4775wuXbqke/fu2W233XL++eev91acr7zySvXf3ap/Xz//+c9z+OGHp1+/fqmtrc0nPvGJNc6ZOXNmzjnnnOy+++7p3r17OnTokH79+uWjH/1oTjnllNxxxx1ZvHhxg2uypr/92Tdnzpxcdtll2W233dKtW7f07t07I0aMyG233dbgHEuWLMmECRPyla98JUOHDk2vXr3Svn379OnTJx/72Mdy+eWXZ/bs2evdx8CBA1NTU5NRo0YlSZ566qmccMIJGTBgQDp16pQBAwbk5JNPzrPPPrtRr+uFF17I+eefnz322CM9e/ZM586d88EPfjCjRo3KlClTGjxv8uTJ1e/H5MmTs2LFitx4440ZMWJEtt1227Rr1666R2AlGYqiZKj3yVCbLxlqTTIUFCM/UZT89D75afMlP61JfoLiZCiKkqHeJ0NtvmSoNclQ61Bhi3TZZZdVklSSVGbOnFkZOHBg9fGqPzfddFOlUqlUfvvb31b69eu31vOr//nMZz5TWbhwYYPrvfzyy5WPfvSj650jSWXSpEkN7nNdHnrooUrPnj0rSSoDBw6svPDCC031LWI1Dz30UKVHjx4b/PubMGHCWuf++te/rmy99dYbPPdvjRw5spKksuOOO6713LJlyyr/8A//sME599577yZ9LZVKpfKtb32r0q5duwbP69ixY+VHP/rROs99+eWXq+NuvPHGyhe/+MW1zj/wwAOr43/6059WOnTosMG9Pv300w38zfG3Vv+Z8tJLL1U+9KEPNfh9/cd//MfK0qVL15pj1bW5vj99+vSpPPLIIw3uY8cdd6wkqYwcObJyww03VGpraxu8nn7605+u9zVdffXVlfbt2ze4l5qamsrXvva1dZ47adKk6rhf/epXlUMOOWSt80eOHFnoewxtnQxFETLU+2SozZsMtSYZCoqRnyhCfnqf/LR5k5/WJD9BcTIURchQ75OhNm8y1JpkqLXVhi3ecccdl9dffz3nnHNOjjrqqPTq1SvPP/98dtxxx7z++us5+OCD8+6771Zbip/73OfSp0+f/PGPf8w111yT3//+9/n5z3+eUaNG5fbbb19r/jfeeCP77bdf/vznPydJDjrooIwcOTJDhgxJTU1NXn755fzmN7/Jf//3fxfa98SJE/MP//APWbhwYXbbbbfcd9996devX5N8T3jf4sWL87nPfS719fXp3r17zjrrrIwYMSJ9+/bNkiVL8vLLL+exxx7LnXfeuda5kyZNyt///d9n2bJl2WqrrfLFL34xRx99dHbYYYcsWrQof/zjH/OrX/0qEyZMKLSnH/zgB9XrZf/9989pp52WD33oQ+natWvefvvtTJ8+Pf/zP/+TuXPnNtlrSZL/+I//yJgxY5Ik22yzTb761a9mv/32y/Lly3P//ffn6quvzvz58zNq1KhsvfXWOfzwwxt8Dd/73vcyffr0HHDAATnrrLMyePDgzJkzp9ocf+ONN3LyySdnyZIl6du3b77yla9k3333zdZbb52FCxfmhRdeyIMPPpi77rqr0PeO9x1//PF5+eWX86UvfSnHHXdcevbsmenTp2fcuHGZOXNmfvrTn6Zfv3757ne/u8Z5y5Ytywc/+MEce+yxGTZsWHbYYYfU1tbm1Vdfzf33358bb7wxb7/9do499tg888wz6du3b4N7+N3vfpef/OQn6du3by655JIMGzYsixYtyj333JPvfe97Wbx4cb7whS9kp512yj777LPW+VdffXUuuuiiJMlHPvKRnHXWWdl5551TV1eX5557Lv/+7/+exx9/PN/4xjey9dZb55/+6Z8a3MtXv/rVTJ8+PUcddVRGjRqVHXfcMW+88cZ6/ysI2NLJUKyPDPU+GaptkaHWJENBMfIT6yM/vU9+alvkpzXJT1CcDMX6yFDvk6HaFhlqTTLU/19rN6xoHau3H9u1a1e599571znuuOOOq44bP378Ws8vWrSoMmLEiOqYe+65Z60xxx57bPX5cePGNbinefPmVd55550G97m6H//4x9Wm5LBhwypvv/32xrxsGuGBBx7YYBu6UqlUli5dWpk7d2718cKFC6v/xUCXLl3WavWv7k9/+tNaX1tf8/uAAw6oJKl87GMfW2czd5W/vS4a+1oqlUrlzTffrHTp0qWSpNKvX7917vm3v/1tpWvXrpUkle23376yZMmSNZ5fvfmdpHLSSSdVVqxYsc493HDDDRvV7F6wYEFlwYIFDT7Pmlb/mZKk8pOf/GStMfX19ZX/7//7/6o/H//2+//CCy80+PdWqVQq06dPr3Tr1q2SpPIv//Iv6xyzqvm96hr/y1/+staY3/zmN9Wfc0OHDl3r+T/84Q/Vxvdll122zj0tX768cuKJJ1aSVLp167bWz9jVm9/r2y/wPhmKjSVDrSRDtQ0ylAwFm0J+YmPJTyvJT22D/CQ/waaSodhYMtRKMlTbIEPJUBui4LSFWv2HwymnnLLOMa+//nplq622qiSpfOpTn2pwrpdffrn6j/fwww9f47lnn322UlNTU0lSOeaYYzZpn6t8//vfr8558MEHV+bNm1d4Xjbej3/84+rfwd+GhfX54Q9/WD3ve9/7XuF11xeMdt5550qSyvnnn19ozsa+lkqlUhk3blz13Ntvv73Bcd/85jer4/72loSrB6O6urpKfX19g/NceeWVlSSVXr16Fdon67f6z5RPf/rTDY773//93+q4s88+u/A65513XiVJZffdd1/n86sHo5/97GcNznPWWWdVxz311FNrPHfKKadUklT22Wef9Qa1d999t9KxY8dKksr111+/xnOrB6PBgwdXli1bVuBVwpZJhmJjyVAryVBtgwwlQ8GmkJ/YWPLTSvJT2yA/yU+wqWQoNpYMtZIM1TbIUDLUhrQLW7wvfOEL6/z65MmTs3z58iTJqaee2uD5AwcOzKGHHrrWOcnK209WKpUkyfnnn7/Je/3mN7+Zs88+O5VKJccee2wmTpyYbt26bfK8NGy77barHt90000bfd4vf/nLJEnXrl1z+umnN8ueJkyYkNmzZxc+Lyn2WpLk/vvvT5LU1dXlM5/5TIPjTjvttLXOWZcjjzwy3bt3b/D5VXt99913c/fddxfaKxvn5JNPbvC5YcOGZbfddkuy/r/HZOXf0Ysvvpg//OEPeeaZZ/LMM8+krq4uSfLHP/4xS5cubfDcXr165eijj27w+VNOOaV6/Lf7WHVL2M9+9rOpqalpcI66urrsscceSZLHH3+8wXHHH398ttpqqwafB9YmQ7E+MtRKMlTbI0OtSYaCYuQn1kd+Wkl+anvkpzXJT1CcDMX6yFAryVBtjwy1JhlqJQUn8pGPfGSdX3/mmWeqxx/72MfWO8eq5xcsWJCXXnqp+vVp06YlSdq3b5999913k/Z5/vnn52tf+1qSlT/Q/vu//zsdO3bcpDnZsP333z8f/OAHkyTnnXdehg0blrFjx+bRRx/NkiVLGjxv1d/93nvvnS5dujTpnkaOHJkkeeGFFzJo0KCccsopue222/Laa6+t97zGvpbk/X8Pe+21V9q3b9/guG233TYDBw5c45x1aejf3SpHHXVU9c312GOPzUEHHZTvfve7mTp16hr/44PGGzp06HqfHzZsWJJk5syZa10fTz/9dE455ZRst9126d27dwYNGpTdd989e+yxR/bYY49cfvnlSZIVK1bk3XffbXCNPffcM7W1tQ0+/9GPfjQdOnSorrnKq6++mrfeeitJcskll6Smpma9f6ZMmZIk+etf/9rgWhu6JoG1yVCsjwy1kgzV9shQa5KhoBj5ifWRn1aSn9oe+WlN8hMUJ0OxPjLUSjJU2yNDrUmGWknBifTq1WudX3/nnXeqx3379l3vHB/4wAfWed6qVm7v3r2r/7gb63vf+16SZPfdd8/48eM1FFtI+/btM2HChOyyyy5JkqeeeipjxozJ/vvvn7q6unzqU5/KT37yk7XerFf93a/etm4qp5xySsaMGZPa2trMnTs3N910Uz7/+c9nwIABGTRoUEaPHr1GQN/U15K8f11v6N9C8v6/h9X/Lfythv7drdKnT5/84he/yPbbb59KpZJJkyblggsuyD777JPevXvnM5/5TLVdT+Ns6O9y2223TZJUKpU1ws0NN9yQvfbaKzfddNN6g8YqCxcubPQeamtr07t37yRrXk9vvvnmBtddlwULFjT43IauSWBtMhTrI0OtJEO1PTLUmmQoKEZ+Yn3kp5Xkp7ZHflqT/ATFyVCsjwy1kgzV9shQa5KhVlJwYqMCxvpum9ZSPvvZzyZZ2aY999xzW3k3W5Zdd901Tz/9dO68886ccsopGTRoUJKVP/DvvffefOELX8jHPvaxRv+wbowrr7wyL7zwQq688socdNBB1Xb5iy++mH/913/NkCFDct111zX5a2mqfwsb8+/ugAMOyAsvvJD/+q//yuc///n0798/SVJfX58777wzRx55ZD71qU+t982OhjXm7/LZZ5/Nl770pSxbtix9+/bN1VdfnalTp+btt9/OkiVLUqlUUqlUcsMNN1TPWXV736baQ5I1wvull16ap59+eqP+rO+Wrv7HJhQnQ7EhMtT7ZKi2Q4ZakwwFxchPbIj89D75qe2Qn9YkP0FxMhQbIkO9T4ZqO2SoNclQKyk40aBVbcMkeeONN9Y7dvX24+rnbb311klWNhY3dOvADbnttttyzDHHJEn+/d//vUk+C5iNt9VWW+WYY47JDTfckOeffz5//vOfc+ONN2bvvf9/7d19aFV1HMfx73Z3b9sUNzWQbdK2wtlkYE0rUWiiPYAmwZT9EziMohyNRSWUQmuEShJMKCOM5sxlmBlo2paWumxbu3vAldIfkaKLHAZrYHPZbnz6Y+xs1/uwc53Tzb1fcOByfw/n97v3nLPP4HfPWWBmZm1tbfbCCy849Qe/+0uXLo3ZmDIzM23jxo323XffWU9PjzU0NFhZWZklJiZaf3+/lZSUOLfYHM1czIaO65HOBbOh82H4uXCjEhMT7ZlnnrFPP/3UOjs77dy5c/bee+9ZTk6OmZl98803tmnTplHvZzIa6bscLI+Li3NWRVdXV1sgEDCPx2P19fX22muvWX5+vs2YMSPolqfRVv3HMoZAIOD0Nfx4mjlzpvPa6/VaXl6eqy07O9vVuACMDhkKw5GhyFB3GjIUgLFAfsJw5Cfy052G/ARgrJChMBwZigx1pyFDIRwWOCGivLw853Vzc3PUun6/38zMkpOTnWejmg0859TMrL+/35qamkY1Hq/Xa/v27bOnnnrKzAZudblhw4ZR9Ykbl5aWZuvWrbOmpibnez58+LBzG7/B91pbW2/JymSv12uLFy+27du32969e81sYMXtF198MWLbkeZiNnQ+tLe3WyAQiNjX5cuX7cKFC0Ftbqbs7Gx76aWXrKWlxVkJ/vnnn9/0/UwGLS0trsrnzJnj3Jr37NmzZmY2f/58u//++yO2HXxW7khOnz4d9Xjq6Ohw/qkcfjzde++9lpKSYmZmDQ0NrvYF4NYhQyEaMlR4ZKiJgwwFYCyQnxAN+Sk88tPEQX4CMFbIUIiGDBUeGWriIEMhHBY4IaKlS5c6tzqrqqqKWO/ixYt27NixkDZmZitXrnRu3Tb47N3R8Pl8duDAAVuxYoWZmb377rv2+uuvj7pf3Div12sFBQVmNrBKtaenx8zMVq1aZWYDzwrduXPnLR3T8uXLndeDzxB2I9JczMwee+wxMzPr6emxL7/8MmIfH3/8sXMrw8E2Y2HatGn20EMPmVlsc8SQ3bt3RyxraWmxM2fOmFnw9zgYYnp7eyO2vXTpkh06dMjVGLq7u+2rr76KWD782jt8HB6Px7kOHj161H755RdX+wNwa5Ch4AYZKhgZauIgQwEYC+QnuEF+CkZ+mjjITwDGChkKbpChgpGhJg4yFMISJqXy8nKZmUY6BNasWePUq66uDim/du2ali9f7tT5+uuvQ+oUFhY65du2bYu4r7///lvd3d2uxvnPP//oySefdMo2bdoUdR64cd9//71+/fXXiOXXrl1Tfn6+zExTp05Vf3+/JKmvr08ZGRkyMyUnJ+vkyZMR++js7Ax5r7i4WGamzMzMkLI9e/Y4+wln//79zrGxdevWUc9Fki5fvqzk5GSZmWbPnq3ff/89pP3p06c1depUmZkyMjL077//BpWfP3/eGdeuXbsijkOS6urq9Mcff0Qs7+npcT7fuXPnRu0LQ4ZfU8xM+/btC6lz5coVPfjggzIzxcfHq6OjwykrLS113m9oaAhp29vbq6VLlwbt4/z58yH1MjMznfKsrCx1dXWF1Dl58qQSEhJkZlqwYEFI+U8//SSPxyMz07x588KeR4MCgYBqampC6pw4ccIZx4kTJyK2BzCEDAW3yFADyFB3BjIUGQoYDfIT3CI/DSA/3RnIT+QnYLTIUHCLDDWADHVnIEORoUbCAqdJym0w6uzs1PTp050LwXPPPadjx46ptbVVNTU1euCBB5x+ioqKwvbR1dWl9PR0p96yZcv0ySefyO/3q6WlRfv371dJSYlmzJgRcmJGG2dfX58ef/xxp/zNN9+84c8DkZWXlys+Pl4FBQXatm2b6urq1NbWph9++EFVVVV6+OGHne+grKwsqO3x48edC3tCQoLWrVungwcPqq2tTY2NjaqqqtKaNWvk8/lC9hstGJmZZs2apfXr12vPnj1qbGxUe3u7amtr9corrygpKckJNxcvXrwpc5GkHTt2OOWzZs1SZWWlmpub1dDQoIqKCicUxcXF6ciRIyHtYwlGxcXF8nq9WrFihbZv365vv/1W7e3tqq+v144dO5Sbm+v0VVlZGbUvDBl+TVm4cKE8Ho9KSkp0/Phxtba2qqqqSnPnznXqlJaWBrX3+/1OWWpqqjZv3qz6+no1Nzfrgw8+0Jw5c2RmWrJkiatgNH/+fHm9XmVkZOj999+X3+/XqVOn9MYbbygxMdE5d3788cew86msrHT2k5KSog0bNqi2tlbt7e1qbGzU3r17VVpaqrS0NJmZfv7556D2BCMgdmQouEWGGkKGmvjIUGQoYDTIT3CL/DSE/DTxkZ/IT8BokaHgFhlqCBlq4iNDkaFGwgKnScptMJKk9vb2oGATbissLFRfX1/EPn777Tfl5eVF7SPciTnSOK9evaply5Y5dd5+++2YPgeM7PqVspG2p59+WlevXg1pX1dX54TraNv1RgpGI20pKSmqra29qXORpM2bNys+Pj5i27vuuku7d+8O2zbWYORmrC+++KL++++/qH1hyPBj4Ny5c8rOzo742a5evTrsLwwqKiqifievvvqqdu3a5SoYFRcX66OPPnL+gbh+8/l8+uyzz6LOaefOnc6vEqJtPp8v5JcPBCMgdmQouEWGCkaGmtjIUGQoYDTIT3CL/BSM/DSxkZ/IT8BokaHgFhkqGBlqYiNDkaFGwgKnSSqWYCQN3Opt69ateuSRR5Samiqfz6f09HQVFhbq0KFDrvoIBAKqrq7WypUrlZaWJq/Xq6SkJOXk5Gjt2rU6ePCgAoFAzOPs7e1VQUGBU2/Lli2uxgN3rly5ogMHDmj9+vVatGiR7rnnHiUmJioxMVFZWVkqKirS4cOHo/bx119/acuWLVq8eLFmzpwpj8ejadOmKT8/Xy+//LL8fn9Im2jB6MyZM3rnnXe0atUqzZs3z+kzNTVVixYtUnl5edhbBd6MuUhSR0eHnn/+ed13331KSkrSlClTlJubq7KysrB/BAfFEoy6u7tVU1OjZ599VgsXLlRGRoZ8Pp9zzhQXF+vUqVMjjhXBrr+mdHd3a+PGjcrNzVVycrJSUlL06KOPqqamJmo/R44c0RNPPKHp06fL5/Np9uzZKiws1NGjRyUppmAkSU1NTSoqKlJ6erp8Pp8yMjK0du1anT171tW8urq6VFFRoSVLlujuu+9WQkKCpkyZopycHK1evVoffvih/vzzz5B2BCMgdmQouEWGCkWGmrjIUMHIUEBsyE9wi/wUivw0cZGfgpGfgNiRoeAWGSoUGWriIkMFI0OFipMkAwBgEnjrrbesoqLCzMxu55+/rKwsu3DhghUXF1t1dfVtGwcAAIAbZCgAAIDYkJ8AAABiR4bCSOJv9wAAAAAAAAAAAAAAAAAAIBIWOAEAAAAAAAAAAAAAAAAYt1jgBAAAAAAAAAAAAAAAAGDcYoETAAAAAAAAAAAAAAAAgHGLBU4AAAAAAAAAAAAAAAAAxq04SbrdgwAAAAAAAAAAAAAAAACAcLiDEwAAAAAAAAAAAAAAAIBxiwVOAAAAAAAAAAAAAAAAAMYtFjgBAAAAAAAAAAAAAAAAGLdY4AQAAAAAAAAAAAAAAABg3GKBEwAAAAAAAAAAAAAAAIBxiwVOAAAAAAAAAAAAAAAAAMYtFjgBAAAAAAAAAAAAAAAAGLdY4AQAAAAAAAAAAAAAAABg3GKBEwAAAAAAAAAAAAAAAIBx638JMM7UYzjluwAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 2400x600 with 4 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# 具体例(じゃんけん)\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "transition_matrix = np.zeros([3,3])\n",
    "\n",
    "# 0: rock, 1:scissors, 2:paper\n",
    "transition_matrix[0,0] = 0.4\n",
    "transition_matrix[0,1] = 0.3\n",
    "transition_matrix[0,2] = 0.3\n",
    "transition_matrix[1,0] = 0.6\n",
    "transition_matrix[1,1] = 0.1\n",
    "transition_matrix[1,2] = 0.3\n",
    "transition_matrix[2,0] = 0.5\n",
    "transition_matrix[2,1] = 0.2\n",
    "transition_matrix[2,2] = 0.3\n",
    "\n",
    "T = 1000 # 試行回数\n",
    "np.random.seed(1)\n",
    "\n",
    "# rockを初期状態として計算\n",
    "h_init = 0\n",
    "h_arr = [h_init]\n",
    "h = h_init\n",
    "for i in range(T):\n",
    "    h = np.random.choice(3,1,p=transition_matrix[h,:])[0]\n",
    "    h_arr.append(h)\n",
    "\n",
    "# グラフに描画\n",
    "hist1, bin1 = np.histogram(h_arr[:10],bins=np.arange(4),density=True)\n",
    "hist2, bin2 = np.histogram(h_arr[:100],bins=np.arange(4),density=True)\n",
    "hist3, bin3 = np.histogram(h_arr[:500],bins=np.arange(4),density=True)\n",
    "hist4, bin4 = np.histogram(h_arr[:T],bins=np.arange(4),density=True)\n",
    "\n",
    "plt.rcParams[\"font.size\"] = 20\n",
    "hands_label = ['rock','scissors','paper']\n",
    "\n",
    "plt.figure(figsize=(24,6))\n",
    "plt.subplot(1,4,1)\n",
    "plt.bar(hands_label,hist1)\n",
    "plt.ylim([0, 1])\n",
    "plt.title('10')\n",
    "\n",
    "plt.subplot(1,4,2)\n",
    "plt.bar(hands_label,hist2)\n",
    "plt.ylim([0, 1])\n",
    "plt.title('100')\n",
    "\n",
    "plt.subplot(1,4,3)\n",
    "plt.bar(hands_label,hist3)\n",
    "plt.ylim([0, 1])\n",
    "plt.title('500')\n",
    "\n",
    "plt.subplot(1,4,4)\n",
    "plt.bar(hands_label,hist4)\n",
    "plt.ylim([0, 1])\n",
    "plt.title('1000')\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "回数を重ねるごとに分布の形の変化が小さくなっている=収束していることがわかる。  \n",
    "十分に時間が経過したとき、とりあえずグーを出し続ければ勝つ確率が高いと言える。"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "定常分布は、固有値1の固有ベクトルの和が1となるように規格化するとを求められる。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "定常分布π=:[0.475 0.225 0.3  ]\n"
     ]
    }
   ],
   "source": [
    "import numpy.linalg as LA\n",
    "\n",
    "w, v = LA.eig(transition_matrix.T)\n",
    "idx = np.argmin(np.abs(np.real(w)-1))\n",
    "w = np.real(v[:, idx]).T \n",
    "w /= np.sum(w)\n",
    "print(f\"定常分布π=:{w}\")"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 有限マルコフ連鎖のパラメータ推定\n",
    "\n",
    "推移行列$Q$が未知パラメータ$\\theta$に依存しているとする:\n",
    "$$\n",
    "Q_\\theta = (p_\\theta(i,j))_{(1 \\leq i,j \\leq N)}\n",
    "$$\n",
    "簡単のために、$\\theta$は1次元のパラメータとする。  \n",
    "$X$に関する観測データ$x_0,x_1,...,x_n$が与えられたとき、その実現確率は\n",
    "$$\n",
    "P(X_0=x_0, X_1=x_1,...,X_n=x_n) = p_0(x_0)\\prod_{j=1}^{n}p_\\theta(x_{j-1}, x_j) \n",
    "$$\n",
    "と書ける。ただし、$p_0$はパラメータによらないとする。  \n",
    "このとき対数尤度関数は、\n",
    "$$\n",
    "l_n(\\theta) = \\sum_{j=1}^{n} \\log p_{\\theta}(x_{j-1}, x_j)\n",
    "$$\n",
    "で与えられる。従って、$\\theta$の最尤推定値$\\hat{\\theta}$は以下の尤度方程式の解として求められる:\n",
    "$$\n",
    "\\frac{\\partial}{\\partial \\theta} l_n(\\hat{\\theta})=0\n",
    "$$\n",
    "独立で同一分布に従うようなマルコフ連鎖が複数観測される場合も同様に上記で求められる。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "最尤推定値(Θ):0.07\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkEAAAHFCAYAAAD1zS3+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABNOUlEQVR4nO3deXhTVf4/8PdN0qRL2nSHFkpZSgFlGQYKKAgUAQcVccVxVJBR9Osyo/4cUXRE3EYYZUZHxwGXkUFHURFR3GAEhCqggOxroQstpXTf26RJ7u+P9KYtdEnam9wk9/16nj6Tltt7P3WifXPO55wjiKIogoiIiEhlNEoXQERERKQEhiAiIiJSJYYgIiIiUiWGICIiIlIlhiAiIiJSJYYgIiIiUiWGICIiIlIlhiAiIiJSJYYgIiIiUiWGICIiIlIlhiAXFRUV4csvv8SiRYswY8YMxMbGQhAECIKAO+64Q+nyZJOXl4dPP/0Ujz/+OKZMmQKTyeT8ORcvXuzSPXbv3o1ly5bht7/9LYYPH46EhAQYDAaEh4dj0KBBmDt3LrZs2eLZH4SIiKgTOqUL8Bc9evRQugSPy83NRd++fbt9n4ceegg//vjjBV+3WCw4ceIETpw4gVWrVuGmm27CqlWrEBwc3O1nEhERuYshqAv69OmDwYMHY+PGjUqXIquWZ+kKgoABAwYgMTER27Ztc+s+BoMBkyZNwqWXXoohQ4YgISEB0dHRKC4uxv79+7F8+XJkZ2fjk08+gUajwerVq+X+UYiIiDrFEOSiRYsWIS0tDWlpaejRowdycnLQr18/pcuSVXh4OJ5//nmMGTMGo0ePRlRUFL7//nukp6e7dZ8NGzZAp2v7rXXFFVfgD3/4A6ZMmYKdO3fio48+whNPPIHhw4fL8SMQERG5jCHIRc8884zSJXhcTEwMnnzyyW7fp70AJAkJCcGDDz6InTt3AgAyMjIYgoiIyOvYGK2AhoYGvP7667j88svRs2dP6PV6xMfHY+rUqXjnnXdgtVqVLtHjwsPDna8bGhoUrISIiNSKI0Fetn//fsyaNQu5ubmtvl5cXIxNmzZh06ZNWLFiBdavXx/Qzdgt+4AGDx6sYCVERKRWDEFedPLkSUyaNAmVlZWIiIjA/fffjzFjxiApKQmlpaX44osvsGLFCuzatQuzZs1CRkYGgoKClC5bFna7HcXFxTh8+DBee+01rFu3DoAjAF1xxRXKFkdERKrEEORFc+fORWVlJUaOHImNGzciNja21Z9Pnz4dV199Na666ir89NNPWLlyJebPn69QtfLo27fvBaNekv79+2Pt2rWd9hARERF5AnuCvCQjIwPbt28HAPznP/+5IABJfvOb3+DGG28EAKxcudJb5XmVTqfD888/j3379mHIkCFKl0NERCrFEOQlX3zxBQBg0KBBGDZsWIfXTpw4EQCwa9cuv2+S3rhxIw4ePIj9+/dj8+bNeO655xAfH49nn30W999/P2pqapQukYiIVIrzEF6ye/duAMDx48chCIJL39PY2IiysjLEx8c7v7Zy5UrMmzev2/W03BjRk1JTU1t9np6ejvvvvx9XXHEF3nvvPezfvx8//vgjjEajV+ohIiKScCTIS4qKirr0fXV1dTJXoryoqCj85z//AQAcOHAAf/nLXxSuiIiI1IgjQV5is9kAACNGjMD777/v8vf16tWr1efXXXcdxo0bJ2ttShgyZAgGDhyIzMxMrFmzhkGIiIi8jiHIS2JiYgAANTU1GDp0aJfvYzKZYDKZ5CpLUXFxccjMzGx39RgREZEncTrMS0aOHAkAyMrKQmFhocLV+IYzZ84AAPuBiIhIEQxBXnLNNdcAcDQkv/rqqwpXo7xdu3Y5R4A6Wy1HRETkCQxBXjJ9+nSMGTMGAPDSSy/h448/7vD6gwcPYv369d4oTVY///wzfvnllw6vOXPmDObOnev8fM6cOZ4ui4iI6ALsCXLRDz/8gJMnTzo/Lykpcb4+efLkBRsb3nHHHRfc44MPPsCYMWNQVlaGm2++Ge+//z5uvvlmDBw4EFqtFkVFRdi7dy/Wr1+PnTt34pFHHsHMmTM99SO16dtvv201XXfs2DHn63379rX6OY1Go3NjR8mRI0cwb948XHrppZg5cyZ+9atfIS4uDoAj/GzZsgXvvvsuKisrAQBTp05t858VERGRpwmitzaM8XN33HGHc1m3K9r7x3rixAnccMMNOHToUKf3eOaZZ7Bo0SKXnymHyZMnY+vWrS5dm5ycjJycnFZfc2cfozvuuAP//Oc/ERoa6m6ZRERE3caRIC9LTU3Fvn378PHHH+PTTz/Frl27UFxcDJvNhpiYGAwaNAgTJkzAddddh1//+tdKl+u2m2++GVFRUdi8eTN++eUXFBQU4Ny5c2hsbITJZEJKSgrGjx+P22+/HcOHD1e6XCIiUjGOBBEREZEqsTGaiIiIVIkhiIiIiFSJPUHtsNvtKCgoQHh4uMsHnhIREZGyRFFEdXU1EhMTodF0PNbDENSOgoICJCUlKV0GERERdUFeXh569+7d4TUMQe0IDw8H4PiHGBERoXA1RERE5IqqqiokJSU5f493hCGoHdIUWEREBEMQERGRn3GllYWN0URERKRKDEFERESkSgxBREREpEoMQURERKRKDEFERESkSgxBREREpEoMQURERKRKDEFERESkSgxBREREpEoMQURERKRKDEFERESkSgxBREREpEoMQeRVRwqq8M4P2WhotCldChERqRxPkSev+fJAAR75eD/MVjtySmrx3LVDlS6JiIhUjCNB5HGiKOKfW07igQ/2wmy1AwDe/ykXe3LLFK6MiIjUjCGIPMpitePRNQfw0objAIB54/vihl/3higCj396EJamUERERORtnA4jj7FY7Zj775+xI6sUGgFYfM3FmHNJX1TUWfD98SJkFtVgxdZT+MPlA5UulYiIVIgjQeQxGw4XYkdWKcL0WrxzRxrmXNIXABAZqseimRcBAF7bfBKnimsUrJKIiNSKIYg8ZuuJYgDAreOSkT4ovtWfXTMiEZNS42Cx2bFw7UHY7aISJRIRkYoxBJFHiKLoDEGTUuMu+HNBEPD8tUMREqTFz9ll+GRPnrdLJCIilWMIIo84crYKxdVmhOq1GN03qs1rkqJD8cj0VADASxuOw8bRICIi8iKGIPIIaRTo0gExMOi07V4399K+CDfoUFJjwaEzld4qj4iIiCGIPGPr8fanwloK0mpwaUoMAGBbU3AiIiLyBoYgkl11QyP25JYDACalxndyNTCxKShty2QIIiIi72EIItn9eLIUVruI/rFh6BMT2un1Ewc6QtAvpytQ1dDo6fKIiIgAMASRB0j9QBM7mQqTJEWHon9cGGx2EdtPlniyNCIiIieGIJKVKIrYerwIADB5kGshCGjuHdrKviAiIvIShiCS1cmiGhRUNsCg02Bc/xiXv8/ZF3SiBKLIpfJEROR5DEEkK2kkZ2z/GAQHtb80/nzj+sVAr9PgTEU9j9EgIiKvYAgiWX3v4tL484XotRjbLxoAsPUE+4KIiMjzGIJINnUWK37OLgPgXj+QhH1BRETkTQxBJJudWaWw2OzoHRWC/rFhbn+/1Bf0U1YpGhptcpdHRETUCkMQyablVJggCG5//8B4IxJMwTBb7c4RJSIiIk9hCCLZ7MwqBeD6/kDnEwTBuXEip8SIiMjTGIJIFlabHdkltQCAixMjunyf5qXyDEFERORZDEEki9yyOjTaRIQEaZFoCunyfSakxEIjAJlFNSioqJexQiIiotYYgkgWJ4sce/sMiA+DRuN+P5DEFBqEXyVFAuBoEBEReRZDEMlCCkED48O7fa9LBjh2mt57uqLb9yIiImoPQxDJQgpBKfHGbt9rWC8TAOBQQWW370VERNQehiCShXM6LK77IWhoUwg6ca4aZiv3CyIiIs9gCKJus9tF53lfcowE9YoMQWRoEBptIk4U8hwxIiLyDIYg6razVQ2os9ig0whIjgnt9v0EQcDQRE6JERGRZzEEUbdJU2F9Y8MQpJXnLSVNiR06wxBERESewRBE3eZsipahH0gytJdjw0WGICIi8hSGIOo2OVeGSaTpsKOF1Wi02WW7LxERkYQhiLrtlAdCUHJMKMKDdbBY7c6QRUREJCeGIOq2kzKuDJMIguA8g+wgp8SIiMgDGIKoW0przCirtQAA+seFyXpvaUrsMEMQERF5AEMQdYs0VdUrMgShep2s9x7WW1omXyXrfYmIiACGIOomT0yFSS5uGgk6UlAFm12U/f5ERKRuDEHULc0Hp8ofgvrFhiFUr0V9ow1ZxWyOJiIieTEEUbd4Ynm8RKtpbo7mztFERCQ3hiDqFk8sj29JmhI7mM++ICIikhdDEHVZrdmKgsoGAJ4LQc7jMzgSREREMmMIoi6TTo6PNeoRGar3yDOk4zOOFFTBzuZoIiKSEUMQdZnUDzRAxjPDzpcSZ4RBp0GN2Yqc0lqPPYeIiNSHIYi6zJNN0RKdVoMhCVJzNPuCiIhIPgxB1GWZXghBQPOUGHeOJiIiOTEEUZd5emWYZFhTczTPECMiIjkxBFGXWKx25JbVAfB8CJKWyR86UwlRZHM0ERHJgyGIuiSntBY2uwijQYeeEcEefVZqj3DoNAKqGqw427Qkn4iIqLsCKgSdPn0aTz/9NEaPHo24uDgEBwcjKSkJl112GRYtWoRDhw4pXWLAOOVcGRYGQRA8+iy9ToO+sY4T6qU+JCIiou6S99hvBb322mtYuHAhamtbL6POz89Hfn4+fvjhB1RVVeGVV15RpsAAk1fumArrExPmlecNjDfiZFENMs9VY1JqnFeeSUREgS0gQtDzzz+Pp556CgCQmpqK+fPnIy0tDSaTCaWlpdi7dy8+++wzaDQBNfClqIIKx7RUr8gQrzxvYI9wfHOo0Lksn4iIqLv8PgRt2rTJGYDmzJmDt99+G0FBQa2uufzyy/GnP/0JFotFiRIDUkFFPQCgV6Rn+4Ek0in1J85Ve+V5REQU+Px6aMRut+Pee+8FAIwYMQLvvPPOBQGoJb3eM0c7qFFBpSMEJXptJMgRgjKLarhCjIiIZOHXIWjjxo3IzMwEADz22GPQ6fx+YMtvSNNhCSbvhKB+sWHQCEB1gxVF1WavPJOIiAKbX4egTz75BAAgCAKuvvpq59fLysqQmZmJsrIypUoLaPUWG8pqHVOL3uoJMui06NvUhJ15jn1BRETUfX4dgnbu3AkA6Nu3L8LDw/HBBx9g2LBhiImJQWpqKmJiYjBo0CC8/PLLMJs5eiAXaSosTK9FRIj3Rt+kTRkzi9gXRERE3ee3Ichut+PYsWMAgNjYWDz44IO49dZbL9gL6MSJE3j00UcxZcoUVFRUKFBp4DnbNBWWGBni8T2CWkrtEQ6AewUREZE8/DYEVVZWwm63AwAOHjyIf/zjH0hISMD777+PsrIy1NXVYevWrRg3bhwAYPv27fj973/f7v3MZjOqqqpafVDbpJVhCV6aCpM4m6O5QoyIiGTgtyGo5aaIDQ0NCA0NxZYtW3DrrbciKioKISEhmDhxIjZv3owRI0YAAD777DP89NNPbd7vxRdfhMlkcn4kJSV55efwR2e8vDxekuJcJs8VYkRE1H1+G4KCg1v/Ar7rrrswaNCgC64LCQnBCy+84Pz8o48+avN+CxcuRGVlpfMjLy9P3oIDiDQSlOillWGSAXFGCAJQWd+Ikhru+URERN3jt2vKw8PDW30+ffr0dq+9/PLLodPpYLVasWvXrjavMRgMMBgMstYYqKRDTL21R5AkOEiLPtGhyC2tQ2ZRNeLC+f8XERF1nd+OBBkMBsTFNZ8h1dH0VXBwMGJjYwEAxcXFHq8t0DlHgrwcgoDmnaN5fAYREXWX34YgALj44oudr202W4fXSn/ODRW7RxRFZ09Qopd7ggDHGWIA9woiIqLu8+sQNHHiROfrrKysdq+rqqpCSUkJAKBXr14eryuQldVaYLY6VuX1NCkQgniGGBERycSvQ9ANN9zgfP3ZZ5+1e91nn33mXE102WWXebyuQCYdlxEXboBBp/X68wfGO0aCOB1GRETd5dchaPjw4ZgxYwYA4MMPP8SmTZsuuKawsBB//vOfATgOUJ03b55Xaww03j449XwD4h1HZ5TWWlBaw13AiYio6/w6BAHAK6+8gsjISNjtdlx99dVYuHAhMjIysHv3brzxxhtIS0tDfn4+AOC5557jdFg3NS+P9/5UGACE6nXoHeUIYBwNIiKi7vD7LuHU1FSsX78eN954I86dO4clS5ZgyZIlra4RBAFPPvkkFixYoFCVgUPJlWGSgfFG5JfXI7OoBmP7xyhWBxER+Te/HwkCgAkTJuDw4cN4+umnMWLECERERCA4OBj9+vXDvHnzsGfPHjz33HNKlxkQCiqU2SOoJekMMY4EERFRd/j9SJAkJiYGixcvxuLFi5UuJaBJPUHePjKjpRSuECMiIhkExEgQeY9PTIfxNHkiIpIBQxC5zGK1o6jasSIrwcvnhrUkjQQVV5tRUcczxIiIqGsYgshl56oaIIqAXqdBTJhesTqMBp1zdRr7goiIqKsYgshlZ1osj9doBEVrSeGUGBERdRNDELnsrMIbJbbE4zOIiKi7GILIZdLyeCX7gSRSX1BWca3ClRARkb9iCCKXSdNhSi6Pl/SLdRyfkVXC6TAiIuoahiBymS8sj5f0j3OEoPzyejQ02hSuhoiI/BFDELnsrA/sFi2JMxoQbtBBFIHc0jqlyyEiIj/EEEQu86WRIEEQnKNBWcWcEiMiIvcxBJFLqhoaUW22AgASfaAnCAD6xzU1R5ewOZqIiNzHEEQukUaBIkODEKr3jSPn+kvN0VwhRkREXcAQRC5xToX5wPJ4SfNIEKfDiIjIfQxB5JICH2qKljT3BNVCFEWFqyEiIn/DEEQuaW6K9o1+IMCxV5AgAJX1jSir5UGqRETkHoYgcokvrQyTBAdpndNzbI4mIiJ3MQSRS3xxOgwAl8kTEVGXMQSRS3zpyIyWuEKMiIi6iiGIOiWKIoqqHSNBPSJ8LAQ1rRA7xRBERERuYgiiTpXXNaLR5lh9FRduULia1qTpsGwukyciIjcxBFGniqvNABwbJRp0WoWraU0aCTpdVgerza5wNURE5E8YgqhT0lRYvI+NAgFAQkQwgoM0aLSJyCuvV7ocIiLyIwxB1KmiKsdIUHy4b/UDAYBGI6BfbNPO0VwhRkREbmAIok4VVUshyPdGgoDWO0cTERG5iiGIOiVNh8VF+GYIGiAtk2dzNBERuYEhiDrVPBLke9NhAJfJExFR1zAEUaeKm3qCfG15vKQfN0wkIqIuYAiiThXX+EdPUEmNGVUNjQpXQ0RE/oIhiDpVVOW7S+QBIDw4yDlKlc3RICIichFDEHWo1mxFrcUGAIj3sSMzWurP5mgiInITQxB1SGqKDtVrYTToFK6mfVJzNPuCiIjIVQxB1CFfnwqTDOBeQURE5CaGIOqQry+Pl0jN0ae4azQREbmIIYg6JIUgX90oUdK/6eiMnNJa2O2iwtUQEZE/YAiiDjl3izb6dgjqHRWCIK2AhkY7Cip5kCoREXWOIYg6VCxNh/n4SJBOq0Gf6FAAQHYJ+4KIiKhzDEHUoWI/6QkCmleIMQQREZErGIKoQ0VVvr1bdEv9eXwGERG5gSGIOiT1BPn6dBjQfIYYR4KIiMgVbu9+9+yzz3qiDixatMgj96Wus1jtKK9znMXlD9Nh/bhrNBERucHtELR48WIIgiB7IQxBvkc6ODVIKyAqNEjhajon9QTll9fDbLXBoNMqXBEREfmyLk2HiaLY4UdXriHfI+0WHWc0eCT4yi3WqEe4QQdRBE6X1ildDhER+Ti3Q5Ddbm/3IysrC2lpaRBFETNmzMAnn3yC3NxcNDQ0oKGhAbm5uVizZg1mzJgBURSRlpaGrKws2O12T/xs1E3Fzo0SfX8qDAAEQUA/6fgM9gUREVEnZDsRs7KyEtOmTUNOTg5WrVqF22677YJrkpKSkJSUhOuvvx7//e9/MXfuXEybNg27d++GyWSSqxSSiXO3aB/fKLGlfrFhOJBfyRViRETUKdlWh/3973/HqVOnMH/+/DYD0PluvfVWzJ8/H6dOncKyZcvkKoNkVOQnGyW2JB2fkc3maCIi6oRsIejTTz+FIAi46aabXP6e2bNnAwDWrl0rVxkko+Jq/zhBviVpOozL5ImIqDOyhaCcnBwAcGtaS7o2NzdXrjJIRs0bJfpHTxDQvGEiQxAREXVGthAUFORYQn3w4EGXv0e6Vvpe8i3O6TA/Ggnq2xSCSmosqKxvVLgaIiLyZbKFoBEjRkAURSxduhR1dZ0vT66rq8PSpUshCAKGDx8uVxkkI3/aLVpiNOjQo6lejgYREVFHZAtB8+fPBwAcP34ckydPxr59+9q9dv/+/UhPT8exY8cAAHfffbdcZZBMbHYRJTUWAP41HQa0PD6DzdFERNQ+2ZbI/+53v8O6deuwZs0a7NmzB6NGjcKwYcOQlpaG+Ph4CIKAc+fOYdeuXa2mzK6//nr87ne/k6sMkkl5nQU2uwhBAGKMeqXLcUu/WCN2ZpUhm8vkiYioA7KFIAD48MMPkZiYiNdffx12ux0HDhxos0dIFEUIgoAHHngAf/vb3+QsgWQiNUVHh+oRpPWvc3YHNK0QO8XpMCIi6oCsv920Wi1eeeUV7Nu3D/fddx9SU1MBtD4eIyUlBffeey/27t2Lf/zjH9DpZM1hJBOpHyjOj5qiJc7pMI4EERFRBzySQIYOHYrXX38dAGCxWFBeXg4AiIyMhMHgf79U1ah5o0T/6gcCWvYE1TpHHYmIiM7n8WEYvV6PHj16ePoxJLNiP1weL0mKDoVOI6C+0YZzVWb0NPlfkCMiIs/zr2YP8hrpBHl/DEFBWg36RIcCALKKuUKMiIja5pGRIJvNhnXr1uG7777DoUOHUFZWBgCIjo7G0KFDMXXqVFx77bXQarWeeDzJwB83SmypX2wYskpqkVVSi0tTYpUuh4iIfJDsIejbb7/F3XffjTNnzji/JooiAEAQBGzfvh1vvvkmevfujTfffBNXXHGF3CWQDIr9uCcIaN0XRERE1BZZp8Pee+89XH311Thz5oxzNVhycjLGjRuHcePGITk5GYAjFOXl5eGqq67Cf//7XzlLIJn4/UgQD1IlIqJOyBaCcnNzcffdd8NutyM0NBTPP/88CgsLkZWVhe3bt2P79u3IyspCYWEhXnjhBRiNRtjtdsyfPx+nT5+WqwySgSiKzUdm+Nlu0ZL+sUYA7AkiIqL2yRaCXn31VZjNZhiNRmRkZOCJJ55AfHz8BdfFxcVh4cKFyMjIgNFohNlsxquvvipXGSSDarMVDY12AP65TxAA9G8aCcorr4fFale4GiIi8kWyhaCNGzdCEAQ8+uij+NWvftXp9SNGjMCf/vQniKKIDRs2yFUGyUDaLTrcoEOI3j+b1+PDDQjVa2Gzi8gr7/xAXyIiUh/ZQpA0pTV16lSXv2fatGmtvpd8g3O3aD86Pf58giA4m6OzuHM0ERG1QbYQZLPZAMCtZe/SkRl2O6crfIk/b5TYUv84R18QT5MnIqK2yBaCevXqBQDYvn27y98jXZuYmChXGSSD5hDkn03REi6TJyKijsgWgtLT0yGKIpYsWYKCgoJOry8oKMCSJUsgCAKmTJkiVxkkAykE+WtTtKR/Uwg6xekwIiJqg2wh6A9/+AM0Gg2Ki4sxduxYrFmzxjlF1pLdbseaNWtwySWX4Ny5c9BoNHjggQfkKoNkUFzjCEGxRj8PQXHsCSIiovbJtmP00KFD8dxzz+HJJ59EQUEBbr75ZkRGRuLXv/414uPjIQgCzp07h19++QUVFRXOXaSfe+45DB06VK4ySAYBMxLU1BNUUmNGZX0jTCFBCldERES+RNZjMxYuXAiTyYQFCxagrq4O5eXl2Lx5c6trpPATGhqKl156Cffee6+cJZAMSmosAIBYo17hSrrHaNChR4QB56rMyCquwcg+UUqXREREPkT2U+Tvu+8+5OTkYOnSpZg2bRp69uwJvV4PvV6Pnj17Ytq0aVi6dClycnIYgHxUSYBMhwHAgKbRIPYFERHR+TxyinxsbCweffRRPProo564PXmQzS6irNYxEuTv02GAoy9o+6lSHp9BREQXkH0kyBc89thjEATB+fH9998rXZLfKK+zwGYXIQhAdJh/T4cBLUeCGIKIiKi1gAtB+/btw9/+9jely/Bb0lRYVKgeQVr/f3tIzdFcIUZEROfzyG+5srIyLFu2DDNmzEBSUhLCwsIQFhaGpKQkzJgxA8uWLUNZWZnsz7Xb7bj77rthtVrbPLyVOldSHRhN0ZIBTcvkc0prYbVxZ3IiImomewhasWIFkpOTsWDBAmzcuBFnzpxBfX096uvrcebMGWzcuBELFixAcnIy3nzzTVmf/Y9//AO7du3C4MGDceedd8p6b7UIpKZoAEg0hSA4SINGm4j88nqlyyEiIh8ia2P0kiVL8OSTTzqXwZtMJowcORI9e/YEABQWFmLv3r2orKxEbW0t7r33XlRUVGDBggXdfvbp06fx1FNPAQCWL1+OLVu2dPueahQoewRJNBoB/WKNOHq2CqeKa9C3aRdpIiIi2UaCDh06hKeeegqiKCIhIQHvv/8+ioqKsHnzZnzwwQf44IMPsHnzZhQXF+O///0vEhMTIYoi/vznP+Pw4cPdfv7999+PmpoazJ07F5MmTZLhJ1KnQBsJArhzNBERtU22EPT666/DZrMhLi4OO3bswO9+9zsEBV24Q69Op8Mtt9yCHTt2ID4+HjabDa+//nq3nv3xxx/jyy+/RHR0NF5++eVu3UvtAuXIjJYGOM8Q4woxIiJqJlsI2rx5MwRBwMKFC9GnT59Or09KSsJjjz0GURSxadOmLj+3oqICDz74IABg6dKliI2N7fK9qHk6LFAaowFgQDxXiBER0YVkC0FnzpwBAFx66aUuf8/48eMBwKVT59uzYMECFBYWYvz48WyGloF0ZEag9AQBQP/YphBUwpEgIiJqJltjtFarBQBYrVaXv0c6ZV6j6VoWy8jIwNtvvw2dTofly5dDEIQu3QcAzGYzzGaz8/Oqqqou38ufBXJPUEmNBZV1jTCF8iBVIiKScSRImgJzZ2pLutaV6bPzWSwW3H333RBFEQ8//HC3T6J/8cUXYTKZnB9JSUndup8/CrQjMyRhBh16RgQDAE5xNIiIiJrIFoKmTZsGURTx8ssv4+DBg51ef+jQIbz00ksQBAHTp093+3l/+ctfcOzYMfTp0wdPP/10V0puZeHChaisrHR+5OXldfue/kY6MgMIjCMzWhoQ39QcXcQQREREDrKFoIceeggGgwE1NTWYMGECXn75ZZSWll5wXWlpKV5++WVcdtllqK6uhsFgwEMPPeTWs44dO4YXX3wRAPDaa68hLKz7e78YDAZERES0+lAbaSosOiwwjsxoqbkviM3RRETkIFtPUHJyMlasWIF58+ahpqYGjz32GB5//HH069cP8fHxEAQB586dQ3Z2NkRRhCiKEAQBK1ascHs67O9//zssFgv69++Puro6rF69+oJrDh065Hy9efNmFBYWAgBmzpwpS2gKRIF2ZEZL0vEZHAkiIiKJrDtGz5kzBzExMbjnnntQUFAAURRx6tQpZGVlAYBzJ2kASExMxJtvvokrr7zS7edIDcxZWVm45ZZbOr3+ueeec77Ozs5mCGpHIDZFS5wHqXIkiIiImsgaggDgqquuQk5ODj777DN89913OHTokPOw1OjoaAwdOhRTp07Ftdde2+ZmiqSc5j2CAjEEOYJvbtNBqroAm+4jIiL3yR6CAMeu0DfddBNuuukmT9weK1euxMqVKzu8ZvHixXjmmWcAAFu2bMHkyZM9UksgkUaCAmllmEQ6SLWh0Y688nr04xliRESqx78Ok1MgHpkhkQ5SBdgXREREDgxB5CTtFh2IjdFAc3M0d44mIiLAQ9NhAFBdXY3s7GxUV1c7d4buyMSJEz1VCrnI2RMUgNNhQHNz9KkiNkcTEZEHQtBbb72FN954AwcPHmy1GqwjgiC4ddwGeYazJygAp8MAjgQREVFrsoUgm82GG264AevXrwcAlwOQpyxevBiLFy9WtAZ/Yg/QIzNaGiCNBPE0eSIigowhaPny5fjiiy8AAD169MC8efMwatQoREdHd/mAVPKeQD4yQyKtCCurtaC81oKoAP05iYjINbKFoFWrVgEALrroImRkZCAqKkquW5MXSCvDokKDAu7IDEmYQYcEUzDOVjYgq6QGo8KilS6JiIgUJNtvu6NHj0IQBDz11FMMQH5IOjIjUKfCJCnxjimxk1wmT0SkerL/lX/QoEFy35K8IJCPzGhJCkEnzjEEERGpnWwhaODAgQDgPCKD/ItaQlBqj3AAQCZHgoiIVE+2EPTb3/4Woijiyy+/lOuW5EWBfG5YSwOl6bBz1QpXQkRESpMtBP3xj3/EiBEj8K9//QsZGRly3Za8pDiAzw1raWC8YySooLIB1Q2NCldDRERKki0EGQwGbNiwAaNGjcK0adOwYMEC7Nu3Dw0NDXI9gjwo0I/MkJhCgxDfFPS4XxARkbq5vUReq9V2eo0oili2bBmWLVvm0j25Y7TySgL8yIyWBvYwoqjajBPnqvGrpEilyyEiIoW4PRIkimKnH65ed/73kHKKA/zIjJakKTEukyciUje3R4KefvppT9RBClLDkRktScvkM9kcTUSkagxBpIojM1riMnkiIgI8sFki+R+pKTqQj8xoSVomn19ej1oze9GIiNQq8H/jUafUskeQJCpM71wFd6qYo0FERGrFEETO3aLV0A8kae4LYggiIlIrt3uCpNPiAWDOnDltfr0rWt6LvEstR2a0lNojHDuzytgXRESkYm6HoDvuuAOCIEAQhFbBRfp6V5x/L/KuYhWGoIFcIUZEpHpuhyAA7e7rw/1+/JOzJyg88FeGSVLiuUKMiEjt3A5B2dnZbn2dfJ+0OkwNGyVKUns4RoLyyutQb7EhRN/5TuhERBRY3A5BycnJbn2dfJ+ajsyQxBgNiA7To6zWglPFNRjay6R0SURE5GVcHUbNq8NUNBIEtFghVsS+ICIiNWIIUjm7XURprXSCvLpC0EAukyciUjWGIJVreWRGjFE9jdEAj88gIlI7t3uCfv/738tehCAIeOedd2S/L3VOWh4fHaZXxZEZLUkjQTxNnohIndwOQStXruzyfkBtEUWRIUhBRVWOEBSvoqZoSUrTCrHc0lo0NNoQHMQVYkREauJ2COrTp4+sIYiUJe0RpKYjMyRxRgNMIUGorG9EVnEtLkqMULokIiLyIrdDUE5OjgfKIKUUq3RlGOCYhk3tYcSunHJkFlUzBBERqYy6mkDoAmoeCQKad45mXxARkfowBKmc2kOQ1Bx9gmeIERGpDkOQyhVVNwBQbwiSlsmf4F5BRESq45EQtGnTJtx+++1ISUmB0WiETqfDkSNHWl2zbds2vPHGG3j//fc9UQK5SO0jQYMTHCEop7QWdRarwtUQEZE3dekU+fbU1dVh7ty5WLt2LYDmU+XbWk2m1WrxwAMPQBAEjB07FgMHDpSzFHKRFILUuEQecOySHWs0oKTGjOOF1RjZJ0rpkoiIyEtkHQmaPXs21q5dC1EUkZaWhj/96U/tXjt+/HgMHToUAPDpp5/KWQa5qKHRhqoGx+hHnDFY4WqUM6RpNOhYIfuCiIjURLYQ9Omnn+Lrr78GALz55pvYuXMn/vrXv3b4Pddffz1EUcTWrVvlKoPcIB2cqtdpEBEi66CgX7kowbE0/ujZKoUrISIib5ItBP3nP/8BANx222246667XPqeUaNGAQCOHj0qVxnkhqLq5j2C1LwBptQXxBBERKQusoWg3bt3QxAE3HzzzS5/T0JCAgCguLhYrjLIDWpvipYMaRoJOna22tnHRkREgU+2EFRaWgoASExMdP3hGsfj7Xa7XGWQGxiCHAbEGRGkFVBttiK/vF7pcoiIyEtkC0EmkwkAUFBQ4PL3ZGdnAwBiY2PlKoPcwBDkEKTVOHeO5pQYEZF6yBaCUlNTAQD79+93+XvWrVsHABg5cqRcZZAb1Hxu2Pm4QoyISH1kC0FXXXUVRFHEa6+9hoaGhk6vz8jIwOrVqyEIAmbOnClXGeSGoqqmPYIiGIK4QoyISH1kC0H3338/oqOjce7cOdx4440oKytr8zqr1Yq33noLV199Nex2O5KSknDHHXfIVQa5gSNBzQb3ZAgiIlIb2TaHiYiIwEcffYQrr7wS33zzDZKSkjBp0iTnny9YsAAWiwW7d+9GZWUlRFFEcHAwPv74YwQFBclVBrmhhD1BTtJ0WG5ZHWrNVoQZ1LtvEhGRWsi6Y/Tll1+OzZs3o0+fPqivr8e3337r3H/mm2++waZNm1BRUQFRFJGUlIQtW7ZgzJgxcpZALhJFkY3RLcQYDYgPN0AUgeM8UZ6ISBVkP0B1/PjxyMzMxKpVq3DjjTciOTkZISEh0Ov1SEhIwFVXXYUVK1YgMzMTY8eOlfvx5KKqeissNsfWBLGcDgMADGZfEBGRqnhkzF+n0+G2227Dbbfd5tL1BQUFbu0vRN1XVO1oXjeFBCE4SKtwNb5hSEI4tp0oxrGzHAkiIlID2UaC/vjHP3bp+/Lz8zF58mS5yiAXcSrsQlwhRkSkLrKFoNdffx2LFy9263vy8vIwadIknDp1Sq4yyEVcGXYhaYXYscJq2O08PoOIKNDJ2hP03HPP4bXXXnPp2uzsbEyaNAnZ2dlcHaYAjgRdqH9cGPRaDWrMVpyp4PEZRESBTrYQdM8990AURTz88MN47733Orw2KysL6enpyMnJgV6vx0cffSRXGeSiIoagCwRpNRjYwwgAOMIpMSKigCdbCHrjjTdw8803w263484778T69evbvO7kyZOYPHkyTp8+DYPBgDVr1mDWrFlylUEukkaC4hmCWuGmiURE6iFbCBIEAe+99x6uuOIKWK1W3Hzzzdi6dWura06cOIHJkycjPz8fwcHBWLt2La6++mq5SiA3cDqsbc4zxLhCjIgo4MnaE6TT6bB27VpceumlaGhowDXXXIM9e/YAAI4ePYr09HQUFBQgJCQE69atw4wZM+R8PLmBIahtzhVihRwJIiIKdLJvlhgSEoKvvvoKw4cPR3V1Na688kp8+umnmDJlCs6ePYvQ0FB8/vnnmD59utyPJjc4V4cxBLUibZiYW1qHGrNV4WqIiMiTZA9BAGAymbBhwwYMGDAAxcXFmD17Ns6dO4fQ0FCsX78eU6dO9cRjyUUWqx1ltRYAXCJ/vugwPXpEOP6ZHC/klBgRUSDzSAgCgB49euC7775DYmIiRFGE0WjEN998g/T0dE89klxUWusYBdJpBESF6hWuxvcMaRoNOlJQqXAlRETkSW4fm/Hss8+6df2ECRPw8ccfY8KECdiyZQu2bNnS5nWLFi1ytxTqIqkfKNZogEYjKFyN7xney4Tvjxdjf34lble6GCIi8hi3Q9DixYudJ8O7ShAEbNiwARs2bGj3GoYg72FTdMeG944EABzIr1C0DiIi8qwuHaAqijxSwJ8xBHVseJIJAHCyqAa1ZivCDB45Z5iIiBTmdk+Q3W73yAd5j3O3aDZFtyk+PBgJpmDYReDQGfYFEREFKo81RpPvcu4WHcEQ1J7hvR2jQQfyGYKIiAIVQ5AKcTqsc1Jf0H72BRERBSyGIBVybpTI6bB2cSSIiCjwMQSpEEeCOje8VyQA4HRZHcqbNpYkIqLA4vayl/79+wNwLHs/derUBV/vivPvRZ4jiiKKqhsAMAR1xBQahL4xocgprcPBM5WYmBqndElERCQzt0NQTk4OAFywV5D09a5wd98h6roasxUNjY7VeAxBHRveOxI5pXU4kF/BEEREFIDcDkFz58516+vkW6SpMKNBh1A997/pyPDeJnyxvwD72RdERBSQ3P4t+O6777r1dfIt7Ady3YikSADcOZqIKFCxMVpluDLMdRcnRkAjAOeqzDhX1aB0OUREJDOGIJUpquJIkKtC9Tqk9ggHAOzPq1C2GCIikp3b02GnT5/2RB3o06ePR+5LrTlHghiCXDKslwnHCqtxIL8S0y/uqXQ5REQkI7dDUL9+/WQvQhAEWK1W2e9LF5KmdRiCXDM8KRKf7MnnztFERAHI7RDkSyfI7969G19//TV++OEHHDlyBMXFxQgKCkJiYiLGjx+PO++8ExMmTFC6TJ9SWOkIQYmRwQpX4h9GNO0cffBMJURR5HYOREQBRLbVYd42ceJEZGRkXPB1i8WCzMxMZGZmYuXKlZgzZw7eeust6PV6Bar0PVII6hkRonAl/mFwzwjotRpU1DXidFkdkmPClC6JiIhkIts+Qd5WUFAAAEhMTMRNN92Eyy67DH369IHNZsOOHTuwbNkynDlzBqtWrUJjYyM++OADhStWniiKKKisBwAkmDgS5Aq9ToMhCeHYn1+J/fmVDEFERAHEb1eHDR48GB999BFOnz6NV155BTfccAPS0tIwbtw4PPzww9i3bx9SU1MBAB9++CG2bdumcMXKq6xvdO4W3ZMhyGXSifIHuEKMiCig+G0I+vLLLzF79mxotdo2/zw2NhbLli1zfr5mzRpvleazzjZNhUWH6REc1PY/N7oQT5QnIgpMfhuCXJGenu58zQNagbNNU2E9IzgK5A5p5+iDZyrRaLMrWwwREckmoEOQ2Wx2vm5vxEhNznJlWJekxBlhCglCfaMNhwuqlC6HiIhkEtAhaOvWrc7XQ4YMUbAS3+BcGcZ+ILdoNALS+kYBAHZllylcDRERySVgQ5DdbseSJUucn8+ePVvBanyDNBKUYOLyeHel9Y0GAPzEEEREFDDcXiLvL/7+97/j559/BgBcf/31GDVqVIfXm83mVtNnVVWBN+1xlsvju2xMP0cI2p1bBrtdhEbDTROJiPxdQI4Ebd26FY8//jgAID4+Hv/61786/Z4XX3wRJpPJ+ZGUlOTpMr3uLKfDumxoLxNCgrSoqGvEyeIapcshIiIZBFwIOnz4MK677jpYrVYEBwfjk08+QXx8fKfft3DhQlRWVjo/8vLyvFCt94ii6OwJ4nSY+4K0GozsEwmAU2JERIEioEJQdnY2pk+fjvLycmi1WqxevRoTJ0506XsNBgMiIiJafQSSqnor6iw2AFwi31XSlBibo4mIAkPAhKCCggJMnToVBQUFEAQB//73vzFr1iyly/IZZ6sc/UBRoUEI0XO7gK4Y09Qc/XN2mU8dJExERF0TECGopKQE06ZNQ1ZWFgDgtddew5w5cxSuyrc09wNxKqyrRvaJgk4joLCqAfnl9UqXQ0RE3eT3IaiyshJXXHEFjhw5AgBYsmQJ7r//foWr8j1nK6R+IE6FdVWIXothTUdo/MwpMSIiv+fXIaiurg5XXXUVfvnlFwDAk08+iccee0zhqnxTIZfHy0KaEtuVwxBEROTv/DYEWSwWXHfddfjxxx8BAA8++CCef/55havyXc0bJTIEdUdai74gIiLyb367WeItt9yCjRs3AgCmTJmCO++8E4cOHWr3er1ej9TUVG+V53MKq9gTJIe0vtEQBCCrpBbF1WbEhRuULomIiLrIb0PQ2rVrna83b96M4cOHd3h9cnIycnJyPFyV7yqocEyHJXIkqFtMoUEY1CMcxwqrsTunDDOGJShdEhERdZHfToeR60RR5G7RMnJOibEviIjIr/ltCBJF0a0PNY8CVZubN0rkbtHdJ22ayL4gIiL/5rchiFwnLY+P5EaJspBC0NGzVahuaFS4GiIi6iqGIBWQTo/ncRny6BERjD7RobCLwJ7ccqXLISKiLmIIUoFCLo+XnTQatONUqcKVEBFRVzEEqUCBFIIi2Q8kl8sGxgIAtp4oVrgSIiLqKoYgFXDuFs3pMNlMHBgHQQCOFVY7pxuJiMi/MASpAJfHyy8qTI8RvSMBANs4GkRE5JcYglSg+cgMTofJafKgOADA98cZgoiI/BFDkAo4G6MjORIkp8mD4gEAP2SWoNFmV7gaIiJyF0NQgKtuaESN2QqAS+TlNqyXCVGhQag2W7H3dIXS5RARkZsYggKcNAoUEaxDmMFvj4rzSVqNgImp0pRYkcLVEBGRuxiCApy0PD6Ry+M9YlJTCOJSeSIi/8MQFOCk5fFcGeYZ0kjQ4YIqFFU3KFwNERG5gyEowJ3lbtEeFWs0YHhvEwBgK1eJERH5FYagACcdnsrl8Z7DKTEiIv/EEBTgzlZxo0RPk/YLysgsgZVL5YmI/AZDUIBzHpnBEOQxI3pHwhQShMr6RuzPr1S6HCIichFDUIDjdJjn6bQaTJAOVOVSeSIiv8EQFMCqGxpRLW2UyJEgj5os7RfEviAiIr/BEBTACiqaN0o0cqNEj5Kaow/kVzo3qCQiIt/GEBTAsktqAQD9YsMUriTwxUcEY1RyFADg64NnFa6GiIhcwRAUwBiCvGvm8AQAwPoDBQpXQkRErmAICmDZJTUAgH6xRoUrUYcrhyVAEIC9pyuQX16ndDlERNQJhqAA5hwJiuNIkDfERwRjbL9oAMBXBzglRkTk6xiCApgUgvpzOsxrZo5IBMApMSIif8AQFKAq6xtRUmMBAPRlCPKaGUMToNUIOHSmyhlCiYjINzEEBSjpF3B8uIHL470oOkyP8SmOjRO/3M/RICIiX8YQFKCam6I5CuRtVzetEvuSfUFERD6NIShAZRc39QOxKdrrrrioJ4K0Ao6fq8aJc9VKl0NERO1gCApQWdwjSDGm0CDnDtKcEiMi8l0MQQGqeWUY9whSwtXDHavEvjxwFqIoKlwNERG1hSEoAImiyD2CFDb1oh4w6DTIKqnF4YIqpcshIqI2MAQFoKJqM+osNmg1ApKiQpUuR5WMBh2mDI4HAHzBKTEiIp/EEBSAspqaopOiQqDX8f9ipVz/694AgE9256Gh0aZwNUREdD7+hgxAPDjVN0wZHI9ekSEor2vkcnkiIh/EEBSAeHCqb9BqBPxubB8AwHs7cpQthoiILsAQFIDYFO07fpuWBL1Wg/35ldifV6F0OURE1AJDUADK4sGpPiPGaMBVTTtIr9qRq3A1RETUEkNQgLHa7DhdWgeAPUG+4vZLkgE4TpYvq7UoXA0REUkYggJMfnk9rHYRwUEa9IwIVrocAjAyKRJDe0XAYrXj4915SpdDRERNGIICjNQP1DcmDBqNoHA1BACCIGDOuL4AgPd35sJm5w7SRES+gCEowDj7gdgU7VOu+VUiIkODkF9ej++PFyldDhERgSEo4DQvj2cI8iXBQVrMHp0EgA3SRES+giEowDRvlMg9gnzNbWOTIQjA1hPFOHGuWulyiIhUjyEowGQXc7doX9UnJhS/ubgnAOBvG08oXA0RETEEBZB6iw0FlQ0AuEeQr3p4WioEAfj2cCEO5lcqXQ4RkaoxBAWQnFLHKFBkaBCiwvQKV0NtSe0Rjmt/1QsA8PLG4wpXQ0SkbgxBAYQHp/qHh6YOhE4jYOuJYuzKKVO6HCIi1WIICiAMQf4hOSYMNzWtFHtpw3GIIvcNIiJSAkNQADlV7Fgez34g3/fHy1Og12nwc3YZMjJLlC6HiEiVGIICyIGmRtvUHuEKV0KdSTCF4LaxjjPFXt7I0SAiIiUwBAWI0hozThY5RoLS+kYrXA254r70AQjVa3EgvxIbj5xTuhwiItVhCAoQu3LKAQCpPYxcGeYnYo0GzBvfFwDwwldHUWexKlsQEZHKMAQFCGmVEUeB/Mv/TRqARFMwTpfV4eUN3ECRiMibGIICxM/ZjhA0ph9DkD8JDw7CC9cPAwC8uz0be3LLFa6IiEg9GIICQI3ZisMFjqZohiD/kz4oHjf8ujdEEViwZj8aGm1Kl0REpAoMQQHgl9xy2EWgd1QIEkwhSpdDXfDU1UMQazTgVHEtXtucqXQ5RESqwBAUADgV5v8iQ/V4/tqhAIDlW7Nw6AzPFSMi8jSGoADwc1NT9Bg2Rfu13wztiauGJcBmF/HomgOwWO1Kl0REFNAYgvyc2WrDvrwKABwJCgSLr7kYUaFBOHq2Ck+tO8RNFImIPIghyM8dyK+ExWpHrFHPM8MCQFy4AX+7+VfQCMBHu/Pwzg/ZSpdERBSwGIL8nNQPlNY3GoIgKFwNySF9UDyevOoiAMBfvj6Kzce4mzQRkScwBPk5aZNEToUFlt+P74tbxiTBLgJ//HAfjhdWK10SEVHAYQjyYza7iD1Nx2Vwp+jAIggCnrlmKMb1j0aN2Yo7/7MLpTVmpcsiIgooDEF+7OjZKlSbrQg36DAkIULpckhmep0G/7p1FJJjQpFfXo873t2FslqL0mUREQUMhiA/Jk2FjeobBa2G/UCBKCpMj3fmjkZUaBAOnqnE7BU7cLayXumyiIgCAkOQH2vZFE2BKyU+HJ/83yVIMAXjZFENbvzXDmQV1yhdFhGR32MI8lOiKDpHgsayKTrgpcSHY829l6J/bBjOVNTjpuU7uKs0EVE3MQT5qZNFNSipsUCv02BYb5PS5ZAX9IoMwcf/dwkuToxAaa0Ft7y5E98cPKt0WUREfoshyE/9+0fHJnrjB8TAoNMqXA15S6zRgA/vHoex/aJRbbbi3v/+ggVr9qPWbFW6NCIiv8MQ5Ifyyurwye58AMD96SkKV0PeFhEchPfuHIt7Jw+AIAAf787Hlf/IwN7T5UqXRkTkVxiC/NA/t5yE1S7isoGxGM2maFXS6zR47DeD8eH8cUg0BSO3tA43Lt+BlzYcQw1HhYiIXBIQISg3NxePPPIIBg8ejLCwMERHRyMtLQ0vvfQS6urqlC5PVnlldVizxzEK9ODlAxWuhpQ2rn8MvnlwImaOSITNLuKfW05h0l+3YOWP2TyFnoioE4Lo58dUr1+/Hrfddhuqqqra/PPU1FR89dVXSElxb9qoqqoKJpMJlZWViIjwnY0IH//0AFbvysNlA2Px3p1jlS6HfIQoithwuBB//fY4skpqAQBJ0SF4ZNogXD08ATptQPx9h4ioU+78/vbrELR3716MHz8e9fX1MBqNWLhwIdLT01FfX4/Vq1fjrbfeAuAIQrt370Z4eLjL9/bFEJRXVof0l7+H1S7i03svwahkToVRa402Oz7enYdXv8tEUbXjmI0EUzBmj07Cb8ckIcEUonCFRESepZoQNHHiRGRkZECn02Hbtm245JJLWv35Sy+9hAULFgAAnn76aSxevNjle/tiCOIoELmqzmLFuz/m4J0fsp1HbWgEYMrgeNw4KgkTU2MRqtcpXCURkfxUEYJ+/vlnjB3rCAL33HMPli9ffsE1drsdQ4cOxdGjRxEZGYmioiIEBQW5dH9fC0EcBaKuMFtt2HD4HD74KRc7s8qcX9frNLikfwwuHxKPKYPj0TsqVMEqiYjk487vb7/9q+C6deucr+fNm9fmNRqNBnPmzMHChQtRUVGBLVu2YPr06V6qUD7ltRb8ed0h54owBiBylUGnxTUjEnHNiEScKq7BR7vy8O2hQpwuq8PWE8XYeqIYiz4/jF6RIRjZJxIj+0RhZJ9IXJQQgeAg7j9FRIHNb0PQDz/8AAAICwvDqFGj2r1u0qRJztc//vijX4UgURTx+b4CPPvlEZTVWqDVCPh/01KVLov81IA4I564cggWzhiMU8U1+O5oETYfLcLu3DKcqajHmYp6fHnAsQO1ViOgT3QoBsSFYUC8ESlxRiRFh6JXZAh6RARDr2OjNRH5P78NQUePHgUApKSkQKdr/8cYPHjwBd/jD/LK6vDkukPYdqIYAJDaw4gXrx+GkX2iFK6M/J0gCEiJD0dKfDj+b9IA1JitOJBXgb15Fdh7uhx7T1egtNaC7JJaZJfU4rujRed9PxBnNKCnKRjRYXrEhBkQY9QjOkwPU0gQwoN1CA9u+l+DDiF6LUL1OoQEaREcpIEgCAr95ERErfllCGpoaEBJSQkAoHfv3h1eGxUVhbCwMNTW1iIvL6/d68xmM8xms/Pz9pbcd1dGZjE2tfilIrVk1VlsKK4xo7ja8VFSY4ZdBPRaDf4wJQX3TBrAv32TRxgNOlyaEotLU2IBON6TxdVmnCyqwcniGpwsqsGp4hoUVDTgTEU9LFY7iqrNztVn7hAEx3vaoNNAr9M2/a8GOo0AnVaDIK2AIK0GWo0AnUaAVvoQBGg0AjQCoHG+bv5cgCPcCQKaXgMCmj4XAMdX4fzzlvU4/lQ47/O2avdeeGNOJLWYkBKLy4f0UOz5fhmCqqurna+NRmOn10shqKampt1rXnzxRTzzzDOy1NeRA/mVWLk9x6Vrx/WPxgvXDcOAuM5/RiK5CIKA+IhgxEcEO4ORRBRFlNVaUFDRgMKqBpTVmlFaa0FZjQVltRZUNTSiqsGK6gYrqhsaUWO2ot5ig7lp40ZRBMxWe9Pn3NmaSO1CgrQMQe5qaGhwvtbr9Z1ebzAYAAD19fXtXrNw4UL8v//3/5yfV1VVISkpqRtVtm1UchTuTx/Q6m+eAgBDkBZxRgPiwh0f8REGxIcHy/58ou4QBAExRgNijAYMg8nl77PZRdQ32lBnscLcaIfFZofF6vgwW+2w2uxotIuO/7XZYbWLsJ33YRcBuyg6PuwibGLzSKrj646QJUKEtOZVFMWmrzX/mcR5zfnFtrFg9vyvuLqmVrzw7kTUwth+MYo+3y9DUHBwcziwWCydXi9Nc4WEtL9RnMFgcIYlTxrXPwbj+iv7fzqRt2k1AowGHYwGv/xPDhEFKL9sMmm583NHU1yS2lrHMQKuTJ0RERGROvhlCAoODkZMjGM0JT8/v8Nry8vLnSHIE9NbRERE5J/8MgQBwEUXXQQAOHnyJKzW9hssjx075nw9ZMgQj9dFRERE/sFvQ9CECRMAOKa69uzZ0+51W7dudb4eP368x+siIiIi/+C3Iejaa691vn733XfbvMZut2PVqlUAgMjISKSnp3ujNCIiIvIDfhuCxowZg8suuwwA8M4772DHjh0XXLNs2TLnLtEPPvigy4enEhERUeDz21PkAWDv3r0YP3486uvrYTQa8cQTTyA9PR319fVYvXo13nzzTQBAamoqdu/e3WpVWWd87RR5IiIi6pw7v7/9OgQBwPr163Hbbbe1e8xFamoqvvrqK6SkpLh1X4YgIiIi/+PO72+/nQ6TzJw5EwcOHMDDDz+M1NRUhIaGIjIyEqNHj8bSpUuxd+9etwMQERERBT6/HwnyFI4EERER+R9VjQQRERERdQVDEBEREakSQxARERGpEo90bofUKtXeqjMiIiLyPdLvbVdanhmC2lFdXQ2Ah64SERH5o+rqaphMpg6v4eqwdtjtdhQUFCA8PByCIMh676qqKiQlJSEvL48rz8jv8f1MgYTvZ/8niiKqq6uRmJgIjabjrh+OBLVDo9Ggd+/eHn1GREQE/yWjgMH3MwUSvp/9W2cjQBI2RhMREZEqMQQRERGRKjEEKcBgMODpp5+GwWBQuhSibuP7mQIJ38/qwsZoIiIiUiWOBBEREZEqMQQRERGRKjEEERERkSoxBBEREZEqMQR1QW5uLh555BEMHjwYYWFhiI6ORlpaGl566SXU1dXJ9pxvvvkG1113HXr37g2DwYDevXvjuuuuwzfffCPbM4gAz76n6+rqsHbtWtx7771IS0tDVFQUgoKCEBMTg0suuQSLFy9GYWGhTD8Jkff+G91SXV0d+vfvD0EQIAgC+vbt65HnkMxEcssXX3whRkREiADa/EhNTRUzMzO79QybzSbeeeed7T4DgHjXXXeJNptNpp+K1MyT7+n9+/eLRqOxw/cyADEiIkJcvXq1zD8ZqZE3/hvdlkceeaTVc5KTk2V/BsmPIcgNv/zyixgSEiICEI1Go/jCCy+I27dvFzdt2iTOnz+/1b9kVVVVXX7O448/7rzXyJEjxQ8//FD8+eefxQ8//FAcOXKk888WLlwo409HauTp93RGRobzHuPHjxdffPFF8X//+5/4yy+/iBs2bBDvueceUaPRiABErVYrfv311x74KUktvPXf6Laeq9VqxeDgYDE8PJwhyI8wBLnhsssuEwGIOp1O3L59+wV//te//tX5L9nTTz/dpWccP35c1Ol0IgBx9OjRYl1dXas/r62tFUePHu2swxN/oyH18PR7+scffxRnz54tHj58uN1r1q1bJwqCIAIQBwwYINrtdrefQySK3vlv9PmsVqs4atQoEYD47LPPisnJyQxBfoQhyEU//fST81+ee+65p81rbDabOGTIEBGAGBkZKVosFrefc++99zqfs2PHjjav2bFjh/Oa++67z+1nEImi997TrrjhhhuctezZs8cjz6DAptT7edmyZSIAcdCgQaLZbGYI8jNsjHbRunXrnK/nzZvX5jUajQZz5swBAFRUVGDLli1uPUMURXz++ecAgMGDB2PcuHFtXjdu3DgMGjQIAPD5559D5Kbf1AXeeE+7Kj093fn61KlTHnkGBTYl3s+5ublYtGgRAGD58uXQ6/Xduh95H0OQi3744QcAQFhYGEaNGtXudZMmTXK+/vHHH916RnZ2NgoKCi64T0fPOXPmDHJyctx6DhHgnfe0q8xms/O1Vqv1yDMosCnxfr7vvvtQW1uL22+/HZMnT+7WvUgZDEEuOnr0KAAgJSUFOp2u3esGDx58wfe46siRI23eR+7nEAHeeU+7auvWrc7XQ4YM8cgzKLB5+/28evVqfP3114iKisKyZcu6fB9SFkOQCxoaGlBSUgIA6N27d4fXRkVFISwsDACQl5fn1nPy8/Odrzt7TlJSkvO1u88h8tZ72hX79+/HV199BQAYNmwYQxC5zdvv5/Lycjz00EMAgCVLliAuLq5L9yHlMQS5oLq62vnaaDR2er30L1hNTY3HniM9oyvPIfLWe7ozZrMZd911F2w2GwDghRdekPX+pA7efj8/+uijOHfuHC655BLMnz+/S/cg38AQ5IKGhgbna1ca3wwGAwCgvr7eY8+RntGV5xB56z3dmQceeAC7d+8GAMydOxczZ86U9f6kDt58P2/btg3//ve/odPpsHz5cgiC4PY9yHcwBLkgODjY+dpisXR6vdTkGRIS4rHntGwkdfc5RN56T3fkxRdfxNtvvw0ASEtLwz//+U/Z7k3q4q33s9lsxt133w1RFPHggw9i+PDh7hVKPochyAXh4eHO164Mn9bW1gJwbVi2q8+RntGV5xB56z3dnhUrVuCJJ54A4GhU/frrr1tN8RK5w1vv5xdeeAHHjx9HUlISnnnmGfeKJJ/Ufgs9OQUHByMmJgalpaWtmpfbUl5e7vwXrGXzsitaNvR19pyWDX3uPofIW+/ptnz44Ye47777AADJycn43//+h9jY2G7fl9TLW+/npUuXAgCmTp2K9evXt3mNdO/a2lqsXr0aABAfH48pU6a49SzyDoYgF1100UXIyMjAyZMnYbVa212CeezYMedrd1e5XHTRRW3eR+7nEAHeeU+f74svvsCcOXNgt9uRkJCATZs2dbqah8gV3ng/S1Nt7777Lt59990Ory0pKcEtt9wCwLE3EUOQb+J0mIsmTJgAwJHu9+zZ0+51Lfc7GT9+vFvP6NevHxITEy+4T1u2bdsGAOjVqxf69u3r1nOIAO+8p1vatGkTZs+eDavVipiYGPzvf//DgAEDunw/opa8/X6mwMAQ5KJrr73W+bq9vwHY7XasWrUKABAZGdnqKABXCIKAWbNmAXD8bWXnzp1tXrdz507n32ZmzZrF1QnUJd54T0u2b9+OWbNmwWw2w2QyYcOGDbj44ou7dC+itnjj/Sw6ztvs8CM5ORmAY6pX+tr333/fpZ+JvECJA8v8VXdPKN6yZYvzz+fOndvmM44fPy5qtdp2T5Gvq6trdYr8iRMn5PjRSKW88Z7eu3evGBkZKQIQw8LCxB9++EHmn4LIwRvv587wAFX/wp4gN7z66qsYP3486uvrMX36dDzxxBNIT09HfX09Vq9ejTfffBMAkJqaikceeaRLz0hNTcWjjz6KJUuWYPfu3Rg/fjwee+wxDBgwAKdOncLSpUuxd+9eAI4NuwYOHCjbz0fq4+n39KlTp3DFFVegoqICAPD888/DZDLh0KFD7X5PfHw84uPju/TzkLp547/RFGCUTmH+5osvvhAjIiKcf1s4/yM1NVXMzMxs83td/VuGzWYTf//737f7DADinXfeKdpsNg/9lKQmnnxPv/vuux2+j9v6aOtv6ESu8sZ/ozvCkSD/wp4gN82cORMHDhzAww8/jNTUVISGhiIyMhKjR492jtKkpKR06xkajQbvvPMOvvrqK8yaNQuJiYnQ6/VITEzErFmz8PXXX+Ptt9+GRsP/+6j7vPGeJvIWvp/JHYIoiqLSRRARERF5G4cSiIiISJUYgoiIiEiVGIKIiIhIlRiCiIiISJUYgoiIiEiVGIKIiIhIlRiCiIiISJUYgoiIiEiVGIKIiIhIlRiCiIiISJUYgoiIiEiVGIKIiIhIlRiCiIiISJUYgoiIiEiVGIKIiIhIlf4/hzJUgdYlqJUAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "theta = 0 # 未知パラメータ\n",
    "\n",
    "q_matrix = np.zeros([3,3])\n",
    "\n",
    "# 0: A, 1:B, 2:C\n",
    "def update_q_matrix(param):\n",
    "    q_matrix[0,0] = 1-param\n",
    "    q_matrix[0,1] = param\n",
    "    q_matrix[0,2] = 0\n",
    "    q_matrix[1,0] = 0.1\n",
    "    q_matrix[1,1] = 0.9-param\n",
    "    q_matrix[1,2] = param\n",
    "    q_matrix[2,0] = 0\n",
    "    q_matrix[2,1] = 0.1\n",
    "    q_matrix[2,2] = 0.9\n",
    "\n",
    "a_num = 30\n",
    "b_num = 50\n",
    "c_num =10\n",
    "\n",
    "a2b = 5 # A→Bへの移動数\n",
    "b2c = 1 # B→Cへの移動数\n",
    "\n",
    "x_min = 0.0\n",
    "x_max = 0.5\n",
    "x_num = 100\n",
    "\n",
    "x = np.linspace(x_min, x_max, x_num) # 0<Θ<0.5\n",
    "retult_list = []\n",
    "for theta in x:\n",
    "    update_q_matrix(theta)\n",
    "    result = \\\n",
    "        (q_matrix[0,0]**(a_num-a2b))*(q_matrix[0,1]**a2b)*\\\n",
    "        (q_matrix[1,1]**(b_num-b2c))*(q_matrix[1,2]**b2c)*\\\n",
    "        (q_matrix[2,2]**c_num)\n",
    "    retult_list.append(result)\n",
    "\n",
    "print(f\"最尤推定値(Θ):{retult_list.index(max(retult_list))*x_max/x_num}\")\n",
    "plt.ylabel(\"likelihood\")\n",
    "plt.plot(x, retult_list)\n",
    "plt.show()"
   ]
  }
 ],
 "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.11.1"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "949777d72b0d2535278d3dc13498b2535136f6dfe0678499012e853ee9abcab1"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}