1632 lines
124 KiB
Plaintext
1632 lines
124 KiB
Plaintext
{
|
||
"cells": [
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "23f67692",
|
||
"metadata": {},
|
||
"source": [
|
||
"Вар. 27 (513020125)\n",
|
||
"1. В результате эксперимента получены данные, приведенные в таблице 1. \n",
|
||
"a) Построить вариационный ряд, эмпирическую функцию распределения и гистограмму частот. \n",
|
||
"b) Вычислить выборочные аналоги следующих числовых характеристик: \n",
|
||
"(i) математического ожидания, (ii) дисперсии, (iii) медианы, (iv) асимметрии, (v) эксцесса, \n",
|
||
"(vi) вероятности P(X ∈ [a, b]). \n",
|
||
"c) В предположении, что исходные наблюдения являются выборкой из распределения Пуассона, построить оценку \n",
|
||
"максимального правдоподобия параметра λ, а также оценку λ по методу моментов. Найти смещение оценок. \n",
|
||
"d) Построить асимптотический доверительный интервал уровня значимости α1 для параметра λ на базе оценки \n",
|
||
"максимального правдоподобия. \n",
|
||
"e) Используя гистограмму частот, построить критерий значимости χ2 проверки простой гипотезы согласия \n",
|
||
"с распределением Пуассона с параметром λ0. Проверить гипотезу на уровне значимости α1. Вычислить \n",
|
||
"наибольшее значение уровня значимости, на котором еще нет оснований отвергнуть данную гипотезу. \n",
|
||
"f) Построить критерий значимости χ2 проверки сложной гипотезы согласия с распределением Пуассона. Проверить \n",
|
||
"гипотезу на уровне значимости α1. Вычислить наибольшее значение уровня значимости, на котором еще нет \n",
|
||
"оснований отвергнуть данную гипотезу. \n",
|
||
"g) Построить наиболее мощный критерий проверки простой гипотезы пуассоновости с параметром λ = λ0 при \n",
|
||
"альтернативе пуассоновости с параметром λ = λ1. Проверить гипотезу на уровне значимости α1. Что получится, \n",
|
||
"если поменять местами основную и альтернативную гипотезы? \n",
|
||
"\n",
|
||
"Таблица 1 α1 = 0.02; a = 0.00; b = 2.49; λ0 = 2.00; λ1 = 4.00. \n",
|
||
"0 1 2 0 0 7 1 0 2 1 0 1 2 2 0 0 1 8 0 0 14 4 3 0 0 3 0 6 2 2 1 0 0 2 0 4 0 0 3 3 1 1 0 0 6 8 1 \n",
|
||
"4 1 1"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 2,
|
||
"id": "57a523dd",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# Данные\n",
|
||
"import numpy as np\n",
|
||
"data = np.array([0, 1, 2, 0, 0, 7, 1, 0, 2, 1, 0, 1, 2, 2, 0, 0, 1, 8, 0, 0, 14, 4, 3, 0, 0, 3, 0, 6, 2, 2, 1, 0, 0,\n",
|
||
" 2, 0, 4, 0, 0, 3, 3, 1, 1, 0, 0, 6, 8, 1, 4, 1, 1])\n",
|
||
"n = len(data)\n",
|
||
"alpha = 0.02\n",
|
||
"a = 0.00\n",
|
||
"b = 2.49\n",
|
||
"lambda0 = 2.00\n",
|
||
"lambda1 = 4.00\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "8b7561a0",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Пункт a)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "b046ad70",
|
||
"metadata": {},
|
||
"source": [
|
||
"### 1. Вариационный ряд"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 3,
|
||
"id": "db7e1a67",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Вариационный ряд: 0^(19), 1^(11), 2^(7), 3^(4), 4^(3), 6^(2), 7^(1), 8^(2), 14^(1)\n",
|
||
"[ 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1\n",
|
||
" 1 1 1 1 1 1 2 2 2 2 2 2 2 3 3 3 3 4 4 4 6 6 7 8\n",
|
||
" 8 14]\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"# Получение уникальных значений и их частот\n",
|
||
"unique_values, counts = np.unique(data, return_counts=True)\n",
|
||
"\n",
|
||
"# Форматирование вариационного ряда\n",
|
||
"variational_series = [f\"{value}^({count})\" for value, count in zip(unique_values, counts)]\n",
|
||
"print(\"Вариационный ряд:\", \", \".join(variational_series))\n",
|
||
"sorted_data = np.sort(data)\n",
|
||
"print(sorted_data)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "93c7e45f",
|
||
"metadata": {},
|
||
"source": [
|
||
"### 2. Эмпирическая функция распределения (ЭФР)\n",
|
||
"$$\n",
|
||
"\\hat{F}_n(x) = \\frac{1}{n} \\sum_{i=1}^{n} \\text{\\textbf{1}}_{\\{X_i \\leq x\\}},\n",
|
||
"$$\n",
|
||
"где $n$ — объем выборки."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 4,
|
||
"id": "261ad18a",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA04AAAIjCAYAAAA0vUuxAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvc2/+5QAAAAlwSFlzAAAPYQAAD2EBqD+naQAAfHZJREFUeJzt3QmYFNW5//F3hgEGFAVBQBFFcQ/iFsA1JmogalySm6iJibtmUeOSGDW5LolJzKLG3Gjct8QkavzH4BUUVKIQQSGI6KgooEFEZBFG1oFZ+v/8DrfGnqZnehp56Toz38/zNF309PKr01XVfeq8VV2WyWQyBgAAAABoVnnzfwIAAAAACB0nAAAAACiAjhMAAAAAFEDHCQAAAAAKoOMEAAAAAAXQcQIAAACAAug4AQAAAEABdJwAAAAAoAA6TgAAAABQAB0nAAAAACiAjhOAjaKqqsrOPfdc23nnna1z58625ZZb2gEHHGC/+93vbM2aNaWOhxTYfPPN7fTTTy91DAAANkjFhj0MAD62fPly22effWzbbbe1k046yXbbbTdbvXq1TZgwwS655BL7y1/+Yk888YRttdVWpY4KAACwQeg4AfjEGhoa7KKLLrKf//znYbQpccEFF4QO0zHHHGNnnnmm/eMf/yhpTgAAgA1FqR6AT0xleddff32TTlPiqKOOshNPPNFGjhxpU6ZMabx9wIABVlZWFjpcuUaMGBH+9sUvfrHxtmeffTbc1twluwTsvvvuC7f95z//adK5Gzx4cLhdf0/ocSohe/vtt8PrbrbZZmHk7Kc//allMpnG++m5ch8r55133nqvr2nNXy7d75prrmly27x580Knsk+fPqH9PvWpT9k999yz3mNramrCY3fddVerrKy0bbbZxr785S/b7Nmzm82nkcD999/fdtxxR5s/f37j7XqvDjroIOvZs6d16dIl3OeRRx5Z7zVXrFhh3//+922nnXayjh07NmnvxYsXW0vU3pdeemlYNtQWTz75ZOPfLrvsMuvWrZvtsssuoWOduPfee8NzT5s2bb3n+8UvfmEdOnQI7SWf/exnbdCgQevdT/OW+97r9XNLBFVWqnbUcpV9v+xlLnH++eeH58ym9yL3NrVX3759w+3Zzyu33npryNu1a9cm7Ziv3fO9zowZM8J6tMUWW4T37cILLwzLRDa13+GHH269e/cOy9Kee+4ZXjcftfthhx0W3gc955AhQ8LIcELt29L6ltu+arexY8eGkWe1q17773//+3qvW11dHdb5/v37h4wq7f3Vr34VlpdcyXqce8m3bql9vvKVr4RRbb3+pz/9aXvsscfyzntz85a7br/44ov2hS98ISzDet/UXs8//3ze9yd3ffj3v/+dd1uTm33u3LlhHcxt07q6OvvZz34W1ne1U3ZOPTeA0qDjBMCdOgaS+0VGX3D+/Oc/W21tbeNt7733nj3zzDPhb/l873vfsz/96U9NLvk6bLl0v1dffTXv3+rr68MXJHVefv3rX4eOxNVXXx0uLZk1a5bdeeedtqEWLFgQjgN7+umnw5dzHQ+mL5JnnXWW3XTTTU3y6YvpT37yk5DthhtuCF+cP/roo3BsWT5q0//6r/+yd99918aMGRM6Wgm9zr777hs6h+qQVFRU2Fe/+lUbNWpUk+dQx+fGG28MX8bvuOOO0IZf+tKXWjVv+jKsTszxxx9vF198cbisXbs2vMZLL70URif1hVGdv3feeSc8Rl98dZuWiVy6TV94+/XrZ5+U3te7777bHnjggfCcG4veF72nuR566CH77ne/a1tvvbX99re/De34ox/9qKjnVqdJHaXrrrvOjj76aPuf//mf0PnLpk7SDjvsEJ5bWdQ50evecsstTe6nL/MaBV6yZIldccUV9stf/jJ0eLI7t7Lddtutt6597Wtfy5tv5syZoUxXO0qUMVmmnnrqqcb7rFq1KnQ+1O6nnnpqmIeDDz44ZFBJb3O0nCavr5y5XnvttbAevfHGG3b55ZeHedcOkBNOOMEeffTRvM+5++67Nz6n3pNc48aNs8985jO2bNmysLxoPVGnT+vC5MmTbWO56qqr1usAi+bhyiuvDJ3tP/zhDyFn7vsNoAQyALCRrFy5MrNo0aL1LjNmzNDQTebLX/5y43132GGHzOc///lMr169Mo888kjj7ddee23moIMOCn8/5phjGm//5z//GZ7jb3/723qvu9lmm2VOO+20xv/fe++94b7vvPNO+H9NTU1m++23zxx11FHhdv09ocfptgsuuKDxtoaGhvDanTp1CvlFz5X72BNPPDEzaNCgTP/+/Zu8/hlnnBFeL5cef/XVVzf+/6yzzspss802mcWLFze538knn5zZcsstM6tWrQr/v+eee8Jjb7zxxvWeU1lz8+m2U045JdO1a9fMiy++uN5jkudNrF27NszH4Ycf3uR2ZRsxYkST25Rfr5O0Sz5q7969e2e+9rWvNd42ffr0TIcOHTJ77713Zs2aNeE2zXe3bt0yF154YeP99Jhtt902U19f33jbSy+9tF7bH3bYYZlPfepT6732b37zmybvvWhZSt6f22+/Pfz997///XqPzV3mEuedd154TL52SCxcuDDMS7KMaXnNnqfu3btnVq9e3arlOd/rHHfccU1u/+53vxtuV7s2976K3r+ddtqp8f/V1dUh57Bhw5rkyV6WNqR9ddv/+3//r/G2jz76KCw/++67b5N1W+vqW2+91eQ5L7/88rBsvPvuu01uv+OOO8Lz/vvf/268Te+PXi/bEUcckdlrr73Ccpc9L9qO7LLLLuvNw8EHH5z53Oc+1/j/3HVbj9Xj1HbZbaL23XHHHcN2q9D6MGXKlLzbmuzsVVVVmfLy8sZlJrtNDzzwwMwee+zR5PWT7ZqeG0BpMOIEYKPRaI32qudetHdXtPc2W6dOneyUU04JJUbZe8PPOOOMjZpLe9w//PDDFkeQNOKTUDmM/q8REo0G5TN16lT729/+Fvaul5c33ZSqVGrhwoXh8c1RP+r//b//Z8cee2yYVqlPclHJoEaTNDIjul+vXr3CMWO5csvFkpEijdA8/PDDNnTo0PX+rlGdxNKlS8NrHXrooY2vl13qp7KwYmlkT/Ov0aSEyiQ1iqgRA73voufWXn2NMCY0EvH+++/bP//5z8bbNC/KrBG0bBqJy243XTSq0RyVi2oERu2T/X5vDNdee20o6dKIaC61o0q9mhtFbQ2VhGZLloXRo0fnfV/1nqo9NMKjMlT9XzQCpDwamcnNk29Zai2Vt2aPRqr8T++lyi4/+OCDcJvWFy1nPXr0aPKeHXnkkeG9HD9+fJPnTEZiWmo3jZppdEgjcpqv5Dm1vms90khYUt6Z0HrZ0ij1yy+/HB739a9/PTxP8pwrV660I444IuTMLS1Ujux5Stq7JRpp22+//cLIXC7Ni9rpk7wnADY+Tg4BYKPRF6VDDjlkvdt17Ie+VOnLVC51klR+pmNw3nrrrXCtL0Gq798Y9AVGZTYqBVIpXj7q+Og4nmw6tkCyjzvIpi+e+hKoErrcL+E6fkilav/93/8dvkjn++K3aNGiUPqjEjhd8lHnQ3Qck85UqPKnQm6//XZ74YUXGjtF+Tz++OOhffUFMftU8blf0g488MBQ6qTjcFRSpeOcWuqYZB+3Ia0pq9N9/vWvfzX+//Of/3woK1RnSV9S9QX1r3/9ayj50/E4uce1qGPeGppXdST1BV1fcjcmlRqq3VUql++9VjuqzXU8jMpW1YlqzRfrbDoeLNvAgQPDcpu9fOr4G+0cmDRp0nrvk15PHbvkmLh8x4d9EioxzV1+stchHfulzsgrr7zS7HuWLO+J5Lgh5W6pXFY7HlTWpktzz5u9LGq9U0ljc5RTTjvttGbvo/ZUxyah9bMYWub/93//N+w0UDltvmXmrrvuCsuVtjHq6Gk7CqC06DgB2GjU+cjtgEhyoPwee+yx3t/23nvvcPnjH/8YjlHQqEK+DtaGUgdGXzA1yqC9xxuDDoLXSJS+oOZz3HHHhS/Iv/nNb8Iln2SP9Te+8Y1mv6BplKZY6jTp+CGdiEPHFenYLY1WJXSKeOXTSI+OnVAnRR0ijfplnxxA1KHTMS359oi3JN8xGy3RqesTOgGE9vTr2DHlU2dAI1Bqp1w60D73GDONauTriE6fPj0cf6POmJYFPd/GOr7pxz/+cejY6H1U++bS+/Dmm2+GUSkdp7Yx5HZS1CHSvGl0V8el6fgmjexpRErH8OQ7+cKmpgzqGP/whz/M+/eko5VQh0vLpkazWnpO+cEPfhBGmJrr1GXTCFhz981+Tq27+Y6pEp1QJptGhbO3W9oJlDtKmE0nSFEGHTOVe1IK0Ui2Rsq+/e1vN/scADY9Ok4A3OlgcNEX9nzUydCXO32h0V7YjUVfuHUiBH0J0WhFcx0nfVFSOVP2Fzd98ZHcs2Bp77ZGmzSCpgPSm6OTD+jAb32hTb6I6UtjQnvdlUkjICpVaolGF3SGL53wQV8kW6K21MkBNO86s5m+tOvA8uwveBoV0QkjssuVssslE5p3vXd77bVXeF4dbK8Obvbz5ZOciEIZCtGXw9wvxhq51MHxWhZ09je1Vb4vujoBQG7baWQpH82DOlUqZ9O1DrTX6McnKZ8TlaI9+OCD4VT76vTlo9dUB0/31eiJRoXUkdOX/dbSKIjOjpg90qLlKlk+1VYaPdQJWLbffvvG+2WXPCbLkuikIrkdik8iGfnJ7tDlrkN6bY2aFFreEzp7nErZckthsyU7arRetOZ5dfIZlcHl24mT20bqCLU2q3ZEZO+g6N69e7P31bKinS65pbHZVMaq9Uxn2dQo/re+9a2ww6a5HTEANg2OcQLwienYJZUhZZ8dL/uLm76U6yxe+Y63EY0w6Au0jg3amGc50959lee1Zq/tzTff3DitL4D6v76MaS9+Nn1J1hdudcYKUTmQ9ijry1fuFzB9ydbomjoy+c6Mp1K+hO6nsqXsjNlZs6l8UNQZ0WibOj76wpX9uvpyqw5b9p79fL+xpVMi6xg0fXlTx1bzkG9EMZdOba3OQvYZzdRmGolSxyY59kslczpeRF86c0fadFGpktrn5JNPblWZYkv0BVwdLX0J1/NqnnW2tk9KnWiVMTa3UyD7eBaVZOn9UDuqPLUYuWfG+/3vfx+uNYomSacte3lQOVluh3j48OGhw67lN3dkMHdZKoY6ydnvt7YJ6mRrxEZleqISXHUY1GnPpfI5LW+J119/PVxUotmSZJuhkrbsU+7nW4+S9Ve0XjZH7406TzorZL7yuNznLIbWO+3Y0DavudGshDr3GjXU8qplRjtCAJQWI04APjF9EdcxM/fff384JbHKlvSlTGVW2ruvU1/nK0dJ6FgBfelJvtRvLOow6FiZ5GQEzdGog07FrFKrYcOGhVEOnTZbX3Byj8fQc55zzjlFH9OQj04DrY6lXlPPqS9G6kxoT7RKAZNjcTQCoy+hOk5Lp0JW50gHqus+OtlBc18u9cVL5XfqOKpzpmNr1IFVKZdK+PTlTcd/6Eu5Rh/UucnteOpEDxopKTTSlU0dFJ0uXfOnDo86LbfddlvotOh9VgZ1NPSFUKMk+UZeNM/J7fnK9D4JHd+jUinlU6csuyRSX4pzT8udHIOi29X5128FZS8Pub/tk0vvU3Ia8paOrSl0HJXaTO+bOh/qgOn9U5lr0iHScq6TjWh0Ql/4NcqljkV2h0KjKMpy9tlnhw6unkPrn0bAdFyU1uENodFanUZfJaLaWaHfItOp2bM7biqR1IiYjtnRbxqpg6LlWMuYjqNTZ1ajNupYJe+9OuDJiLVoB4seo9uS5ULLr0ZlNKqo9Uide7222kkjTJo3/V8jfVrm9J4nJ6zJJ+lcq1OqnQY6DlPHSOm1tb6qDTd0ZFx5khLKlmjEWh1RvV5Lx3gB2MRKdDY/AG3Ma6+9lvnWt74VTn2s03jrlMdDhgwJp9DOPe1xS6d+bu7vG3I68n322afJ6XzznVJcj9PjZ8+enRk+fHg4hXefPn3CaYazT4mdPLZLly6ZefPmrZc1+/Wbk3s6clmwYEE43bVOad6xY8dM3759w+mVdSrmbDoV8o9//ONwOuTkfl/5yldC7ubmTd58881MZWVl5uKLL2687e677w6nW+7cuXNm9913D4/JPb32hAkTwimidfruYk9HLrW1tZmLLrooLAc6NfuTTz7Z+D5ddtllmc033zwsK4899ljex8+fPz+8/q677pr37xt6OvKETl2tedcyWldX13g/PbalS3Ka8aQdjj/++CbPmyynyf10ynWdXj371Owbcjry119/Pbzfas8ePXpkzj///PXWK7Xl4MGDw/s9YMCAzK9+9avGU9lnt0dyX52uW8vzFltskRk6dGjmr3/96wa3r9bVMWPGhNdPlqt887Z8+fLMFVdckdl5553DdkI/R6Ac119/fTgtfvLahd6H3K8vWg9OPfXUsF5o/ejXr1/mi1/8YuNPHTz//PPhNa+55prG0+Enmlt3pk2bFn5CoWfPnmGeNJ/6CYJnnnlmg09HrtuyT7+f7+cTZs6cGdYVtVO++3E6cqB0yvTPpu6sAUBaaM+39nZzxip/OqBeP3Lb0uhjQqWJOlZKx4k1d7a0TU2joRoB2JjlpIWoBFYjfxoJyz6GJk10DJNG8XTmwI1B7auL5j0fjUzpeC++vgDY1DjGCQCQOupc6XiQb37zm6WOAgBAwDFOAIDU0I+Z6qQAOqW6zuKXe1bDUtKZ/bKPb4IPnX2ypbPeaeRSJy0BgE2NjhMAIDV0pruJEyeGM9UlZ45Li9yTRsCHfherJSpZzD5hBABsKhzjBAAAAAAFcIwTAAAAABRAxwkAAAAACmh3xzg1NDSEXzjXL6dvzB/aBAAAABAXHbW0fPly23bbbcMPYLek3XWc1Gnq379/qWMAAAAASIm5c+fadttt1+J92l3HSSNNSeNsscUWpY4DAAAAoESWLVsWBlWSPkJL2l3HKSnPU6eJjlPxamtrbezYsTZ8+HDr2LGjpR15fZHXF3l9kdcXeX2R1xd525+yVhzC0+5OR65e5ZZbbmkfffQRHadPUAcayzFi5PVFXl/k9UVeX+T1RV5f5G0/lhXRN6DjBAAAAKBdWlZE34DTkaPooeCRI0eG6xiQ1xd5fZHXF3l9kdcXeX2RF/kw4oSiaHGpqamxysrKKIaCyeuLvL7I64u8vsjri7y+yNt+LGPECZ4qKuI6pwh5fZHXF3l9kdcXeX2R1xd5kYuOE4pSV1dno0ePDtcxIK8v8voiry/y+iKvL/L6Ii/yoVQPRdHiopVSezViGAomry/y+iKvL/L6Iq8v8voib/uxjFI9eIptbwZ5fZHXF3l9kdcXeX2R1xd5kYuOE4peKfUDa7GsnOT1RV5f5PVFXl/k9UVeX+RFPpTqAQAAAGiXllGqBy/qZ2sBi6W/TV5f5PVFXl/k9UVeX+T1RV7kQ8cJRdEQ8IQJE6IZCiavL/L6Iq8v8voiry/y+iIv8qFUDwAAAEC7tIxSPXhpaGiwJUuWhOsYkNcXeX2R1xd5fZHXF3l9kRep6ziNHz/ejj32WNt2223DOef/8Y9/FHzMs88+a/vtt5917tzZdt55Z7vvvvs2SVasU19fb1OmTAnXMSCvL/L6Iq8v8voiry/y+iIvUleq98QTT9jzzz9v+++/v335y1+2Rx991E444YRm7//OO+/YoEGD7Nvf/radffbZ9swzz9hFF11ko0aNshEjRrTqNSnVAwAAAFBs36DCSuioo44Kl9a67bbbbMcdd7Qbbrgh/H+PPfawf/3rX/bb3/621R0nfDIaAl68eLH16tXLysvTX+lJXl/k9UVeX+T1RV5f5PUVW94JExps1qzlNnRoNxswIH/eDh3MKis//v/Klc0/n2a5S5cNu++qVTrLX/77lpWZde1q0Sppx6lYkyZNsiOPPLLJbeowadSpOWvWrAmX7F6lzJs3z1auXNk4pLnZZpuFXqbORqIaUV2rfLBDhw5hul+/fmHFWbBgQfi/ppPrrbbayioqKqympsZWrFgRbtfj9Hj9vXfv3mEFfP/998PtotfVY7beeutw6sjly5fb6tWrw/10u66VqXv37iGnMiW361JZWRleV49dtGhR43zo+TWtFV33+fDDD8PzJvOhTJrPrl27Wm1trVVXVzfOhy56vr59+4ZptZEk86rnUI/81VdftX333TfMr14rmadOnTqF1127dm1op+y8mu7Tp0+YXrp0abhP9rwqU7du3UJvX++Rbk/mSVn1uvq/5id7XpVLeTt27GgLFy4M85T93my++eYh76c//WlbtWpVk3lVWyiTrufOndvkvdG03hvRe6o2TN4z5VWmHj16hNu1Yc2eV5WRqh00rbzKpHZN5qlnz57WpUuX0A56b5P3THS78g4bNiwsE9l5dR+VtWp6/vz54fmz51XPq+fKXQ4lWdb0HHqs/q9MyTwly6jaXo/Pfm/0vuj90XMqczIferze86qqKjv44IPD33KXQz2v7qNlVOth8p4pr95TLaNaFvS+N7ccvvfee43vVzJPal/9X+2n9zWZDz1O7a+20HzodbPfGz3PjBkzQl7Nq96bZF71HHpPtd5pvdD8Zr83pdhGiPIecsghYVlL+zZC3njjDTv00EPD+5psY9O6jdC8vP766/bZz342PEbPm+ZthC6vvfaaHXHEEeH++nuatxG6aHv2+c9/PrxOvuUwTdsIvQ/Kq+8Vej6tU2neRqi9lHf48OHhtrR/j9C8vfLKK+F7nN5zvS9p3kaofadPnx62D2pH/T/N24jf/rbBXnihPHRaPvpoS+vYsdZ69Vrc5DvxQQdl7OGH1703Ws8PPjhjNTVljX9furSH1dRU2uabr7DDDltut9768Tbi0EMr7e23t7IOHbTeLmzyvHvsYTZ27MfbiOHD19j8+R8/r/KsWtXVunZdZfvsU22PPFLW+D1C7aBltLnlcFNsI7SstcmO0wcffBAW0Gz6v94kLXhaiHJdd9119pOf/GS92//4xz+Gxk7stddeYTRLDTp69Oj17n/66afbDjvsYH/5y18aO1+JL33pS2EB1IL+3HPPNfmbVuLvfOc74c2/55571nveb33rWzZ16tTwnG+99VaTvynPqaeeai+99JKNGzeuyd+0MOjDUyvNmDFj1qtp1QfVQQcdZI899pi9++67Tf6mL41aULUQ6bHZtHCeddZZYUFVG2lhzHbaaaeFDcDkyZPtxRdfbPI3bYAvuOCCsID+6U9/avI3LaTnnntu+JI9c+bM8F5mGzx4cGhHHfem+c226667hvdHK8zjjz++XhueeOKJYfTx4YcfDhuYbBrR1IqhDxeVdmbThlCZtLHN996cc845jfXCenw2fQCqXFRfYlQqmk0bwuOPPz58YKgUVR8O2fQlWO+dlrNZs2Y1+Zs6eGp7fZlTKWtuG55xxhnhtR944IH1nvfkk08Oz6flUK+bTRuX8847Lywv+eb1u9/9bmh3Lcdz5sxp8rfddtstPLd2XEycOLHJ37ROHHDAAeF9zW0H0TGMOibx73//+3rv+eGHH974QfDUU081+ZvWZbWvlql77713vYNdzzzzzJBXy3Du8qJtgsp5Z8+ebY888kiTv2m5//rXvx6W35dffjm0R7YhQ4bY0UcfbU8//XR4b7OVchuh9T+mbYR+vV7zGss2Qsu15iuWbYS+fGrdiWUboTJ7fW7Eso1QHv0tlm3E0KFDo/oesffee4cOSSzbiG222SYsv2nfRuy997rL5MmfttGjjwmdpm99644mj6ur62jz55/WuI047bSmz/vXv55sb765m+2zzzQbMmScZc/SAQfsZm+/fbJtttnK9Z5Xliz5eBtx3HFNtxGPPXasvfTSfrb77jNs+PD/tTvuWP97RCm3EbnbrShOR64vt4WOcdKCrw+FK664ovE2LTTHHHNMWKjydZzyjTj1798/fDFVD5cRp+JGnLQy6/W0IU/2nqRxT1EyT3pevZ5yJ+9NGvcUZS+HWhb0+Ny9J2nbm6zHaxnT/dWOue9N2vYmJ/Ok51BeLftp3pus103ed+XVvKZ9G6HH6jFaTjWvad6bnIyQabncfvvtw+1p3pucrHNa3vVlOGmnNG8j9Dfdtssuu4Q2TPuIk+6j5X3PPfcMGdM+4pSsn/rCr/ul/XuE3kctH+poK0faR5x0rfvqRGRaH9M+4nT++Q02aVK5/fd/V9opp3QLOT/6SJUges/W5TWrtwEDPt5GrFzZdBvRo8e6bcSaNTW2atUKq6j4eBtRW1thPXuu20YsWNB0G9GxY4Vtv/3H24jq6nXbiA4d1s3T5pt3s27dtrCVK1eETJtvXpG6Eafdd9+9dec/yKSEojz66KMt3ufQQw/NXHjhhU1uu+eeezJbbLFFq1/no48+Cq+laxSvtrY289xzz4XrGJDXF3l9kdcXeX2R1xd5fcWWd8SI+oy+1d9zT12po0SnmL5BVCNOl112WRhhyh7uVNmN9qI8+eSTrXodzqoHAACAtkTnWtNX4fvvNzv11FKniUs0P4CroUvVEOsiqoPWdFJLq5I81eYmVJP89ttv2w9/+MNwwPQf/vCHUI968cUXl2we2hsNuaq+PbdWNK3I64u8vsjri7y+yOuLvL5iy5uMg8SSN1Yl7Tj9+9//Dmdn00UuueSSMH3VVVeF/6v+MvuARNV168AxHQCmgwt1WvK77rqLU5FvQlohVbMcy4pJXl/k9UVeX+T1RV5f5PUVW16zdR2nlBSStVmpKdXbVCjVAwAAQFtCqV47+AFcxEdnOFFJpUb/kjOtpBl5fZHXF3l9xZRXZ3Q+8siMVVXpTFUqFvn4N1LSS2fcIq8f8vqKK291tcZByv5vhCz9P9gbKzpOKEryo4sDBgywGJDXF3l9kddXTHlff91swgR9eUt3B68p8voir6/48nbokLHddqPj5IlSPQAAUk7nUNLhwL16meX8PioABNo+9O5d6hTxoVQPrqUt+sVu/aBh2ktbhLy+yOuLvL5iyytlZbW2227lUeSNrX3J64u8myZvz55x5I0VY3komn7lOibk9UVeX+T1FVve2GpEYmtf8voir6/Y8saIUj0AACIp1dtmG7P33y91GgBoO6L5AVzER0PBVVVV4ToG5PVFXl/k9RVbXqmrq40mb2ztS15f5PUVW95Y0XECAAAAgAIo1QMAIOUo1QMAH5TqwY2GgKdNmxbNUDB5fZHXF3l9xZZXamvjKtWLqX3J64u8vmLLGys6Tihaly5dLCbk9UVeX+T1FVveMv0mZ0Ria1/y+iKvr9jyxohSPQAAUo5SPQDwQake3NTV1dmUKVPCdQzI64u8vsjrK7a8Ulu7Npq8sbUveX2R11dseWNFxwlFKSsrsx49eoTrGJDXF3l9kddXbHmlvLw8mryxtS95fZHXV2x5Y0WpHgAAKUepHgD4oFQPbjQEPHHixGiGgsnri7y+yOsrtryydm1cpXoxtS95fZHXV2x5Y1VR6gCIi8pE+vXrF65jQF5f5PVFXm/llslsb3PnllvaIyejTB06KGvKw0a6PJDXF3l9xZY3VpTqAQDapS9+0WzUKIsKpXoAsHFRqgc3GgIeP358NEPB5PVFXl/k9TVlyrr9hp06Zayy0lJ/6dIlY4cc8l407Rvb8kBeX+T1FVveWFGqh6JoCHjgwIHRDAWT1xd5fZF303WgBg9O/5moGhoyNn9+WTTtG9vyQF5f5PUVW95YUaoHAGiX+vQxW7jQ7NVXzQYNKnUaAEApUKoHNxoCHjduXDRDweT1RV5f5PW2br9hLHlja1/y+iKvL/IiHzpOKIqGgAcNGhTNUDB5fZHXF3k3jVjyxta+5PVFXl/kRT6U6gEA2iVK9QAAyyjVg5fa2lobM2ZMuI4BeX2R1xd5va3bbxhL3tjal7y+yOuLvMiHEScUpaGhwaqrq6179+5RDAeT1xd5fcWWd8qUBrv66jqrr+9oZWXpP0vduHEZq60ts+nTG2zw4PS3b2zLA3l9kdcXeduPZUX0Deg4AQA2ijPPNLv3XouKvl/Mm2fWt2+pkwAA0t434HecUBQNAY8dO9aGDx9uHTt2tLQjry/y+oot75o1DaEC/KSTGuyYY9K/x1Nnn6qufsF69hxmZulv39iWB/L6Iq8v8iIfRpxQFC0uy5cvt27dukVRikNeX+T1FVveb34zYw88UGbXX5+x738//Xlja1/y+iKvL/L6ii1vmjDiBDdaGWPqcJLXF3l9xZbXbN2HdSwf2rG1L3l9kdcXeX3FljdW6a+lQOqGgkeOHBnNWVvI64u8vmLLq4OTpb6+3mIQW/uS1xd5fZHXV2x5Y0WpHoqixaWmpsYqKyuj2KtMXl/k9RVb3hhL9WJqX/L6Iq8v8vqKLW+a8DtOcFVREVeFJ3l9kddXbHljE1v7ktcXeX2R11dseWNExwlFn4Vq9OjR4ToG5PVFXl+x5W1oyDQp2Uu72NqXvL7I64u8vmLLGytK9VAULS5aKbVXI4ahYPL6Iq+v2PLGWKoXU/uS1xd5fZHXV2x504RSPbiKbW8GeX2R11dseWMTW/uS1xd5fZHXV2x5Y0THCUWvlPqBtVhWTvL6Iq+v2PLGWKoXU/uS1xd5fZHXV2x5Y0WpHgBgo/jmN80eeMDshhvMLrmk1GkAACiMUj24UT9bC1gs/W3y+iKvr9jymq3LGUve2NqXvL7I64u8vmLLGys6TiiKhoAnTJgQzVAweX2R11dseWMs1Yupfcnri7y+yOsrtryxolQPALBRUKoHAIgNpXpwoz3JS5YsiWaPMnl9kddXbHmT/XCZTBx5Y2tf8voiry/y+ootb6z4iWEUpb6+3qZMmWKHH364lZenv99NXl+x5f3d7zJ27bWbWYcOsfzGRZnV1m5uHTvGkXfZsqYle2kX2/JLXl/k9UVeX7HljRWlegDajSFDzP7971KnaNv0u4tPPGE2YkSpkwAAsHH7Bow4oSgaAl68eLH16tUrij0a5PUVW951Z30rs1tuabDPfjaO9l26dKn16NEjivZV3oaGJTZo0FZRVILHtvyS1xd5fZHXV2x5Y0XHCUWvmFVVVfaZz3wmihWTvL5iy5uMr/fvn7E997TUq6trsPHjp9vuu3/GKirKI8n7ijU0xLE8xLb8ktcXeX2R11dseWNFqR6AdleqN2qU2dFHlzoNAAAoNc6qB9c9GvPmzYvmrC3k9RVb3uQHWmPJG1v7ktcXeX2R1xd5fcWWN1Z0nFAUrZCzZ8+OZsUkr6/Y8ibj67EMtMfWvuT1RV5f5PVFXl+x5Y0VpXoA2g1K9QAAQDZK9eBGezLmzJkTzR4N8vqKLW+MpXoxtS95fZHXF3l9kddXbHljRccJbbqGlry+YssbY6leTO1LXl/k9UVeX+T1FVveWFGqB6DdoFQPAABko1QPburr623WrFnhOgbk9RVb3qRUL5a8sbUveX2R1xd5fZHXV2x5Y0XHCUXRAOXSpUujKXUir6/48lpU4mtf8noiry/y+iKvr9jyxopSPQDtBqV6AAAgG6V6cKMh4BkzZkQzFExeX7HljbFUL6b2Ja8v8voiry/y+ootb6zoOKFoq1evtpiQ11dseWMTW/uS1xd5fZHXF3l9xZY3RpTqAWg3KNUDAADZKNWDGw0BV1VVRTMUTF5fseWNsVQvpvYlry/y+iKvL/L6ii1vrOg4AQAAAEABlOoBaDco1QMAANko1YMbDQFPmzYtmqFg8vqKLW+MpXoxtS95fZHXF3l9kddXbHljRccJRevSpYvFhLy+Yssbm9jal7y+yOuLvL7I6yu2vDGiVA9Au0GpHgAAyEapHtzU1dXZlClTwnUMyOsrtrzJfqJYShlia1/y+iKvL/L6Iq+v2PLGio4TilJWVmY9evQI1zEgr6/48lpU4mtf8noiry/y+iKvr9jyxopSPQAbRIM2F11kNnOmRWPiRLPlyynVAwAAxfcNKlr8K5BDQ8CTJ0+2oUOHWkVF+hcf8vqZNs3s5pstSr16qZQh3e0b2/Ig5PVFXl/k9UVeX7HljRUti6KUl5dbv379wnUMyOuntnbddc+e9XbDDWVRZG5oaLDOnRfa/vv3thjEtDwIeX2R1xd5fZHXV2x5Y0WpHoANMmmS2UEHmQ0caDZrVqnTAAAAFI+z6sF1KHj8+PHRnLWFvP5qalZHkze29iWvL/L6Iq8v8voiL/Kh44SiaAh44MCB0QwFk9dfRUXHaPLG1r7k9UVeX+T1RV5f5EU+JW/dW265xQYMGGCVlZU2bNiwcGBbS2666Sbbbbfdwq8j9+/f3y6++GKrqanZZHnbu9hqaMnrTwehxpI3tvYlry/y+iKvL/L6Ii/yKWnrPvTQQ3bJJZfY1VdfbS+99JLtvffeNmLECFu4cGHe+//lL3+xyy+/PNz/jTfesLvvvjs8x49+9KNNnr290hDwuHHjohkKJq+/1atXRZM3tvYlry/y+iKvL/L6Ii9S13G68cYb7ZxzzrEzzjjD9txzT7vtttusa9euds899+S9/8SJE+3ggw+2r3/962GUavjw4fa1r32t4CgVNh7tyRg0aFA0ezTI669Tp07R5I2tfcnri7y+yOuLvL7Ii3xK1rpr1661qVOn2pFHHvlxmPLy8P9JOl1XHgcddFB4TNJRevvtt2306NF2dAu/ZLlmzZpwtozsi9Tr1zv/7zrftHrs2dM6jXFL07W1tU2mk5MVJtO65E5L9rQenz2d7DVoblr5sqc3xTzpPUp+mTqGeVLerbbaqnF5SPv7lORNcsWw7HXo0KHxNdO+Pql9e/Xq1fjcpV6fCs2T8m699daNr1Pq9anQPGm70Lt37/AcaVifCs2TKK9uS8v61NI86TWSvGlYnwrNk/6vvHqeNKxPheYpyZu9fKR5G6HnVl6td2lYnwrNkx6n7Zm2a2lYnwrNk671eaG8aVifCs2Tnr9nz56NeUu9PtWlaNlrzTylvuO0ePHi0AB9+vRpcrv+/8EHH+R9jEaafvrTn9ohhxxiHTt2DAfBffazn22xVO+6664LpxhMLjouSqqqqsK1Sv50kVdeecVmzpwZpqdNm2bvvPNOmFZHbe7cuY2jXvPnzw/TOnuJ5kM0PFpdXR2mx44da8uXLw/T6tjpGCy9mZrWtf6vadH9dH/R4/U8Sfvo+UWvp9cV5Ug6jsqnnKLcyu89T1rIlF3TMcxTkjdZptL+PinvE088YbP+7/zeaV723nrrrTC9atUqe/nll6NYn9S+Tz75pD3//POpWJ8KzVOSNy3rU6F5Wrp0qY0ZMyY161OheXrvvfdC3ueee67k61Nr5mnGjBkhr3YgpmF9KjRPyqm8yp2G9anQPGk5UF4tF2lYn1ozT8qr9S4N61Nr5knbs+RzudTrU6F50udEkjcN61OhedLnsL4/KG8a1qfxKVv2Cs1T6n/H6f333w8HsWmGDzzwwMbbf/jDH4aN14svvrjeY5599lk7+eST7Wc/+1k4kYS+XF544YWh3O/KK69sdsRJl4RGnNR5WrJkSRg5SXqv2muePa2G1F6cZFo9+KQXn29aC6rum0zrgPlkL1DyC866f/a0On9q/mQ62YuYTOui+zc3neyBTKbzzcfGnie9nhbK7L0aaZ4n+fDDD8N7redM+/ukLMrbvXv38PxpXvYmTszYZz5TYQMG1NtbbzWEx6Z9fdLzad3X7zSoxLDU61OhedL99aVo8803t86dO5d8fSo0T7rW72BsttlmoX1LvT4VmqfkM0Htm5zkJM3bCN1fH/7Kq9tKvT4VmiddVqxYYd26dQuvXer1qdA86XblTX7HpdTrU6F5UuXOypUrw05hZSz1+lRonvRdTO2bfPcq9fpUaJ7Uvto+qAok+fxI8zZCr6/tb5K31OtTQ4qWvULzpPdZ37ta8ztOJes4aYHU8UyPPPKInXDCCY23n3baaaGHOHLkyPUec+ihh9oBBxxgv/nNbxpve+CBB+zcc88NK2Nr6jr5AVxg4+AHcAEAQOyi+AFc7Y3cf//97Zlnnmm8Tb1Q/T97BCqbSoJyO0fqTUqJ+n/tjnrqo0aNaqxjTTvy+lu1amU0eWNrX/L6Iq8v8voiry/yIlUjTqJTiWuE6fbbb7ehQ4eG32h6+OGHQz20jnU69dRTQzmfjlOSa665JpyJ74477mgs1fvOd74TOmB6rtZgxOmT0eKiUhGVXmjIM+3I6z/itOOODTZ7dlnq88bWvkJeX+T1RV5f5PVF3vZjWRF9g3VFfiVy0kkn2aJFi+yqq64KB+/vs88+4UC85IQR7777bpMRpv/+7/8OC4Ou582bF87Ocuyxx9rPf/7zEs5F+6L2j6nDSV5/Wkdj2UbH1r7k9UVeX+T1RV5f5EU+JT/Z+/nnn29z5swJBw3qhBAaSco+GcR9993X+H8dwKUfv9VI0+rVq0PH6pZbbgkHdGHT0BCwjj+LZSiYvP5WrlwRTd7Y2pe8vsjri7y+yOuLvEhdqV4pUKr3yWhx0WkbKysroxgKJq9/qd5OOzXYrFnxlOrF0r5CXl/k9UVeX+T1Rd72Y1kMJ4dAvJLTOMaCvL5i20DH1r7k9UVeX+T1RV5f5EUuOk4oSvYPhsWAvP70OyKx5I2tfcnri7y+yOuLvL7Ii3wo1UNRkh87S35ILO3Iuyl+xylj+mHvtOeNrX2FvL7I64u8vsjri7ztxzJK9eAptr0Z5PUV276X2NqXvL7I64u8vsjri7zIRccJRa+UY8eOjWblJK8//TB1LHlja1/y+iKvL/L6Iq8v8iIfSvUAfMJSPbNZs0qdBgAAoHiU6sGN+tlawGLpb5PXX0NDQzR5Y2tf8voiry/y+iKvL/IiHzpOKIqGgCdMmBDNUDB5/dXUrI4mb2ztS15f5PVFXl/k9UVe5EOpHoANQqkeAACIHaV6cC3LWrJkSbiOAXn9NTTUR5M3tvYlry/y+iKvL/L6Ii/yoeOEotTX19uUKVPCdQzI66+mZk00eWNrX/L6Iq8v8voiry/yIh9K9QBsEEr1AABA7CjVgxsNAS9cuDCaoWDy+quvr4smb2ztS15f5PVFXl/k9UVe5EPHCUXRCllVVRXNiklef2vXro0mb2ztS15f5PVFXl/k9UVe5EOpHoANQqkeAABoT32Dik2WCm2C9mTMnz/fttlmGysvT/+AZWx5a2sb7KWXFlrv3r1Tn3f+/HXX+s2Ihoby1OeNcXkgry/y+iKvL/L6Ii/yoeOEolfM2bNnW58+faJYMWPLO3y42bPP9rWY1NXVWkNDxyjaN7blgby+yOuLvL7I64u8yIdSPSBFunUzW7HCrFMnsxi2e2VlZhdfbPbzn5c6CQAAQPE4qx5c92jMmTMnmoMPY8trtm4/xmuvNdjq1Zb6y4oVDXbuufG0b2zLA3l9kdcXeX2R1xd5kQ8dJxRFK+S8efOiWTFjy5uIJW9s7UteX+T1RV5f5PVFXl+x5Y0VpXpACkv1Zs8222mnUqcBAABo25ZRqgcv9fX1NmvWrHAdg9jyJqV6seSNrX3J64u8vsjri7y+yOsrtryxouOEomiAcunSpeE6BrHlTcSSN7b2Ja8v8voiry/y+iKvr9jyxopSPSBFKNUDAADYdCjVgxsNAc+YMSOaoeDY8sZYqhdT+5LXF3l9kdcXeX2R11dseWNFxwlFW63zUEcktryxia19yeuLvL7I64u8vsjrK7a8MaJUD0gRSvUAAAA2HUr14EZDwFVVVdEMBceWN8ZSvZjal7y+yOuLvL7I64u8vmLLGys6TgAAAABQAKV6QIpQqgcAALDpUKoHNxoCnjZtWjRDwbHljbFUL6b2Ja8v8voiry/y+iKvr9jyxoqOE4rWpUsXi0lseWMTW/uS1xd5fZHXF3l9kddXbHljRKkekCKU6gEAAGw6lOrBTV1dnU2ZMiVcxyC2vEmpXix5Y2tf8voiry/y+iKvL/L6ii1vrOg4oShlZWXWo0ePcB2D2PImYskbW/uS1xd5fZHXF3l9kddXbHljRakekCKU6gEAAGw6lOrBjYaAJ06cGM1QcGx5YyzVi6l9yeuLvL7I64u8vsjrK7a8saLjhKKUl5dbv379wnUMYsubiCVvbO1LXl/k9UVeX+T1RV5fseWNFaV6QIpQqgcAALDpUKoHNxoCHj9+fDRDwbHljbFUL6b2Ja8v8voiry/y+iKvr9jyxoqOE4qiIeCBAwdGMxQcW95ELHlja1/y+iKvL/L6Iq8v8vqKLW+sKNUDUoRSPQAAgE2HUj240RDwuHHjohkKji1vjKV6MbUveX2R1xd5fZHXF3l9xZY3VnScUBQNAQ8aNCiaoeDY8iZiyRtb+5LXF3l9kdcXeX2R11dseWNFqR6QIpTqAQAAbDqU6sFNbW2tjRkzJlzHILa8SaleLHlja1/y+iKvL/L6Iq8v8vqKLW+sGHFCURoaGqy6utq6d+8exXBwbHm7dcvYihVlNnNmg+28c/rzxta+5PVFXl/k9UVeX+T1FVveWPsGdJyAFKFUDwAAYNOhVA9uNAQ8atSoaIaCY8sbY6leTO1LXl/k9UVeX+T1RV5fseWNFSNOKIoWl+XLl1u3bt2srKzM0i62vEmp3qxZGRs4MP15Y2tf8voiry/y+iKvL/L6ii1vmlCq1wI6TkgzSvUAAAA2HUr14EZDwCNHjoxmKDi2vDGW6sXUvuT1RV5f5PVFXl/k9RVb3lgx4oSiaHGpqamxysrKKIaCY8sbY6leTO1LXl/k9UVeX+T1RV5fseVNE0ac4KqiosJiElve2MTWvuT1RV5f5PVFXl/k9RVb3hjRcUJR6urqbPTo0eE6BrHlTcSSN7b2Ja8v8voiry/y+iKvr9jyxopSPRRFi4tWSu3ViGEoOLa8MZbqxdS+5PVFXl/k9UVeX+T1FVveNKFUD65i25sRW97YxNa+5PVFXl/k9UVeX+T1FVveGNFxQtEr5dixY6NZOWPLm4glb2ztS15f5PVFXl/k9UVeX7HljRWlekCK8DtOAAAAmw6lenCjfrYWsFj627HlTX7HKZa8sbUveX2R1xd5fZHXF3l9xZY3VnScUBQNAU+YMCGaoeDY8iZiyRtb+5LXF3l9kdcXeX2R11dseWNFqR6QIpTqAQAAbDqU6sFNQ0ODLVmyJFzHILa8SaleLHlja1/y+iKvL/L6Iq8v8vqKLW+s6DihKPX19TZlypRwHYPY8iZiyRtb+5LXF3l9kdcXeX2R11dseWNFqR6QIpTqAQAAbDqU6sGNhoAXLlwYzVBwbHljLNWLqX3J64u8vsjri7y+yOsrtryxouOEomiFrKqqimbFjC1vIpa8sbUveX2R1xd5fZHXF3l9xZY3VpTqASlCqR4AAMCmQ6ke3GhPxrx586LZoxFb3hhL9WJqX/L6Iq8v8voiry/y+ootb6zoOKEoWiFnz54dzYoZW95ELHlja1/y+iKvL/L6Iq8v8vqKLW+sKNUDUoRSPQAAgE2HUj240Z6MOXPmRLNHI7a8MZbqxdS+5PVFXl/k9UVeX+T1FVveWJW843TLLbfYgAEDrLKy0oYNG2aTJ09u8f7V1dV23nnn2TbbbGOdO3e2XXfd1UaPHr3J8rZ3sdXQxpY3EUve2NqXvL7I64u8vsjri7y+Yssbq5KW6j300EN26qmn2m233RY6TTfddJP97W9/szfffNN69+693v3Xrl1rBx98cPjbj370I+vXr1/oXXfv3t323nvvVr0mpXpIM0r1AAAANp1i+gYVVkI33nijnXPOOXbGGWeE/6sDNWrUKLvnnnvs8ssvX+/+un3JkiU2ceJE69ixY7hNo1XYdOrr6+2dd96xHXfc0Tp06GBppl0CEybU26uvLrQ+fXpbeXm680ptrfZjlIV2Nkt/3piWByGvL/L6Iq8v8voir6/Y8saqZB0njR5NnTrVrrjiisbbysvL7cgjj7RJkyblfcxjjz1mBx54YCjVGzlypG299db29a9/3S677LJmF5I1a9aES3avUtZ9Mf34Wo/Pnq6rq7OysrLGaWXTpbnp2tracN9kuqKiIjw+mRbdP3tanT8N+CXTGl5VhmRaF92/uWndV49PpvPNx8aeJ73ehx9+aDvssEPj/dM6T08/XW5HHaXlYhuLR9n/XdeFjlPalz3Rzoz+/fu7L3sbY56SvNttt114/lKvT4XmSX9XXo2ud+nSJfXbCFm6dGnIq/Jrz2VvY8yTnlN5tTzoPqVenwrNk55febW+Zc9fWrcR+fKmeRuhxyrv9ttv3/h+pHkboe9RyqvPY2Us9fpUaJ6UV9sz7fBOw/pUaJ6UMcmbhvWp0DzpubPzlnp9akjRsldonoopvivZMU6LFy8ODdKnT58mt+v/H3zwQd7HvP322/bII4+Ex+m4piuvvNJuuOEG+9nPftbs61x33XVh+C25JBtw/bqyvPHGG+Eir7zyis2cOTNMT5s2LfTcRcddzZ07N0xrtGv+/Plhevz48WE+ZNy4ceH4Kxk7dqwtX748TCtnTU1NeHM0rWv9PzkuS/fT/UWP1/Mk7aPnF72eXleUIzkOTPmUU5Rb+b3nSQuY3h/NR9rn6fXXPwrTm29ea0OHqszTbM89l9qwYbX/N73EDjigLkzvsceHduCB9XbQQQ1hWtf6v6b1d91P99e0Hq/n0bSed9Cg6jA9ZMga22uvj8L0pz9dY4MHLwvT+++/2vbee3mY3m+/VbbPPivC9L77rgwXTes2/U3TJ520QK0SxbKn5UHD2q+//noU65Py6kt9Mh+lXp8KzZPy7r777o3zkfZtxOrVq23IkCFhOobt3qJFi0JeTadhfSo0T7qf8ibTpV6fCs2TcipvMl3q9anQPOm1lFfLRRrWp0LzpPsor9a7NKxPheZJ09qeabuWhvWp0DxpWp8XypuG9anQPOlzWJ/HypuG9Wl8ipa91sxTq2VKZN68eereZSZOnNjk9ksvvTQzdOjQvI/ZZZddMv3798/U1dU13nbDDTdk+vbt2+zr1NTUZD766KPGy9y5c8PrLlmyJPxdz5U8X/Z0bW1tk+n6+voWp9euXdtkuqGhocm0LrnTkj2tx2dP6/lbmla+7Ol887Gx50mPr6qqCvdL+zzddVd9Rkv4YYcta8yQ9vdJl9deey2zZs2aKJa93LxpX5+U4fXXXw/bhTSsT4XmKcm7evXqkq9PrZkn/e2NN94IedOwPhWaJ72u8mp5SMP6VGietJ4pr67TsD4Vmqd8edO8jdByoLzJslzq9anQPGk9U95kWS71+lRonpRX2zM9Pg3rU6F50vKQ5E3D+lRonrSe6fM4+/O5lOtTbYqWvULzVF1dHfoG6icUUrJSvV69eoUhtAULtHf9Y/p/37598z5GZ9LTEGB2Wd4ee+wRRkA0BNypU6f1HqMz7+mSK3mO7OfKnk6G8Vo7nRxzVey0hg2T6WR4sbXTzWX3nie1tXJnZ0/jPP3fU1p9fUPj86f9fdJoqkpLk/+nfdnLzZv29Ul5tWcpea00rE8tTSd507A+tWaelFd7vzWt96K5+UvLNiLJu257Ecc2Qnk1Hcs2IjdvmrcRen7l1XSSLe3bCOXVPOXbpqVxG5Hs2U/L+lRoPpK8aVmfCk0nh6akYX0qT9myV2g+Ul+qp07O/vvvb88880zjbap71P91HFM+OqPerFmzmpxq8a233godqnydJmx8Woj33XffJgtz2qlEM5a8sbUveX2R1xd5fZHXF3l9kRep+x2nSy65xO688067//77Qz3jd77zHVu5cmXjWfZ0qvLsk0fo7zrw7cILLwwdJp2B7xe/+EU4WQQ2De2h1fFhyQF7MVi+fFk0eWNrX/L6Iq8v8voiry/y+iIvUnc68pNOOikchHnVVVeFcrt99tnHnnzyycYTRrz77ruNQ3qiEzuMGTPGLr74Yhs8eHA4aE+dKJ1VDwAAAADa5A/glgI/gNt+3Huv2Zlnmh1zjNnjj5c6DQAAAGLuG5S0VA/x0RCwTgkZ01CwVoRY8sbWvuT1RV5f5PVFXl/k9UVe5EPHCUXTD3HGJLYDJWNrX/L6Iq8v8voiry/y+iIvclGqhzaLUj0AAAC0hFI9uNEvLE+ZMiVcx0K/HB1L3tjal7y+yOuLvL7I64u8vsiLfOg4oSj6kbAePXoU9WNhpZb9Y5xpF1v7ktcXeX2R1xd5fZHXF3mRD6V6aLMo1QMAAEBLKNWDGw0BT5w4Maqh4KVLl0STN7b2Ja8v8voiry/y+iKvL/IiHzpOKIp+kFg/PJz9w8RpV1lZGU3e2NqXvL7I64u8vsjri7y+yIt8KNVDm0WpHgAAAFpCqR7caAh4/PjxUQ0FL1nyYTR5Y2tf8voiry/y+iKvL/L6Ii/yoeOEomgIeODAgVENBXftulk0eWNrX/L6Iq8v8voiry/y+iIv8qFUD20WpXoAAABoCaV6cKMh4HHjxkU1FLx48eJo8sbWvuT1RV5f5PVFXl/k9UVe5EPHCUXREPCgQYOiGgru1q1bNHlja1/y+iKvL/L6Iq8v8voiL/KpKHUAxEUrZO/evS0mnTt3tli2I7G1L3l9kdcXeX2R1xd5fZEX+UTydRJpUVtba2PGjAnXsVi0aGE0eWNrX/L6Iq8v8voiry/y+iIv8qHjhKJ06NDBhgwZEq5j0b1792jyxta+5PVFXl/k9UVeX+T1RV7kQ6keih4K3mqrrSwmHTt2iqpUL6b2Ja8v8voiry/y+iKvL/Iin0i+TiItNAQ8atSoqIaCFy5cEE3e2NqXvL7I64u8vsjri7y+yIt8+B0nFEWLy/Lly8OZ6srKyiyG33EaMaLWnniiIvV5Y2tfIa8v8voiry/y+iKvL/K2H8uK6BtQqoeiaGWMrcNZUdHRYtmGxNa+5PVFXl/k9UVeX+T1RV7kQ6keiqIh4JEjR0Y1FLxgwQfR5I2tfcnri7y+yOuLvL7I64u8yIdSPRRFi0tNTY1VVlZGU6r3hS/U2+jR5anPG1v7Cnl9kdcXeX2R1xd5fZG3/VhWRN+AEScUraIirgrP8vK4NiCxtS95fZHXF3l9kdcXeX2RF7noOKEodXV1Nnr06HAdi4ULF0aTN7b2Ja8v8voiry/y+iKvL/IiH0r1UBQtLloptVcjllK9o49usMcfL0t93tjaV8jri7y+yOuLvL7I64u87ccy77PqvfHGG/bggw/ahAkTbM6cObZq1Srbeuutbd9997URI0bYf/3Xf1nnzp03ND9SLlkxY9HQoH0D8WxEYmtf8voiry/y+iKvL/L6Ii8+UaneSy+9ZEceeWToIP3rX/+yYcOG2UUXXWTXXnutfeMb3wi93R//+Me27bbb2q9+9Stbs2ZNMU+PSFbKsWPHRjUUvHjxomjyxta+5PVFXl/k9UVeX+T1RV584lK9HXfc0S699FL7+te/bt27d2/2fpMmTbLf/e53NnjwYPvRj35kaUKpXvuRlOodc4zZ44+XOg0AAADaTaneW2+9ZR07dix4vwMPPDBcOJd82xPjL1PX1dVaJhNHzW9s7UteX+T1RV5f5PVFXl/kxScu1WtNp0l0zFMx90c8NASsY9tiGgpesmRJNHlja1/y+iKvL/L6Iq8v8voiLzbqWfWOOOII++Mf/2j9+vVrcvvkyZPD8U4anUojSvXaD0r1AAAAUPIfwNUvE+sYpoceeij8v6Ghwa655ho75JBD7Oijj97Qp0XK6X3WCI6uY1FbuzaavLG1L3l9kdcXeX2R1xd5fZEXG7XjNGrUKPvpT39qZ555ZjhZhDpMd955pz3++ON20003bejTIuXq6+ttypQp4ToW1dXV0eSNrX3J64u8vsjri7y+yOuLvHD5AdwrrrginHpc541/9tln7aCDDrI0o1Sv/aBUDwAAACUv1Vu6dGn4odtbb73Vbr/9djvxxBNt+PDh9oc//GFDnxIR0BDwwoULoxoK1u+JxZI3tvYlry/y+iKvL/L6Iq8v8mKjdpwGDRpkCxYssGnTptk555xjDzzwgN1999125ZVX2jHaxY82SStkVVVVVCumTs8ZS97Y2pe8vsjri7y+yOuLvL7Ii41aqnfttdfaj3/8Yysvb9r3eu+99+yMM86wp556ytKIUr32g1I9AAAAlLxUTyNLuZ0m2W677VLbacInpz0Z8+bNi2qPRk1NTTR5Y2tf8voiry/y+iKvL/L6Ii8+ccfp3XffLebu4Q1E26IVcvbs2VGtmKtWrYwmb2ztS15f5PVFXl/k9UVeX+TFJy7V69Onj51wwgl29tln25AhQ/LeR8NcDz/8sP3ud7+zc8891773ve9ZmlCq135QqgcAAICSlOq9/vrrttlmm9nnP/9569u3bzgJhE4MccEFF9g3vvEN22+//ax37952zz332K9//evUdZrwyWlPxpw5c6Lao7F69apo8sbWvuT1RV5f5PVFXl/k9UVefOKOU8+ePe3GG2+0+fPn280332y77LKLLV682GbOnBn+fsopp9jUqVNt0qRJdvTRRxfz1IhEjDW0HOPkh7y+yOuLvL7I64u8vsiLjXJWvbffftt23HFHKysrsxhRqtd+UKoHAACAkp1VT6NMixYtavz/SSedFH7PCe1DfX29zZo1K1zHYuXKldHkja19yeuLvL7I64u8vsjri7zYKB2n3AGq0aNHhy+maB/0/i9dunS95SDNamtro8kbW/uS1xd5fZHXF3l9kdcXebFRSvX0200ffPBBOAmEdOvWzaZPn2477bSTxYBSvfaDUj0AAACUrFRPxzblHt8U6/FOKJ6GgGfMmBHVUPCKFSuiyRtb+5LXF3l9kdcXeX2R1xd5kU+FFUkDVKeffrp17ty58Yxl3/72t8NpyrP9/e9/L/apEYnVq1dbTGLbiMTWvuT1RV5f5PVFXl/k9UVefOJSvTPOOKNV97tXdVIpRKle+0GpHgAAADZW36DoEae0doiw6UZv3njjDdtjjz2sQ4cOFoPly5dZff1mUeSNrX3J64u8vsjri7y+yOuLvNgoxzgBAAAAQHtTdKle7CjVaz8o1QMAAEDJzqqH9k1DwdOmTYvqhAtaEWLJG1v7ktcXeX2R1xd5fZHXF3mRDx0nFK1Lly4Wk9hqfWNrX/L6Iq8v8voiry/y+iIvclGqhzaLUj0AAAC0hFI9uKmrq7MpU6aE61hUV1dHkze29iWvL/L6Iq8v8voiry/yIh86TihKWVmZ9ejRI1zHomPHjtHkja19yeuLvL7I64u8vsjri7zIh1I9tFmU6gEAAKAllOrBjYaAJ06cGNVQ8NKlS6LJG1v7ktcXeX2R1xd5fZHXF3mRDx0nFKW8vNz69esXrmNRWVkZTd7Y2pe8vsjri7y+yOuLvL7Ii3wo1UObRakeAAAAWkKpHtxoCHj8+PFRDQUvWfJhNHlja1/y+iKvL/L6Iq8v8voiL/Kh44SiaAh44MCBUQ0Fd+26WTR5Y2tf8voiry/y+iKvL/L6Ii/yoVQPbRalegAAAGgJpXpwoyHgcePGRTUUvHjx4mjyxta+5PVFXl/k9UVeX+T1RV7kQ8cJRdEQ8KBBg6IaCu7WrVs0eWNrX/L6Iq8v8voiry/y+iIv8qnIeyvQDK2QvXv3tph07tzZYtmOxNa+5PVFXl/k9UVeX+T1RV7kE8nXSaRFbW2tjRkzJlzHYtGihdHkja19yeuLvL7I64u8vsjri7zIh44TitKhQwcbMmRIuI5F9+7do8kbW/uS1xd5fZHXF3l9kdcXeZEPpXooeih4q622sph07NgpqlK9mNqXvL7I64u8vsjri7y+yIt8Ivk6ibTQEPCoUaOiGgpeuHBBNHlja1/y+iKvL/L6Iq8v8voiL/Lhd5xQFC0uy5cvD2eqKysrsxh+x2nEiFp74omK1OeNrX2FvL7I64u8vsjri7y+yNt+LIvtd5xuueUWGzBggFVWVtqwYcNs8uTJrXrcgw8+GBaOE044wT0j1lF7a6GKaaWsqOgYTd7Y2pe8vsjri7y+yOuLvL7Ii1R2nB566CG75JJL7Oqrr7aXXnrJ9t57bxsxYoQtXLiwxcf95z//sR/84Ad26KGHbrKsWDcUPHLkyKiGghcs+CCavLG1L3l9kdcXeX2R1xd5fZEXqSzV0wiTzgJy8803h/83NDRY//797YILLrDLL78872Pq6+vtM5/5jJ155pk2YcIEq66utn/84x+tej1K9T4ZLS41NTVhdDCWUr0vfKHeRo8uT33e2NpXyOuLvL7I64u8vsjri7ztx7JYSvXWrl1rU6dOtSOPPPLjQOXl4f+TJk1q9nE//elPw498nXXWWQVfY82aNaFBsi9J5yu5zjddV1fXZFodupam1cPPnk76o8m0LrnTkj2tx2dP6/lbmla+7OlNNU/JPMQyT8lrxPI+aYMX07KXnTeG9UnbmDStT4XmSaeWTdP6VGieKioqUrU+FZon5U3T+lRonpK8aVmfCs2T8qZpfSo0T8qbpvWp0Dwpb5rWp0LzlJwqOy3rU6F50udFmtanQvOUdJjSsj6ladkrNE+tVdKO0+LFi0Mj9OnTp8nt+v8HH3yQ9zH/+te/7O6777Y777yzVa9x3XXXhV5kctFollRVVYXrN954I1zklVdesZkzZ4bpadOm2TvvvBOmdczV3Llzw/TEiRNt/vz5YXr8+PFhHmTcuHFh5EvGjh0bDtCT0aNHhz0AejM1rWv9X9Oi++n+osfreZK20fOLXk+vK8qRHAOmfMopyq383vOk/E899ZStWLEi9fOkPQfrXneRLViwIIr3Sa+rH7CbPXt2FMteknf69OlRrE+6/5NPPpma9anQPGn6iSeeCOtcqden1syTppVbmdOwPhWap3nz5oVMaVmfCs3Tm2++GV5PZe1pWJ8KzZNy6rmUOw3rU2vmSffTcpGG9anQPGk90//1mDSsT4XmSdsxZU7mo9TrU2vmSZ8XypaG9anQPOlzWJ/Het20rE9pWfZaM09RlOq9//771q9fvzDTBx54YOPtP/zhD+25556zF198scn9NdODBw+2P/zhD3bUUUeF204//fQWS/U04qRLQiNO6jwtWbLEevTo0dh71V6Q7Gk1pHruybT2OiR7qvNNJ3tSkmntBdLjk2lJ9mYl0x07dgy93GRaPWNlSKaTvaHNTeu+yR7e5uZjY8+TrF69unEoOM3zdP/95Xb22eWhVO/xxz/ec5/m90m5tLzqvsme5TQve3rd7LxpX5/0nBrp1vPmzmsp1qdC85SM3kinTp1Sv41Ink+vqYylXp8KzZNuT/ZI5s5HGrcR2XtGla/U61OheUraVq+VvAdp3kZkb4eT9yPN2wjdnt2upV6fCs2Ttr2ibPm2dWnbRiQjFtr26v6lXp8KzZOeW5fOnTs3jj6lfRuRScn6pL5B9+7dW1WqV9KOk1airl272iOPPNLkzHinnXZa6AzpILdsL7/8su27775NfhU52TCrYbRXa+DAgS2+Jsc4tZ8aWo5x8kdeX+T1RV5f5PVFXl/kbT+WxXKMk3rx+++/vz3zzDNNOkL6f/YIVGL33Xe3V199NXSgkstxxx1nn/vc58J0UoYHP+qda8hT17FQqV4seWNrX/L6Iq8v8voiry/y+iIvUnlWPZ2OXCNMt99+uw0dOtRuuukme/jhh23GjBnhWKdTTz01lPPpWKV8CpXq5WLEqf1IRpyOOcZCqR4AAAAQ5YiTnHTSSXb99dfbVVddZfvss08YOdLBeMkJI959993Gg75QeupnawErcX+7KHV1H59NJe1ia1/y+iKvL/L6Iq8v8voiL1LZcZLzzz/f5syZEw4y1wkh9NtOiWeffdbuu+++Zh+rv7V2tAmfnIaA9dtZMQ0F60QgseSNrX3J64u8vsjri7y+yOuLvEhlqd6mRqnehtN5OI4/3uyFFywKOrvkihWU6gEAAOCT9w3WnY8PaAVVTMbYAdl559XW0NC58Yfs0kwnR9ExezotJnk3PvL6Iq8v8voiry/y+ootb6zoOKFoHTo02NSp687Vn3ZlZbU2Z854q68/PIoNiX6zYMqUKXb44eT1QF5f5PVFXl/k9UVeX7HljRWlemi1efPMtttOP16n3+AqdRoAAACgHZ1VDzHKNP7wcNop58KFC8nrhLy+yOuLvL7I64u8vsiLfOg4oWgapIxlxVTOqqoq8johry/y+iKvL/L6Iq8v8iIfSvXQapTqAQAAoC2hVA/O4hpxmjdvHnmdkNcXeX2R1xd5fZHXF3mRDx0ntPlSvdmzZ5PXCXl9kdcXeX2R1xd5fZEX+VCqh1ajVA8AAABtCaV6cBbXiNOcOXPI64S8vsjri7y+yOuLvL7Ii3zoOKHNl+rFVPNLXl/k9UVeX+T1RV5f5PUVW95YUaqHVqNUDwAAAG0JpXpwlrH6+nqLgXLOmjWLvE7I64u8vsjri7y+yOuLvMiHjhOKpkHKWAYqlXPp0qXkdUJeX+T1RV5f5PVFXl/kRT6U6qHVKNUDAABAW0KpHpzFVao3Y8YM8johry/y+iKvL/L6Iq8v8iIfOk4oWmyDlKtXr7aYkNcXeX2R1xd5fZHXF3l9xZY3RpTqodUo1QMAAEBbQqkeXGUyDdEMBStnVVUVeZ2Q1xd5fZHXF3l9kdcXeZEPHScAAAAAKIBSPbQapXoAAABoSyjVg6vYSvWmTZtGXifk9UVeX+T1RV5f5PVFXuRDxwlFKysrs5h06dLFYkJeX+T1RV5f5PVFXl/k9RVb3hhRqodWo1QPAAAAbQmlenAv1aurq7MYKOeUKVPI64S8vsjri7y+yOuLvL7Ii3zoOGGDSvViKddTzh49epDXCXl9kdcXeX2R1xd5fZEX+VCqh1ajVA8AAABtCaV6cBVbqd7EiRPJ64S8vsjri7y+yOuLvL7Ii3zoOKFoGgYuL49j0VHOfv36kdcJeX2R1xd5fZHXF3l9kRf5UKqHVqNUDwAAAG0JpXpwFVup3vjx48nrhLy+yOuLvL7I64u8vsiLfOg4oc2X6g0cOJC8Tsjri7y+yOuLvL7I64u8yIdSPbQapXoAAABoSyjVg6vYSvXGjRtHXifk9UVeX+T1RV5f5PVFXuRDxwltvlRv0KBB5HVCXl/k9UVeX+T1RV5f5EU+lOqh1SjVAwAAQFtCqR7cS/Vqa2stBso5ZswY8johry/y+iKvL/L6Iq8v8iIfRpywASNOGaupyUQxHNzQ0GDV1dXWvXt38jogry/y+iKvL/L6Iq8v8rYfy4roG9BxQqtRqgcAAIC2hFI9uIqtVG/UqFHkdUJeX+T1RV5f5PVFXl/kRT6MOGGDSvXWrFl3dr200+K9fPly69atG3kdkNcXeX2R1xd5fZHXF3nbj2WU6jWPjtOGo1QPAAAAbQmlenAVW6neyJEjyeuEvL7I64u8vsjri7y+yIt8GHFCmy/Vq6mpscrKSvI6IK8v8voiry/y+iKvL/K2H8sYcQI+VlFRYTEhry/y+iKvL/L6Iq8v8vqKLW+M6Dhhg/Zq1NXVWQyUc/To0eR1Ql5f5PVFXl/k9UVeX+RFPpTqoc2X6mkjor0w5N34yOuLvL7I64u8vsjri7ztxzJK9YCPxbb3hby+yOuLvL7I64u8vsjrK7a8MaLjhDZfqjd27FjyOiGvL/L6Iq8v8voiry/yIh9K9dBq/I4TAAAA2hJK9eAsE0adYqCcWiHI64O8vsjri7y+yOuLvL7Ii3zoOKHNl+pNmDCBvE7I64u8vsjri7y+yOuLvMiHUj20GqV6AAAAaEso1YOzjDU0NFgMlHPJkiXkdUJeX+T1RV5f5PVFXl/kRT50nFA0DVLW19dbDJRzypQp5HVCXl/k9UVeX+T1RV5f5EU+lOqh1SjVAwAAQFtCqR6cxVWqt3DhQvI6Ia8v8voiry/y+iKvL/IiHzpOKJoGKWNZMZWzqqqKvE7I64u8vsjri7y+yOuLvMiHUj20GqV6AAAAaEso1YOzuEac5s2bR14n5PVFXl/k9UVeX+T1RV7kQ8cJbb5Ub/bs2eR1Ql5f5PVFXl/k9UVeX+RFPpTqodUo1QMAAEBbQqkenMU14jRnzhzyOiGvL/L6Iq8v8voiry/yIh86TmjzpXox1fyS1xd5fZHXF3l9kdcXeX3FljdWlOqh1SjVAwAAQFtCqR6cZay+vt5ioJyzZs0irxPy+iKvL/L6Iq8v8voiL/Kh44SiaZAyloFK5Vy6dCl5nZDXF3l9kdcXeX2R1xd5kQ+lemg1SvUAAADQllCqB2dxlerNmDGDvE7I64u8vsjri7y+yOuLvMiHjhOKFtsg5erVqy0m5PVFXl/k9UVeX+T1RV5fseWNEaV6aDVK9QAAANCWUKoHV5lMQzRDwcpZVVVFXifk9UVeX+T1RV5f5PVFXqS243TLLbfYgAEDrLKy0oYNG2aTJ09u9r533nmnHXroodajR49wOfLII1u8PwAAAABEX6r30EMP2amnnmq33XZb6DTddNNN9re//c3efPNN692793r3P+WUU+zggw+2gw46KHS0fvWrX9mjjz5qr732mvXr16/g61Gqt+Eo1QMAAEBbElWp3o033mjnnHOOnXHGGbbnnnuGDlTXrl3tnnvuyXv/P//5z/bd737X9tlnH9t9993trrvusoaGBnvmmWc2efb2KrZSvWnTppHXCXl9kdcXeX2R1xd5fZEXqes4rV271qZOnRrK7RoDlZeH/0+aNKlVz7Fq1Sqrra21rbbaKu/f16xZE3qS2RdJFixd55uuq6trMq3OWUvTypA9nQzkJdO65E5L9rQenz2t529pWvmypzfVPCW5Y5mnzp07R/U+KW9My1523hjWJ41Up2l9KjRPypum9anQPHXp0iVV61OheVLeNK1PheZJeXPnL63biHx5S70+FZon5U3T+lRonpQ3TetToXnS9ixN61OheUrypmF9as086fM4TetTmpa9QvMURcdp8eLFoRH69OnT5Hb9/4MPPmjVc1x22WW27bbbNul8ZbvuuuvC8Fty6d+/f7hdB9DJG2+8ES7yyiuv2MyZM8O0eu3vvPNOmNYxVHPnzg3TEydOtPnz54fp8ePHh3mQcePGWXV1dZgeO3asLV++PEyPHj3aampqwpupaV3r/5oW3U/3Fz1ez5O0jZ5f9Hp6XVGO5Jgu5VNOUW7l95ynFStWNLarFrYY5qlDhw5hWdKvacfwPinvypUr7d13341i2VNe7exQqWwM65Pyahj+xRdfLPn61Jp5Ul6VICcj6qVenwrNk3ZkqRJgzJgxqVifCs3TwoULQ97nn38+FetToXl6++23Q97p06enYn0qNE/KqbzKnYb1qdA8aTlQXi0XaVifCs2T1jPl1XqXhvWp0DxpO6btmbZraVifCs2TPif0eaG8aVifCs2TPof1eay8aVifxqdo2WvNPLVapoTmzZunLl5m4sSJTW6/9NJLM0OHDi34+Ouuuy7To0ePzPTp05u9T01NTeajjz5qvMydOze85pIlS8Lf6+rqwiV3ura2tsl0fX19i9Nr165tMt3Q0NBkWpfcacme1uOzp/X8LU0rX/Z0vvnYmPM0d25DRktMRcW6nDHMk65feOGFzJo1a6J4n5K8Wm5jWPZ0efHFFxvzpn19SvKuXr265OtTa+Ypybtq1apUrE+F5kmXyZMnh7xpWJ8KzZO2C8qr5SEN61OhedJ6pry6TsP6VGie8uVN8zZCy4HyarlIw/pUaJ60nilvsu6Ven0qNE/Kq+2Z/p6G9anQPGl5SPKmYX0qNE9az7Lzlnp9qk3Rsldonqqrq0PfQP2EQiqshHr16hV6xgsWLGhyu/7ft2/fFh97/fXX2y9/+Ut7+umnbfDgwc3eT8OWydBlNr1u9nXudEVFRVHTHXXGhA2YLisra5zWngJdWjvdXHaveSors8bMypCdPa3zpOHXnj17Nj4m7e+TRmCVN8mf9mVPeVUmm+RJ+/qU5M13n029PrVmOsnbqVOnZueppelNPU/Kq7OdKq/ei+bmLy3biCSvsiTPn/ZthPIqe/L/tG8jcvOmeRuh25VXuZJsad5GaD1T3mSUIe3bCOXV9kzvQaH5TsM2Qq+f5M2+f1q3EcqeL2/atxFpWJ+Sz6sozqqnM+kNHTrUfv/734f/qx5x++23t/PPP98uv/zyvI/59a9/bT//+c/DMPUBBxxQ1OtxVr0Nx1n1AAAA0JZEdVa9Sy65JPw20/333x9qGr/zne+EYzx0lj3RqcqvuOKKxvvr9ONXXnllOOuefvtJx6/okn38DfzPqpcc0Jd2yqn6V/L6IK8v8voiry/y+iKvL/Iin5KW6slJJ51kixYtsquuuip0gHSa8SeffLLxhBE6SD4Z1pNbb701nI3vK1/5SpPnufrqq+2aa67Z5Pnbo6RULwbKqYNRyeuDvL7I64u8vsjri7y+yItUluptapTqbThK9QAAANCWRFWqh/jEVqqn01OS1wd5fZHXF3l9kdcXeX2RF/nQcUKbL9UbOHAgeZ2Q1xd5fZHXF3l9kdcXeZEPpXpoNUr1AAAA0JZQqgdXsZXq6dejyeuDvL7I64u8vsjri7y+yIt86DihzZfqDRo0iLxOyOuLvL7I64u8vsjri7zIh1I9tBqlegAAAGhLKNWDe6lebW2txUA5x4wZQ14n5PVFXl/k9UVeX+T1RV7kw4gTNmDEKWM1NZkohoMbGhqsurraunfvTl4H5PVFXl/k9UVeX+T1Rd72Y1kRfQM6Tmg1SvUAAADQllCqB1exleqNGjWKvE7I64u8vsjri7y+yOuLvMiHESdsUKnemjXrzq6Xdlq8ly9fbt26dSOvA/L6Iq8v8voiry/y+iJv+7GMUr3m0XHacJTqAQAAoC2hVA+uYivVGzlyJHmdkNcXeX2R1xd5fZHXF3mRDyNOaPOlejU1NVZZWUleB+T1RV5f5PVFXl/k9UXe9mMZI07AxyoqKiwm5PVFXl/k9UVeX+T1RV5fseWNER0nbNBejbq6OouBco4ePZq8Tsjri7y+yOuLvL7I64u8yIdSPbT5Uj1tRLQXhrwbH3l9kdcXeX2R1xd5fZG3/VhGqR7wsdj2vpDXF3l9kdcXeX2R1xd5fcWWN0Z0nNDmS/XGjh1LXifk9UVeX+T1RV5f5PVFXuRDqR5ajd9xAgAAQFtCqR6cZcKoUwyUUysEeX2Q1xd5fZHXF3l9kdcXeZEPHSe0+VK9CRMmkNcJeX2R1xd5fZHXF3l9kRf5UKqHVqNUDwAAAG0JpXpwlrGGhgaLgXIuWbKEvE7I64u8vsjri7y+yOuLvMiHjhOKpkHK+vp6i4FyTpkyhbxOyOuLvL7I64u8vsjri7zIh1I9tBqlegAAAGhLKNWDs7hK9RYuXEheJ+T1RV5f5PVFXl/k9UVe5EPHCUXTIGUsK6ZyVlVVkdcJeX2R1xd5fZHXF3l9kRf5UKqHVqNUDwAAAG0JpXpwFteI07x588jrhLy+yOuLvL7I64u8vsiLfOg4oc2X6s2ePZu8Tsjri7y+yOuLvL7I64u8yIdSPbQapXoAAABoSyjVg7O4RpzmzJlDXifk9UVeX+T1RV5f5PVFXuRDxwltvlQvpppf8voiry/y+iKvL/L6Iq+v2PLGilI9tBqlegAAAGhLKNWDs4zV19dbDJRz1qxZ5HVCXl/k9UVeX+T1RV5f5EU+dJxQNA1SxjJQqZxLly4lrxPy+iKvL/L6Iq8v8voiL/KhVA+tRqkeAAAA2hJK9eAsrlK9GTNmkNcJeX2R1xd5fZHXF3l9kRf50HFC0WIbpFy9erXFhLy+yOuLvL7I64u8vsjrK7a8MaJUD61GqR4AAADaEkr14CqTaYhmKFg5q6qqyOuEvL7I64u8vsjri7y+yIt86DgBAAAAQAGU6qHVKNUDAABAW0KpHlzFVqo3bdo08johry/y+iKvL/L6Iq8v8iIfOk4oWllZmcWkS5cuFhPy+iKvL/L6Iq8v8voir6/Y8saIUj20GqV6AAAAaEso1YN7qV5dXZ3FQDmnTJlCXifk9UVeX+T1RV5f5PVFXuRDxwkbVKoXS7mecvbo0YO8Tsjri7y+yOuLvL7I64u8yIdSPbQapXoAAABoSyjVg6vYSvUmTpxIXifk9UVeX+T1RV5f5PVFXuRDxwlF0zBweXkci45y9uvXj7xOyOuLvL7I64u8vsjri7zIh1I9tBqlegAAAGhLKNWDq9hK9caPH09eJ+T1RV5f5PVFXl/k9UVe5EPHCW2+VG/gwIHkdUJeX+T1RV5f5PVFXl/kRT6U6qHVKNUDAABAW0KpHlzFVqo3btw48johry/y+iKvL/L6Iq8v8iIfOk5o86V6gwYNIq8T8voiry/y+iKvL/L6Ii/yoVQPrUapHgAAANoSSvXgXqpXW1trMVDOMWPGkNcJeX2R1xd5fZHXF3l9kRf5MOKEDRhxylhNTSaK4eCGhgarrq627t27k9cBeX2R1xd5fZHXF3l9kbf9WFZE34COE1qNUj0AAAC0JZTqwVVspXqjRo0irxPy+iKvL/L6Iq8v8voiL/JhxAkbVKq3Zs26s+ulnRbv5cuXW7du3cjrgLy+yOuLvL7I64u8vsjbfiyjVK95dJw2HKV6AAAAaEso1YOr2Er1Ro4cSV4n5PVFXl/k9UVeX+T1RV7kw4gT2nypXk1NjVVWVpLXAXl9kdcXeX2R1xd5fZG3/VjGiBPwsYqKCosJeX2R1xd5fZHXF3l9kddXbHljRMcJG7RXo66uzmKgnKNHjyavE/L6Iq8v8voiry/y+iIv8qFUD22+VE8bEe2FIe/GR15f5PVFXl/k9UVeX+RtP5ZRqgd8LLa9L+T1RV5f5PVFXl/k9UVeX7HljREdJ7T5Ur2xY8eS1wl5fZHXF3l9kdcXeX2RF/lQqodW43ecAAAA0JZEV6p3yy232IABA8IpFIcNG2aTJ09u8f5/+9vfbPfddw/332uvvcLBcNiUMmHUKQbKqRWCvD7I64u8vsjri7y+yOuLvEhlx+mhhx6ySy65xK6++mp76aWXbO+997YRI0bYwoUL895/4sSJ9rWvfc3OOussmzZtmp1wwgnhUlVVtcmzt1exlepNmDCBvE7I64u8vsjri7y+yOuLvEhlqZ5GmIYMGWI333xz+H9DQ4P179/fLrjgArv88svXu/9JJ51kK1eutMcff7zxtgMOOMD22Wcfu+2226Is1Vu5svm/dehgVlnZuvuWl5t16bJh9121Sh2i/PfVyVm6djX74x/NTjuNUj0AAAC0DcX0DUr6S1lr1661qVOn2hVXXNF4W3l5uR155JE2adKkvI/R7RqhyqYRqn/84x95779mzZpwyW4cmTdvXuiA1dfXh/9vttlmobHUU1+yZEm41ukcO3ToEKb79esXsi1YsCD8X9PJ9VZbbRVO/6hfbF6xYkW4XY/T4/X33r17hw7h+++/H24Xva4es/XWW9vmm1dYjx5LrLLy45zrsnazww7b3B5+eGXIpPsffLB+GbrMamsrbPHircP9+vadHzo3++6bsbvuKgvP3atXLxswoNJqa6uta9fVTZ53xYrNbNddu9mECTVWXV0d8p5wQgebP1+PLbeFC/uE+/XuvcA6dGiwHXc0e/DBOrvxxg7WqVNP69Chky1dWm2rV68Or5XMU6dOncLr6n1VOymv5lsXTffp0ydML126NNwnuV3Xavtu3bqFhVbvkW5P3puuXbuGBVr///DDDxtvT96bvn37WseOHcMoZW1tbZP3Rs+b/JJ28t7odl30/iiTrufOndvkvdG03pt17bWicV6TvMrUo0ePcPvixYubzGvnzp1DO2haeZVJ+yeSeerZs6d16dIltIOWQb1W9nKojLpWW2Tn1X223XbbMD1//vzw/NnzqufVc+Uuh5Isa3oOPVb/V6ZknpJlVG2vx2e/N3pf1I56TmVO5kOPT36hXH/PfW80refVcrFo0aKwHibvmfLqPdXjtSzkzqueW++rpt97773G9yuZJ7Wv/q/2W7VqVeN86HFqf7WF5kOvm7scallRXi37em+S2/Ucek/V9vqb5jf3vdnU2wj9TXnVVppXLW/Z740yde/ePfwt2UYk86q21euqTdQOue+N2lD30fum503mI3k/tYyrfZJtRL73RtvRZLudPFYZ1I6a12Qbm9ZtRJJd64Zu1/OmeRuh++gx22yzTWM7pXkbob9p3rbffvvQhvmWwzRtI3Qf3bbTTjuFjFqn0ryNUDYtX8qr+y1fvjzV2wi9j3ovdGiGcmg6zdsIXSvfdtttF9ZHLZ9p3kYol15H86J20P/Tvo3IZDKhHbSMNrccbopthJa11ippx0kLihpFC102/X/GjBl5H/PBBx/kvb9uz+e6666zn/zkJ+vd/sc//jE0dkLHSu24446hQfMdM3X66afbDjvsYH/5y18aO1+JL33pS2Eh0oL+3HPPNfmbFpbvfOc7YT7vueee9Z73W9/6lro+9oUvjLHddnuryd/GjBluZgeGEsZx48aF2zTiI/Pn97Xbb9djzc4++26rqFi3oN1xx7q/f/7znzezg+yww8bbfvtNa/K8EyYcbPPmHWDTp1fZmDFjwm3HHfdxZ+3GG9d1TL/xjT/bFlssD9OK/qUv6eQQ37TevVfZpEnv2pQpU5o8rzbAGinUAvqnP/2pyd+0kJ577rmhpHLmzJnrvV+DBw8O7Th+/Pgwv9l23XXX8P5ohckeaUyceOKJtscee9jDDz8cNjDZhg8fHlY2rVzPPvtsk7/pA0WZtDLne2/OOeecMI9671599dUmf9MH4Nlnn22vvfaajRo1qsnftCE8/vjjwwfG888/Hz4csh1yyCF2xBFHhOVs1qxZTf623377hWut0Do7Tm4bnnHGGeG1H3jggfWe9+STTw7Pp+VQr5tN83/eeeeFDVa+ef3ud78b2l3L8Zw5c5r8bbfddgvPrZ0WKpXNpi9E2pBqw/bEE0+s97zHHntsmKe///3v673nhx9+eOMHwVNPPdXkb/pAUPtqmbr33nvDhjfbmWeeGfJqQ5i7vGh78O1vf9tmz55tjzzyyHrPO2jQoPDeT58+PbRHNo1+H3300fb000+H9zZbqbYR++67b3jv9ZxvvdV0G6E8p556apNtREIfGFrONI9az5MPo4S2EQcddJA99thj9u677zb528EHHxw+zNS+yTYioQ8wlUrrPdd2VB9Y2bQuaj51efHFF1O/jdCXOH0xUt5nnnkm9dsILd9aJvSBn7vOpXEboeVbOzf15TC3HdK4jdBrffOb37T//Oc/9vLLL0exjVA7KGss24hTTjkltF0s2wh9j9Dyq/mKYRuh9UnvrdadGLYRO+ywQ6gcU6eqVNuI3EypLdXTnhO9eQp84IEHNt7+wx/+MGw4clco0YJ7//33h+OcEn/4wx9C5yhfjzHfiJNKAV9//fXQw03DiNOKFRlbuXL9PUVdu25mW22ljB/vKVq5ssEyGe2NqLTu3dftKfrww3U9dJXfde368Z6i+vpKW7r04z1F9fXaw1pm3bptYZtt1tUqKj7eU7R2rXri63rovXuv66HPn79uT1GHDuVWUbFunrQ3Qn9L+56i5L3Re6wNQNr3FCXzpDZIRp3YU7TxR5yURe2g6bTvTc7eRmi+0r43OZknthFsI9hGsI1gG8E2oiGibYSWNZ10rjWleiXtOGlG1WDa46MTPCROO+20sCKOHDlyvcdoD7dK9S666KLG23RiCZXqaQ9yjMc4xUQrgFZwLXBaCNOOvL7I64u8vsjri7y+yOuLvO3HslhOR64e5P77799k6FNvvP6fPQKVTbfnDpVqeK65+2Pj0vujYfLcIc+0Iq8v8voiry/y+iKvL/L6Ii9SeVY9nY5cI0y33367DR061G666aZQY6pjnDT8qdpcDW/rWCVRWd9hhx1mv/zlL+2YY46xBx980H7xi1+EOkUdu1AII04AAAAAohpxSk4vfv3119tVV10VTimuAxyffPLJxhNA6IBE1VImdKCiDqy84447wm8+qcxPZXqt6TThk9OeDNUsx7JHg7y+yOuLvL7I64u8vsjri7xIZcdJzj///HAGDh30pRNC6LedEjoT2n333dfk/l/96lftzTffDPfXsKTOcINNQyukzkQUy4pJXl/k9UVeX+T1RV5f5PVFXqSyVG9To1QPAAAAQHSleoiL9mRodDCWPRrk9UVeX+T1RV5f5PVFXl/kRT50nNCma2jJ64u8vsjri7y+yOuLvL7Ii3wo1QMAAADQLlGqBzf6FedZs2Y1/rJz2pHXF3l9kdcXeX2R1xd5fZEX+dBxQlE0QLl06dJwHQPy+iKvL/L6Iq8v8voiry/yIh9K9QAAAAC0S8so1YMXDQHPmDEjmqFg8voiry/y+iKvL/L6Iq8v8iIfOk4o2urVqy0m5PVFXl/k9UVeX+T1RV5f5EUuSvUAAAAAtEvLKNWDFw0BV1VVRTMUTF5f5PVFXl/k9UVeX+T1RV7kQ8cJAAAAAAqgVA8AAABAu7SsiL5BhbUzST9RjYQNHwoeNGiQdejQwdKOvL7I64u8vsjri7y+yOuLvO3Hsv/rE7RmLKnddZyWL18ervv371/qKAAAAABS0kfQyFNL2l2pXkNDg73//vvWrVs3KysrK3WcKHvl6nTOnTs3ilJH8voiry/y+iKvL/L6Iq8v8rYfmUwmdJq23XZbKy9v+fQP7W7ESQ2y3XbblTpG9LRSxrRiktcXeX2R1xd5fZHXF3l9kbd92LLASFOCs+oBAAAAQAF0nAAAAACgADpOKErnzp3t6quvDtcxIK8v8voiry/y+iKvL/L6Ii/yaXcnhwAAAACAYjHiBAAAAAAF0HECAAAAgALoOAEAAABAAXScAAAAAKAAOk5otVtuucUGDBhglZWVNmzYMJs8ebKl1fjx4+3YY48NvwJdVlZm//jHPyytrrvuOhsyZIh169bNevfubSeccIK9+eablla33nqrDR48uPFH9g488EB74oknLBa//OUvwzJx0UUXWRpdc801IV/2Zffdd7c0mzdvnn3jG9+wnj17WpcuXWyvvfayf//735ZW2o7ltrEu5513nqVNfX29XXnllbbjjjuGth04cKBde+214Zfu02r58uVh/dphhx1C5oMOOsimTJlisXw+qG2vuuoq22abbUL+I4880mbOnJnavH//+99t+PDhYf3T319++WUrpZby1tbW2mWXXRa2EZtttlm4z6mnnmrvv/9+KvMm22Rtg5W3R48eYXl48cUXU5s327e//e1wn5tuummTZmzL6DihVR566CG75JJLwqkuX3rpJdt7771txIgRtnDhQkujlStXhozq7KXdc889F76wvfDCC/bUU0+FDxZ9CGoe0mi77bYLnY+pU6eGL8eHH364HX/88fbaa69Z2unL2+233x46fmn2qU99yubPn994+de//mVptXTpUjv44IOtY8eOoQP9+uuv2w033BC+YKR5OchuX6138tWvftXS5le/+lXYWXHzzTfbG2+8Ef7/61//2n7/+99bWp199tmhTf/0pz/Zq6++GrZn+rKpDnYMnw9q3//5n/+x2267LXxB1hdmfd7V1NRYGvPq74ccckhYNtKgpbyrVq0K3yG0M0DX6vRpR+Fxxx1npVKofXfdddew/mlZ1rZYO160TC9atMjS/P3m0UcfDd8r1MHCRqTTkQOFDB06NHPeeec1/r++vj6z7bbbZq677rpM2mkxf/TRRzOxWLhwYcj83HPPZWLRo0ePzF133ZVJs+XLl2d22WWXzFNPPZU57LDDMhdeeGEmja6++urM3nvvnYnFZZddljnkkEMyMdOyMHDgwExDQ0MmbY455pjMmWee2eS2L3/5y5lTTjklk0arVq3KdOjQIfP44483uX2//fbL/PjHP86k/fNBy0Dfvn0zv/nNbxpvq66uznTu3Dnz17/+NZPmz7N33nkn/H3atGmZmD5/J0+eHO43Z86cTAx5P/roo3C/p59+OpPWvO+9916mX79+maqqqswOO+yQ+e1vf1uSfG0RI04oaO3atWF0QXsME+Xl5eH/kyZNKmm2tuijjz4K11tttZWlncqIHnzwwbAHTCV7aaZRvWOOOabJcpxWKgvSXsKddtrJTjnlFHv33XctrR577DH79Kc/HUZrVGq677772p133mkxbd8eeOABO/PMM0NJS9qozO2ZZ56xt956K/x/+vTpYa/3UUcdZWlUV1cXtgsq6c6mkrc0j5wm3nnnHfvggw+abCe23HLLUJ7O553fZ57Wve7du1sM24s77rgjLBMa9UmjhoYG++Y3v2mXXnppqF7AxlWxkZ8PbdDixYvDB2GfPn2a3K7/z5gxo2S52iJt8HRsgEqfBg0aZGmlkgV1lFS6svnmm4eSgD333NPSSp07lYWk6TiL5ugL2n333We77bZbKCP7yU9+YoceeqhVVVWF4+DS5u233w6lZCrl/dGPfhTa+Hvf+5516tTJTjvtNEs7HR9QXV1tp59+uqXR5ZdfbsuWLQvHWHTo0CFsi3/+85+HDnUaaRnVtkHHYe2xxx7hc+Kvf/1r6HTsvPPOlnbqNEm+z7vkb9h49BmiY56+9rWvhWNm0+rxxx+3k08+OZQa6tg3laL26tXL0kglmxUVFWE7jI2PjhOQslERfUFO+55ZfanXAcjaU/jII4+EL8g6ViuNnae5c+fahRdeGD7ocveCp1H2SIKOxVJHSgfZP/zww3bWWWdZGjv7GnH6xS9+Ef6vESctwzo+JIaO09133x3aPK3HAeh9//Of/2x/+ctfwt5jrXfauaK8aW1fHdukEbx+/fqFzt5+++0XvhircgFI6HjeE088MZyMQztf0uxzn/tcWPe0I1kj6sqt4980yp4mWsd+97vfhR2FaRxBbwso1UNB2quiD78FCxY0uV3/79u3b8lytTXnn39+2Kv1z3/+M5yAIc00mqC9x/vvv384K6BKFrSxTiN9kOgkJvrypr1wuqiTp4O/Na09+Gmm8hUdnDxr1ixLI+19ze0wa6QhzeWFiTlz5tjTTz8dTmaQViq30aiT9nbrTGQqwbn44ovDepdWOvOf1rEVK1aEHRc6A6u+JKv0NO2SzzQ+7zZNp0nroHZqpXm0SXSCEH3mHXDAAWFniz47dJ02EyZMCJ9322+/fePnndr4+9//fjipBT45Ok5o1ZdkfUFWnX32Xmb9P+3HtcRAe9vUaVK527hx48Jph2Oj5WHNmjWWRkcccUQoLdTewuSiERKVOmlaOwXSTF8+Z8+eHTooaaSy0tzT5+t4HI2Spd29994b9hjr2Le0UmmQjinNpmVW61za6cumlludeXHMmDHh7Jtpp+2vOkjZn3cqldToAp93G7fTpGM5teNCp1GPTVo/87Rj5ZVXXmnyeafRae2A0TqIT45SPbSKjl9QWYi+cA4dOjT8JoBOCHDGGWdYWr9sZu+h1wG/2oDohAvaE5O28jyV4YwcOTIcH5DU0evgUx1QnTZXXHFFKG1SO+r3WpT92WefTe1GWW2ae7yYvtDpwzqNx5H94Ac/CL/RoY6HfttEPwGgL8oqdUojjX7oBAYq1dOXIY0u6OBpXdL+xUcdJ23XtFc2rbQs6JgmrW8q1Zs2bZrdeOONoRQurbQt0A4hlfRqO6wvbTpGKy2fF4U+H1QK+bOf/cx22WWX0JHSqbP15VO/sZfGvEuWLAkjvMlvISU7MtQBLMUoWUt51ZH+yle+EkrJVGGhEf/kM09/147aNOXV54TWP50uXdlVqqfTgOvU+qX6+YJCy0NuR1Q/FaHlQOsjNoJSn9YP8fj973+f2X777TOdOnUKpyd/4YUXMmn1z3/+M5ymM/dy2mmnZdImX05d7r333kwa6dTIOr2ploOtt946c8QRR2TGjh2biUmaT0d+0kknZbbZZpvQvjqdrP4/a9asTJr97//+b2bQoEHhlM2777575o477sik3ZgxY8J69uabb2bSbNmyZWFZ1ba3srIys9NOO4XTeq9ZsyaTVg899FDIqWVYp/bWT1nolN6xfD7olORXXnllpk+fPmGZ1jaulMtJobz6rMj3d/20QdryJqdMz3fR49KWd/Xq1ZkvfelL4edXtDxr23zccceFU6jH8v2G05FvXGX6Z2N0wAAAAACgreIYJwAAAAAogI4TAAAAABRAxwkAAAAACqDjBAAAAAAF0HECAAAAgALoOAEAAABAAXScAAAAAKAAOk4AAAAAUAAdJwAAAAAogI4TAAAAABRAxwkAAAAACqDjBABoNxYtWmR9+/a1X/ziF423TZw40Tp16mTPPPNMSbMBANKtLJPJZEodAgCATWX06NF2wgknhA7TbrvtZvvss48df/zxduONN5Y6GgAgxeg4AQDanfPOO8+efvpp+/SnP22vvvqqTZkyxTp37lzqWACAFKPjBABod1avXm2DBg2yuXPn2tSpU22vvfYqdSQAQMpxjBMAoN2ZPXu2vf/++9bQ0GD/+c9/Sh0HABABRpwAAO3K2rVrbejQoeHYJh3jdNNNN4Vyvd69e5c6GgAgxeg4AQDalUsvvdQeeeQRmz59um2++eZ22GGH2ZZbbmmPP/54qaMBAFKMUj0AQLvx7LPPhhGmP/3pT7bFFltYeXl5mJ4wYYLdeuutpY4HAEgxRpwAAAAAoABGnAAAAACgADpOAAAAAFAAHScAAAAAKICOEwAAAAAUQMcJAAAAAAqg4wQAAAAABdBxAgAAAIAC6DgBAAAAQAF0nAAAAACgADpOAAAAAFAAHScAAAAAsJb9f5883IOM+pnEAAAAAElFTkSuQmCC",
|
||
"text/plain": [
|
||
"<Figure size 1000x600 with 1 Axes>"
|
||
]
|
||
},
|
||
"metadata": {},
|
||
"output_type": "display_data"
|
||
}
|
||
],
|
||
"source": [
|
||
"import matplotlib.pyplot as plt\n",
|
||
"from statsmodels.distributions.empirical_distribution import ECDF\n",
|
||
"\n",
|
||
"ecdf = ECDF(data)\n",
|
||
"x = np.linspace(min(data) - 1, max(data) + 1, 1000)\n",
|
||
"y = ecdf(x)\n",
|
||
"\n",
|
||
"# Находим точки, где F(x) переходит от 0 к основному росту и от роста к 1\n",
|
||
"x_left = x[y == 0][-1] # Последняя точка, где F(x)=0\n",
|
||
"x_right = x[y == 1][0] # Первая точка, где F(x)=1\n",
|
||
"\n",
|
||
"# Разделяем данные на 3 части\n",
|
||
"mask_left = (x < x_left) # F(x) = 0\n",
|
||
"mask_mid = (x >= x_left) & (x <= x_right) # Основной рост\n",
|
||
"mask_right = (x > x_right) # F(x) = 1\n",
|
||
"\n",
|
||
"# Рисуем каждую часть своим стилем\n",
|
||
"plt.figure(figsize=(10, 6))\n",
|
||
"plt.step(x[mask_left], y[mask_left], '--', color='blue', where='post', label='F(x)=0') # Пунктир слева\n",
|
||
"plt.step(x[mask_mid], y[mask_mid], '-', color='blue', where='post', label='ЭФР') # Сплошная основная часть\n",
|
||
"plt.step(x[mask_right], y[mask_right], '--', color='blue', where='post', label='F(x)=1') # Пунктир справа\n",
|
||
"\n",
|
||
"# Настройки графика\n",
|
||
"plt.title(\"Эмпирическая функция распределения\")\n",
|
||
"plt.xlabel(\"x\")\n",
|
||
"plt.ylabel(\"F(x)\")\n",
|
||
"# Добавление пунктирных линий для F(x) = 0 и F(x) = 1\n",
|
||
"plt.axhline(y=0, color='gray', linestyle='--', linewidth=1, label='F(x) = 0')\n",
|
||
"plt.axhline(y=1, color='gray', linestyle='--', linewidth=1, label='F(x) = 1')\n",
|
||
"\n",
|
||
"plt.grid(True, linestyle=':')\n",
|
||
"plt.xticks(np.arange(np.floor(min(data)), np.ceil(max(data)) + 1))\n",
|
||
"plt.show()"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "639c228f",
|
||
"metadata": {},
|
||
"source": [
|
||
"### 3. Гистограмма частот"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 5,
|
||
"id": "09541433",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAHHCAYAAABXx+fLAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvc2/+5QAAAAlwSFlzAAAPYQAAD2EBqD+naQAAO2FJREFUeJzt3QlcVPXex/EfO2LgBook4pL7rpVplltp5rVcbotWrll59VVmeU3T1KxstSx9tDS1W7lkLu2a5ZZXve6VlaaGuCtogoICwTyv3/95hgs4IOjMYeB83q/XEebMmfmfMzPO+fL//845Pg6HwyEAAAA24lvUKwAAAGA1AhAAALAdAhAAALAdAhAAALAdAhAAALAdAhAAALAdAhAAALAdAhAAALAdAhAAALAdAhAAALAdAhDghebNmyc+Pj55TkeOHCnqVQSAYs2/qFcAQN6ef/55qV69+iXzy5cvXyTrAwAlBQEI8GJdunSR66+/vqhXAwBKHIbAgBIwVHbw4MGseZmZmdK4cWMzX+/Pbs+ePXLvvfdKRESElCpVSurUqSPPPvusuW/ChAn5DrvptHbt2qznWrx4sbRo0cI8T3h4uDz44INy9OjRHO3179/f5fNcd911WctUq1ZN/va3v8m3334rTZs2leDgYKlfv74sXbo0x3OdOXNGnn76aWnUqJFcc801EhYWZgLijz/+mGM5XUdnO7t27cpxn66fn5+fue/TTz+9ZD21/dwmT55s7tM2s5s7d6506NBBKlasKEFBQWadZ8yYcZl37L/brG26Wu/sr/EPP/wg99xzj1StWtW0ER0dLU8++aRcuHDhkue08r3N633NPmX/TALeiB4goIT58MMP5eeff75k/k8//SS33HKLBAQEyCOPPGJ2wgcOHJAvvvhCXnzxRenZs2eOYKI72nr16pllnfS20mA1YMAAueGGG0xAOHnypEydOlX+/e9/y86dO6Vs2bJZj9Ed9+zZs3OsS2hoaI7b+/btk/vuu08ee+wx6devnwkXuuNfsWKF3H777WaZP/74Q5YvX27m67Cgtvnuu+9K27Zt5ddff5WoqKgcz6lBSp9H18vpgw8+kMDAQLl48eIlr4+/v7/88ssvZv2bNWuWNV+3VZ8rNw07DRo0kLvuuss8Vl/Hf/zjHyaADh06VNxBg0hKSooMGTJEKlSoIFu2bJF33nnH1IDpfUX13j766KNy2223ZT32oYcekh49eph2nDSIAV7NAcDrzJ0716H/Pbdu3Vqg5WJjY83tixcvOqpWrero0qWLma/3O916662O0NBQR1xcXI7nyMzMdPncMTExjn79+l0yPy0tzVGxYkVHw4YNHRcuXMia/+WXX5o2n3vuuax5+vjSpUvnuw3ajj5uyZIlWfMSExMdlStXdjRr1ixrnm5bRkZGjsfqdgcFBTmef/75rHlr1qwxz9e7d29HhQoVHKmpqVn31apVy9GnTx9z/+LFiy9Zz27dujmGDRuWNf+HH35wlCpVytG9e/dLtiMlJeWSbencubOjRo0ajsupXr26o2/fvjnmOddbf+bXxuTJkx0+Pj453seieG+z0/vGjx+f7zYD3oYhMKAEmT59upw+fVrGjx+fY358fLysX79eBg4caIZTstPhisLYtm2bnDp1yvR2ZO8Z6dq1q9StW1e++uqrQq+39t5oD4KTDm/17dvX9DicOHEiqyfJ1/f/vrIyMjLMduqwlA717Nix45Ln7Natm9m2zz//PGs4SXtOtKcpL/r6zJ8/X1JTU81t7UHSXo0yZcpcsqwODzklJiZKQkKC6Y3Sniq9nR8dNivIkXzZ20hOTjZttG7dWv9wNa9NcXhvAW9FAAJKCN3pvvTSSzJixAipVKlSjvt0p6waNmx41e3ExcWZnxo8ctOdpPP+wtDhmdw769q1a5ufzloSHVp68803pVatWiYMaW2KDrPo8I+rwKHDQVq7MmfOHHNbf/bq1cuEq7zojl6Hsz777DMTOD755BMzHOSKDgnpMFDp0qXNsJCuy5gxY8x9lwtAGmLWrVsnCxcuNIFDg42rxxw6dMjU2+hRfxr2tA0NWdnb8Pb3FvBW1AABJcQrr7xiekhGjhxpekdKGg1348aNMz0dkyZNMqFAt3f48OEmHLmiy2o9z969e03NjLM3KC/O0KQ9P1p7o3U3WuisdVXZaX1Nx44dTSiYMmWKKU7W2qKvv/7ahLS81sdJg5IGqN69e+e5jPZyaf2TFn+PGjXKtKVhS4uRNRRdrg0A+SMAASXAsWPHTKGqFq1qgXHuAFSjRg3zc/fu3VfdVkxMjPmpoULDQXY6z3l/Yezfv98M62TvBfr999/NTy3oVXrUVvv27eX999/P8dizZ8+a3iBX9IgxDUDOo6P08drzkh8NTU2aNJHDhw+bgmxXw0haXKzDZBqosg87rVmzpkDbq+u7adMmU7ztHOLTo9n0KDcnLWTX10ALt3U40GnVqlXF6r0FvBVDYEAJMHHiRDPspUdRuaI7/1tvvdUMA+mwSnb/V8NacHpeIq1hmTlzZlatjPrmm2/kt99+M8NIVxLgli1blnU7KSlJ/vWvf5nD0iMjI808PXw997pqr07uQ+9dBRodJnMeun05emSXHgKu4ST3oepOui4q+/rokJT2HBWU9l7psJUOo+mkbV6uDf09+1FtxeG9BbwVPUBACaDn0Pn444/NMExe3n77bWnTpo00b97cHP6sh5JrfY0WtuY+X87lhol0uE1rY7QeRYdxnIdKa2+NHmJdWFrvM2jQINm6dasJcroz1+fMHij0XEF6ZmxtV2totIdEt9nZA5KXwYMHm0PnXRUy52X16tUmAOR1xu1OnTqZ11oLrfWQ8PPnz8usWbNMeDh+/Li4gw551axZ0/QKacjT2qUlS5bIn3/+WazeW8BbEYCAEkB7SvKrJ1E6rLN582ZTR6PnsNFz4eiQhg4PFZb2jISEhMjLL79s6lO0NkWP4tKdZ/ZzABWUFjbr+W20fkmHWnQHvmjRIuncuXOOuhktTNajtPQ+3dnrDv6ZZ57J97m1qDmvIbK86PbolBctEtYhubFjx5qAor1Ueq4e7Y3RHid30DCiQ22PP/64GdrUo7L0NR42bJh5L4vLewt4Kx89Fr6oVwKAfWnPgg4Fffnll0W9KgBshBogAABgOwQgAABgOwQgAABgO9QAAQAA26EHCAAA2A4BCAAA2A7nAXJBr7GjZ6bVSwoU9mrKAACgaGhVz7lz5yQqKsqcbT0/BCAXNPzoxQ0BAEDxo9fyq1KlSr7LEIBc0J4f5wuop58HAADeT68jqB0Yzv14fghALjiHvTT8EIAAACheClK+QhE0AACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHQIQAACwHf+iXgE7io+Pl6SkJMvbDQsLk4iICMvbBQDA2xCAiiD8PDjgYTlzLsXytsuHhshHc2cTggAAtkcAspj2/Gj4iWjVS0qXr2RZu8lnTkr8piWmfQIQAMDuCEBFRMNPWMUqlrYZb2lrAAB4L4qgAQCA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RCAAACA7RRpAFq/fr1069ZNoqKixMfHR5YvX57jfp3nanrttdfyfM4JEyZcsnzdunUt2BoAAFBcFGkASk5OliZNmsj06dNd3n/8+PEc05w5c0yg6dWrV77P26BBgxyP27Bhg4e2AAAAFEf+Rdl4ly5dzJSXyMjIHLc/++wzad++vdSoUSPf5/X397/ksQAAAMWuBujkyZPy1VdfyaBBgy677L59+8ywmgalBx54QA4dOmTJOgIAgOKhSHuACuODDz6Q0NBQ6dmzZ77LtWzZUubNmyd16tQxw18TJ06UW265RXbv3m0e70pqaqqZnJKSkszP9PR0M7lTRkaGBAYGSICvvviZYhVtT9vV9t29TQAAeIPC7N98HA6HQ7yA1vYsW7ZMunfv7vJ+LWS+/fbb5Z133inU8549e1ZiYmJkypQpefYeaeG0BqXc5s+fLyEhIYVqDwAAFI2UlBTp06ePJCYmSlhYWPHvAfrhhx9k7969smjRokI/tmzZslK7dm3Zv39/nsuMHj1aRowYkaMHKDo6Wjp16nTZF7CwYmNjZcDQ4RLTebCEhkeJVc4lHJO4lbNk7vS3pHr16pa1CwCAVZwjOAVRLALQ+++/Ly1atDBHjBXW+fPn5cCBA/LQQw/luUxQUJCZcgsICDCTO/n5+UlaWrqkZ4r8ZWEJlran7Wr77t4mAAC8QWH2b0VaBK3hZNeuXWZy9o7o79mLljXNLV68WB5++GGXz9GxY0eZNm1a1u2nn35a1q1bJwcPHpSNGzdKjx49zE6/d+/eFmwRAAAoDoq0B2jbtm3msHYn5zBUv379TCGzWrhwoWiZUl4BRnt3EhISsm4fOXLELHv69GmJiIiQNm3ayObNm83vAAAARR6A2rVrZ8JNfh555BEz5UV7erLTwAQAAFAizgMEAADgLgQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgO0UagNavXy/dunWTqKgo8fHxkeXLl+e4v3///mZ+9umOO+647PNOnz5dqlWrJsHBwdKyZUvZsmWLB7cCAAAUN0UagJKTk6VJkyYmsORFA8/x48ezpgULFuT7nIsWLZIRI0bI+PHjZceOHeb5O3fuLKdOnfLAFgAAgOLIvygb79Kli5nyExQUJJGRkQV+zilTpsjgwYNlwIAB5vbMmTPlq6++kjlz5sgzzzxz1esMAACKvyINQAWxdu1aqVixopQrV046dOggL7zwglSoUMHlsmlpabJ9+3YZPXp01jxfX1+57bbbZNOmTXm2kZqaaianpKQk8zM9Pd1M7pSRkSGBgQES4KsvfqZYRdvTdrV9d28TAADeoDD7N68OQDr81bNnT6levbocOHBAxowZY3qMNMz4+fldsnxCQoLZwVeqVCnHfL29Z8+ePNuZPHmyTJw48ZL53377rYSEhIi7/XPYo///2wmxTDlfkVqPym+//WYmAABKmpSUlJIRgO6///6s3xs1aiSNGzeWmjVrml6hjh07uq0d7THSuqHsPUDR0dHSqVMnCQsLE3eKjY2VAUOHS0znwRIaHiVWOZdwTOJWzpK5098ygRIAgJLGOYJT7ANQbjVq1JDw8HDZv3+/ywCk92nP0MmTJ3PM19v51RFpnZFOuQUEBJjJnXT90tLSJT1T5C8La9C1PW1X23f3NgEA4A0Ks38rVucBOnLkiJw+fVoqV67s8v7AwEBp0aKFfP/991nzMjMzze1WrVpZuKYAAMCbFWkAOn/+vOzatctMzuEh/f3QoUPmvpEjR8rmzZvl4MGDJsTcfffdct1115nD2p20J2jatGlZt3Uoa9asWfLBBx+YWpchQ4aYw+2dR4UBAAAU6RDYtm3bpH379lm3nXU4/fr1kxkzZshPP/1kgszZs2fNyRK1JmfSpEk5hqu0OFqLn53uu+8+iY+Pl+eee05OnDghTZs2lRUrVlxSGA0AAOyrSANQu3btxOFw5Hn/ypUrL/sc2juU27Bhw8wEAABQ7GuAAAAA3IEABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbIcABAAAbKdIA9D69eulW7duEhUVJT4+PrJ8+fKs+9LT02XUqFHSqFEjKV26tFmmb9++cuzYsXyfc8KECea5sk9169a1YGsAAEBxUaQBKDk5WZo0aSLTp0+/5L6UlBTZsWOHjBs3zvxcunSp7N27V+66667LPm+DBg3k+PHjWdOGDRs8tAUAAKA48i/Kxrt06WImV8qUKSOrVq3KMW/atGly4403yqFDh6Rq1ap5Pq+/v79ERka6fX0BAEDJUKQBqLASExPNkFbZsmXzXW7fvn1myCw4OFhatWolkydPzjcwpaammskpKSkpaxhOJ3fKyMiQwMAACfDVFz9TrKLtabvavru3CQAAb1CY/ZuPw+FwiBfQYLNs2TLp3r27y/svXrwoN998s6nn+fjjj/N8nm+++UbOnz8vderUMcNfEydOlKNHj8ru3bslNDQ0z7ohXS63+fPnS0hIyFVsFQAAsIqWz/Tp08d0mISFhRX/AKSJrlevXnLkyBFZu3btZTcqu7Nnz0pMTIxMmTJFBg0aVOAeoOjoaElISChUWwURGxsrA4YOl5jOgyU0PEqsci7hmMStnCVzp78l1atXt6xdAACsovvv8PDwAgUgrx8C0/Bz7733SlxcnKxevbrQgUSHy2rXri379+/Pc5mgoCAz5RYQEGAmd/Lz85O0tHRJzxT5y8IadG1P29X23b1NAAB4g8Ls33yLQ/jRmp7vvvtOKlSoUOjn0OGwAwcOSOXKlT2yjgAAoPgp0gCk4WTXrl1mcg4P6e96lJeGn7///e+ybds2U/OjxbsnTpwwU1paWtZzdOzY0Rwd5vT000/LunXr5ODBg7Jx40bp0aOH6fXo3bt3kWwjAADwPkU6BKbhpn379lm3R4wYYX7269fPFCZ//vnn5nbTpk1zPG7NmjXSrl0787v27mitjpPWCWnYOX36tEREREibNm1k8+bN5ncAAIAiD0AaYvKrwS5Ifbb29GS3cOFCt6wbAAAouby6BggAAMATCEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2/K/0gcnJybJu3To5dOiQpKWl5bjv8ccfd8e6AQAAeE8A2rlzp9x5552SkpJiglD58uUlISFBQkJCpGLFigQgAABQ8obAnnzySenWrZv8+eefUqpUKdm8ebPExcVJixYt5PXXX3f/WgIAABR1ANq1a5c89dRT4uvrK35+fpKamirR0dHy6quvypgxY9y5fgAAAN4RgAICAkz4UTrkpXVAqkyZMnL48GH3riEAAIA31AA1a9ZMtm7dKrVq1ZK2bdvKc889Z2qAPvzwQ2nYsKG71xFukp6WZoYqrRYWFiYRERGWtwsAgFsD0EsvvSTnzp0zv7/44ovSt29fGTJkiAlEc+bMuZKnhIelnk+Ug7F/yPAxEyQoKMjStsuHhshHc2cTggAAxTsAXX/99Vm/6xDYihUr3LlO8ID01AuS6eMv4Tf1lApRMZa1m3zmpMRvWiJJSUkEIABA8Q5AHTp0kKVLl0rZsmXdv0bwqJByERJWsYqlbcZb2hoAAB4qgl67du0lJz8EAAAo8ZfC8PHxce+aAAAAePulMHr06CGBgYEu71u9evXVrBMAAIB3BqBWrVrJNddc4961AQAA8NYApMNfI0eONEeAAQAA2KIGyOFwuH9NAAAAvDkAjR8/nuEvAABgryEwDUAqPj5e9u7da36vU6cOJ7oDAAAltwcoJSVFBg4cKFFRUXLrrbeaSX8fNGiQuQ8AAKDEBaAnn3xS1q1bJ59//rmcPXvWTJ999pmZ99RTT7l/LQEAAIp6CGzJkiXy6aefSrt27bLm3XnnnVKqVCm59957ZcaMGe5cRwAAAO8YAqtUqdIl8/Ww+MIMga1fv166detmhs/00Prly5dfcrTZc889J5UrVzbh6rbbbpN9+/Zd9nmnT58u1apVk+DgYGnZsqVs2bKlwOsEAABKPt8rPQmiFkJfvHgxa96FCxdk4sSJ5r6CSk5OliZNmpjA4sqrr74qb7/9tsycOVP+85//SOnSpaVz58452s1t0aJFMmLECLN+O3bsMM+vjzl16lQhtxIAAJRUVzQE9tZbb8kdd9whVapUMQFD/fjjj6bHZeXKlQV+ni5dupjJFe390XbGjh0rd999t5n3r3/9y/Q8aU/R/fff7/JxU6ZMkcGDB8uAAQPMbQ1PX331lcyZM0eeeeaZK9haAABQ0lxRAGrUqJEZivr4449lz549Zl7v3r3lgQceMENV7hAbGysnTpwww15OZcqUMUNamzZtchmA9Ar127dvl9GjR2fN8/X1Nc+hj8lLamqqmZySkpLMz/T0dDO5U0ZGhgQGBkiAr774mWKVAD8fCQ4Osr5dXzHbq9vt7tcSAIDsCrOfuaIApLU7rVu3Nj0tnqLhR+WuNdLbzvtyS0hIMDtaV49xBjVXJk+ebIbvcvv2228lJCRE3O2fwx79/99cb4dHtKwsg1q+Yn275XxFaj0qv/32m5kAAPCUwtQhX1EAat++vRw/frzEXAtMe4y0bih7D1B0dLR06tRJwsLC3NqW9mwNGDpcYjoPltDwKLHK8d93yuYFb8nNA8dJxejrLGv3XMIxiVs5S+ZOf0uqV69uWbsAAPtJ+v8RHI8FICuuBRYZGWl+njx50hwF5qS3mzZt6vIx4eHh4ufnZ5bJTm87n8+VoKAgM+UWEBBgJnfS9UtLS5f0TJG/rqwG/YqkZzjk4sVU69vN1KHJdLPd7n4tAQDIrjD7mSsKQEprasqVK+fyPj0z9NXS3gINLd9//31W4NFkp0eDDRkyxOVjAgMDpUWLFuYx3bt3N/MyMzPN7WHDhl31OgEAgJLhigNQjx49XM7X8/loHU5BnD9/Xvbv359jeGjXrl1Svnx5qVq1qgwfPlxeeOEFqVWrlglE48aNM+cMcoYb1bFjR7MuzoCjQ1n9+vWT66+/Xm688UZzJJkebu88KgwAAOCKA5AWIl9tDdC2bdtMPZGTsw5HA8y8efPkn//8pwkvjzzyiLncRps2bWTFihXmcHunAwcOmOJnp/vuu89cpFVPoKjrqL1H+hhXJ24EAAD2dEUBSHt53EEvpZFfPZG28/zzz5spLwcPHrxknvYGMeQFAADy4uutRdAAAABe1QOkhcUAAAC26gHSEwfqpSVy03mvvOI82R4AAEAJCkDvvvuu1K1b95L5DRo0MNfeAgAAKHEBSI+uyn5yQqeIiAhzhmgAAIASF4D0MhH//ve/L5mv8/Q8PQAAACWuCFovgqonKdSrrnbo0MHM07Mt63l7nnrqKXevIwAAQNEHoJEjR8rp06flH//4h6SlpZl5enLCUaNGmQuLAgAAlMgTIerRXnppit9++01KlSplLlfh6oKiAAAAJeZSGOqaa66RG264wX1rAwAA4M0BSK/j9cknn8ihQ4eyhsGcli5d6o51AwAA8J6jwBYuXCitW7c2w1/Lli0zxdC//PKLrF69WsqUKeP+tQQAACjqAPTSSy/Jm2++KV988YUEBgbK1KlTZc+ePXLvvfdK1apV3bl+AAAA3hGADhw4IF27djW/awBKTk42hdFPPvmkvPfee+5eRwAAgKIPQOXKlZNz586Z36+99lrZvXu3+f3s2bOSkpLi3jUEAADwhiLoW2+9VVatWiWNGjWSe+65R5544glT/6PzOnbs6O51BAAAKPoANG3aNLl48aL5/dlnn5WAgADZuHGj9OrVS8aOHeveNQQAACjKAJSUlPR/D/L3N+cAct7WM0LrBAAAUOICUNmyZU2x8+VkZGRczToBAAB4TwBas2ZNjtsOh0PuvPNOmT17timGBgAAKHEBqG3btpfM8/Pzk5tuuklq1KjhzvUCAADwrsPgAQAAbBuADh8+bM77U6FCBfetEQAAgDcNgb399ttZvyckJMiCBQukQ4cOXP8LAACU3ACk1/9SeiRYeHi4dOvWjfP+AACAkh2AYmNjPbcmAAAAFqEIGgAA2A4BCAAA2A4BCAAA2A4BCAAA2A4BCAAA2A4BCAAA2A4BCAAA2A4BCAAA2A4BCAAA2I7XB6Bq1aqZS2/knoYOHepy+Xnz5l2ybHBwsOXrDQAASsilMIrC1q1bJSMjI+v27t275fbbb5d77rknz8eEhYXJ3r17s25rCAIAACg2ASgiIiLH7Zdffllq1qwpbdu2zfMxGngiIyMtWDsAAFAcef0QWHZpaWny0UcfycCBA/Pt1Tl//rzExMRIdHS03H333fLLL79Yup4AAMC7eX0PUHbLly+Xs2fPSv/+/fNcpk6dOjJnzhxp3LixJCYmyuuvvy6tW7c2IahKlSouH5Oammomp6SkJPMzPT3dTO6kw3mBgQES4KsvfqZYJcBPa6GCrG/XV8z26na7+7UEACC7wuxnfBwOh0OKic6dO0tgYKB88cUXhXox6tWrJ71795ZJkya5XGbChAkyceLES+bPnz9fQkJCrmqdAQCANVJSUqRPnz6mA0TrgUtEAIqLi5MaNWrI0qVLzbBWYWjBtL+/vyxYsKDAPUA6fJaQkHDZF7CwYmNjZcDQ4RLTebCEhkeJVY7/vlM2L3hLbh44TipGX2dZu+cSjkncylkyd/pbUr16dcvaBQDYT1JSkoSHhxcoABWbIbC5c+dKxYoVpWvXroV6nA69/Pzzz3LnnXfmuUxQUJCZcgsICDCTO/n5+UlaWrqkZ4r8ZWEJVnqGQy5eTLW+3Uyt3Uo32+3u1xIAgOwKs58pFkXQmZmZJgD169fP9ORk17dvXxk9enTW7eeff16+/fZb+eOPP2THjh3y4IMPmt6jhx9+uAjWHAAAeKNi0QP03XffyaFDh8zRX7npfF/f/+a4P//8UwYPHiwnTpyQcuXKSYsWLWTjxo1Sv359i9caAAB4q2IRgDp16iR5lSqtXbs2x+0333zTTAAAAMV6CAwAAMCdCEAAAMB2isUQGIq39LQ0U4huNT0EMvelVAAAUAQgeFTq+UQ5GPuHDB8zweWpBjypfGiIfDR3NiEIAHAJAhA8Kj31gmT6+Ev4TT2lQlSMZe0mnzkp8ZuWmJNiEYAAALkRgGCJkHIRElbR9bXYPCXe0tYAAMUJRdAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2vDoATZgwQXx8fHJMdevWzfcxixcvNssEBwdLo0aN5Ouvv7ZsfQEAQPHg1QFINWjQQI4fP541bdiwIc9lN27cKL1795ZBgwbJzp07pXv37mbavXu3pesMAAC8m9cHIH9/f4mMjMyawsPD81x26tSpcscdd8jIkSOlXr16MmnSJGnevLlMmzbN0nUGAADezV+83L59+yQqKsoMabVq1UomT54sVatWdbnspk2bZMSIETnmde7cWZYvX55vG6mpqWZySkpKMj/T09PN5E4ZGRkSGBggAb764meKVQL8fCQ4OMg+7fqKeZ319Xb3ewgA8E6F+b73cTgcDvFS33zzjZw/f17q1Kljhr8mTpwoR48eNUNaoaGhlywfGBgoH3zwgRkGc/qf//kf87iTJ0/mW2uky+Q2f/58CQkJceMWAQAAT0lJSZE+ffpIYmKihIWFFd8eoC5dumT93rhxY2nZsqXExMTIJ598Yup83GX06NE5eo60Byg6Olo6dep02RewsGJjY2XA0OES03mwhIZHiVWO/75TNi94S24eOE4qRl9X4ts9l3BM4lbOkrnT35Lq1atb1i4AoOg4R3AKwqsDUG5ly5aV2rVry/79+13erzVCuXt69LbOz09QUJCZcgsICDCTO/n5+UlaWrqkZ4r8ZWEJVnqGQy5eTLVPu5liXmd9vd39HgIAvFNhvu+9vgg6Ox0OO3DggFSuXNnl/Voj9P333+eYt2rVKjMfAACgWASgp59+WtatWycHDx40h7j36NHD/EXvrPHp27evGb5yeuKJJ2TFihXyxhtvyJ49e0xtz7Zt22TYsGFFuBUAAMDbePUQ2JEjR0zYOX36tEREREibNm1k8+bN5nd16NAh8fX9b4Zr3bq1KVweO3asjBkzRmrVqmWOAGvYsGERbgUAAPA2Xh2AFi5cmO/9a9euvWTePffcYyYAAIBiOQQGAADgCQQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgOwQgAABgO159LTDgaqSnpUlcXJzl7YaFhWVdsBcA4J0IQCiRUs8nysHYP2T4mAkSFBRkadvlQ0Pko7mzCUEA4MUIQCiR0lMvSKaPv4Tf1FMqRMVY1m7ymZMSv2mJJCUlEYAAwIsRgFCihZSLkLCKVSxtM97S1gAAV4IiaAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDsEIAAAYDv+Rb0CQEmTnpYmcXFxlrcbFhYmERERlrcLAMURAQhwo9TziXIw9g8ZPmaCBAUFWdp2+dAQ+WjubEIQABQAAQhwo/TUC5Lp4y/hN/WUClExlrWbfOakxG9aIklJSQQgACgAAhDgASHlIiSsYhVL24y3tDUAKN4oggYAALZDAAIAALbj1QFo8uTJcsMNN0hoaKhUrFhRunfvLnv37s33MfPmzRMfH58cU3BwsGXrDAAAvJ9XB6B169bJ0KFDZfPmzbJq1SpJT0+XTp06SXJy8mUPBz5+/HjWVBSHJAMAAO/l1UXQK1asuKR3R3uCtm/fLrfeemuej9Nen8jISAvWEAAAFEdeHYByS0xMND/Lly+f73Lnz5+XmJgYyczMlObNm8tLL70kDRo0yHP51NRUMznpocRKe5x0cqeMjAwJDAyQAF998TPFKgF+OhQYRLsltV1fMZ8r/Xy5+zMLAMVFYb7/fBwOh0OKAQ0zd911l5w9e1Y2bNiQ53KbNm2Sffv2SePGjU1gev3112X9+vXyyy+/SJUqrg9LnjBhgkycOPGS+fPnz5eQkBC3bgcAAPCMlJQU6dOnj9n/azlMiQhAQ4YMkW+++caEn7yCTF5psF69etK7d2+ZNGlSgXuAoqOjJSEh4bIvYGHFxsbKgKHDJabzYAkNjxKrHP99p2xe8JbcPHCcVIy+jnZLWLvnEo5J3MpZMnf6W1K9enXL2gUAb6L77/Dw8AIFoGIxBDZs2DD58ssvTU9OYcKPCggIkGbNmsn+/fvzXEYvWeDqsgX6WJ3cyc/PT9LS0iU9U+QvC2vQ0zMccvFiKu2W1HYzxXyu9PPl7s8sABQXhfn+8+qjwLRzSsPPsmXLZPXq1Vf0l63WRPz8889SuXJlj6wjAAAofry6B0gPgdc6nM8++8ycC+jEiRNmfpkyZaRUqVLm9759+8q1115rzhmknn/+ebnpppvkuuuuM/VCr732mjkM/uGHHy7SbQEAAN7DqwPQjBkzzM927drlmD937lzp37+/+f3QoUPi6/vfjqw///xTBg8ebMJSuXLlpEWLFrJx40apX7++xWsPAAC8lVcHoILUZ69duzbH7TfffNNMAAAAxbIGCAAAwBMIQAAAwHa8eggMgPeLj4/POnu6lfQcHxEREZa3a7ftBUoqAhCAqwoDDw54WM6cS7G87fKhIfLR3NmWhgK7bS9QkhGAAFwx7QnRMBDRqpeULl/JsnaTz5yU+E1LTPtWBgK7bS9QkhGAAFw1DQNhFQt3lvarFS9Fx27bC5REFEEDAADbIQABAADbIQABAADbIQABAADbIQABAADbIQABAADbIQABAADbIQABAADbIQABAADbIQABAADbIQABAADbIQABAADbIQABAADbIQABAADb8S/qFQCAK5GeliZxcXGWtqnt/ZX+l6VtAp4UHx8vSUlJlrcbFhYmERERUpQIQACKndTziXIw9g8ZPmaCBAUFWdbuxQspcuTocamanm5Zm4Anw8+DAx6WM+dSLG+7fGiIfDR3dpGGIAIQgGInPfWCZPr4S/hNPaVCVIxl7Z46sFviDs+RjL8IQCj+kpKSTPiJaNVLSpevZFm7yWdOSvymJaZ9AhAAXIGQchESVrGKZe2dP33CsrYAq5QuX8nS/0cqXooeRdAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2ikUAmj59ulSrVk2Cg4OlZcuWsmXLlnyXX7x4sdStW9cs36hRI/n6668tW1cAAOD9vD4ALVq0SEaMGCHjx4+XHTt2SJMmTaRz585y6tQpl8tv3LhRevfuLYMGDZKdO3dK9+7dzbR7927L1x0AAHgnrw9AU6ZMkcGDB8uAAQOkfv36MnPmTAkJCZE5c+a4XH7q1Klyxx13yMiRI6VevXoyadIkad68uUybNs3ydQcAAN7JqwNQWlqabN++XW677baseb6+vub2pk2bXD5G52dfXmmPUV7LAwAA+/EXL5aQkCAZGRlSqVKlHPP19p49e1w+5sSJEy6X1/l5SU1NNZNTYmKi+XnmzBlJT08Xd9Ln9vUVST4VJ5KWIla5eOaYBAb6y8X4I5LkZ1mztGuR5LPx4sjIkF9++SXr82uFI0eOiCMzk89zCX1/UbIdKaL/v/p51v2gfpZPnz7t1uc+d+6c+elwOC6/sMOLHT16VLfAsXHjxhzzR44c6bjxxhtdPiYgIMAxf/78HPOmT5/uqFixYp7tjB8/3rTDxMTExMTEJMV+Onz48GUzhlf3AIWHh4ufn5+cPHkyx3y9HRkZ6fIxOr8wy6vRo0ebQmunzMxM0/tToUIF8fHxEW+RlJQk0dHRcvjwYQkLC6Nd2qVd2qVd2rVVu5ejPT/aCxQVFXXZZb06AAUGBkqLFi3k+++/N0dyOcOJ3h42bJjLx7Rq1crcP3z48Kx5q1atMvPzEhQUZKbsypYtK95KP2xF8YGjXdqlXdqlXdr1hnbzU6ZMmQIt59UBSGnPTL9+/eT666+XG2+8Ud566y1JTk42R4Wpvn37yrXXXiuTJ082t5944glp27atvPHGG9K1a1dZuHChbNu2Td57770i3hIAAOAtvD4A3XfffRIfHy/PPfecKWRu2rSprFixIqvQ+dChQ+bIMKfWrVvL/PnzZezYsTJmzBipVauWLF++XBo2bFiEWwEAALyJ1wcgpcNdeQ15rV279pJ599xzj5lKGh2m0xNC5h6uo13apV3apV3atUO77uSjldBufUYAAAAv59UnQgQAAPAEAhAAALAdAhAAALAdAhAAALAdAlAxMX36dKlWrZoEBwdLy5YtZcuWLR5vc/369dKtWzdzRk09I7aeTsDT9HxON9xwg4SGhkrFihXNCTD37t3r8XZnzJghjRs3zjqpl54485tvvhGrvfzyy+a1zn4iT0+ZMGGCaSv7VLduXbHC0aNH5cEHHzRnWy9VqpQ0atTInK/Lk/T/T+7t1Wno0KEebVevZzhu3DipXr262daaNWvKpEmTCnatoqukZ8TVz1JMTIxpW08TsnXrVku/J3Q79TQmlStXNuugF6vet2+fx9tdunSpdOrUKeuM/rt27brqNi/Xrl47ctSoUebzXLp0abOMnqvu2LFjHm3X+f9Z//9qu+XKlTOv83/+8x+Pt5vdY489ZpbR8/UVBwSgYmDRokXmhJB6yOGOHTukSZMm5gr3p06d8mi7esJJbUvDl1XWrVtndkibN282Z/DWLxT9EtN18aQqVaqY8LF9+3azI+7QoYPcfffd5uKTVtEd07vvvmuCmFUaNGggx48fz5o2bNjg8Tb//PNPufnmmyUgIMCEzF9//dWcuFS/tD39+mbfVv18KU+fMuOVV14xAXvatGny22+/mduvvvqqvPPOO+JpDz/8sNnODz/8UH7++Wfzf0l3jBpArfqe0G19++23ZebMmWaHrDto/f66ePGiR9vV+9u0aWNeb3fKr92UlBTzHa2BV39qCNM/4O666y6Ptqtq165tPmP6Puv/Yw38+n7refQ82a7TsmXLzPd2QS5B4TUKemFSFB298OvQoUOzbmdkZDiioqIckydPtmwd9KOybNkyh9VOnTpl2l63bp3lbZcrV84xe/ZsS9o6d+6co1atWo5Vq1Y52rZt63jiiSc83qZeBLhJkyYOq40aNcrRpk0bR1HT17hmzZqOzMxMj7bTtWtXx8CBA3PM69mzp+OBBx7waLspKSkOPz8/x5dffpljfvPmzR3PPvusJd8T+tpGRkY6Xnvttax5Z8+edQQFBTkWLFjgsXazi42NNffv3LnTbe0VpF2nLVu2mOXi4uIsbTcxMdEs991333m83SNHjjiuvfZax+7dux0xMTGON99801Ec0APk5dLS0kyvhP7V5qRnvtbbmzZtkpIuMTHR/CxfvrxlbeqQhV5CRf/yye8acu6kvV566Zbs77MVdChC/2KrUaOGPPDAA+bM6p72+eefm0vbaM+LDnM2a9ZMZs2aJVb/v/roo49k4MCBHr/gsQ476fUJf//9d3P7xx9/NH+hd+nSxaPt/vXXX+azrMPm2ekwlBU9fSo2NtacwT/751qv06TD+Hb4/nJ+h+lnzMrrS+rnWy//pK+19t54UmZmpjz00EMycuRI06NcnBSLM0HbWUJCgvkSc176w0lv79mzR0oy/Y+l9Qs6XGLFpUy061gDj3bNX3PNNaZLt379+h5vV8OWdpe7uzbjcnQnNG/ePKlTp44ZEpo4caLccsstsnv3blOD5Sl//PGHGRLSYV29XI1u9+OPP24ufqzX/bOC1jGcPXtW+vfv7/G2nnnmGXPlbK3P8PPzM/+fX3zxRRM4PUnfQ/08a71RvXr1zHfGggULTPC47rrrxAoafpSr7y/nfSWZfpdoTVDv3r0tuWDol19+Kffff78ZitOaKx3+DA8P92ibr7zyivj7+5v/w8UNAQheS3tFdGds1V+rGgS0UFL/Yvv000/NzlhrkjwZgg4fPmwu4KtfVLn/Uve07D0QWnekgUiLZT/55BMZNGiQR4Ot9gC99NJL5rb2AOn7rDUiVgWg999/32y/FfUK+np+/PHH5hqF+heyfsY02Gvbnt5erf3RXi69YLSGr+bNm5udsfYqw7O0fvHee+81ReAa+K3Qvn178/nSP5y1V1Xb17or7Wn1hO3bt8vUqVPNH3Ce7kn1BIbAvJymd/3iOnnyZI75ejsyMlJKKr32m/41s2bNGlOgbAXtgdC/jFu0aGGORtOuY/3P7Un6BaLF7Lpj0r+idNLQpUWj+rv2FlhFu+i1kHL//v0ebUf/Ms0dKrWHworhNxUXFyffffedKRC2gg4NaC+Q/mWuRwfpcMGTTz5pPmOepkec6efp/PnzJmzr0aO6Y9YhTys4v6Ps9v3lDD/6WdM/bqzo/VFaYK7fYTfddJMJ+fodoj895YcffjDfX1WrVs36/tJtfuqpp0wRtrcjAHk53SnrDllrCLL/Ba23rapPsZL+taThR4efVq9ebQ4dLir6Oqempnq0jY4dO5qhN/2rzTlp74gOj+jvGn6tojvJAwcOmIDiSTqkmfvUBlofo71PVpg7d675i1hrrqygwxFat5edvq/6+bKK7hj1fdUj8FauXGmOcLSC/v/VoJP9+0uHA7VXoiR+f2UPP1pfp0FbD8Mvqd9hDz30kPz00085vr+0Z1NDv37OvB1DYMWA1kpoV7nuGG+88UZzjgUt0B0wYIDHd4jZewO0oFE/4FqQrInfU8NeOlTw2WefmRoGZ52AFvNp8aanjB492gyJ6HbpuVN0HdauXevx/8S6jbnrm3RnpV+anq57evrpp835PTR46HlK9DQLumPWIRJP0t4PLQzWITDdUWivhBZs6mTFDkEDkP5/0r9WraCvsdb86GdLh8B27twpU6ZMMUNTnqafX/2jQod39f+y7pi0Fsmd3x2X+57Q4b4XXnhBatWqZQKRHiKuO0k9x5cn2z1z5ozpVXSeg8cZujWQXU3vU37tasj8+9//boaEtAdbe3Cd32F6v/5B64l29ftCP2N6uL2ugw6B6WHrerqDqz3Nw/nLvM65A56e3kJfX/3Meb2iPgwNBfPOO+84qlat6ggMDDSHxW/evNnjba5Zs8Yc9ph76tevn8fadNWeTnPnznV4kh6mrIdv6usbERHh6Nixo+Pbb791FAWrDoO/7777HJUrVzbbrIew6u39+/c7rPDFF184GjZsaA6Hrlu3ruO9996zpN2VK1eaz9PevXsdVklKSjLvp/7/DQ4OdtSoUcMchp6amurxthctWmTa0/dYD0fX02noYehWfk/oofDjxo1zVKpUybzf+n/LHa//5drV7wxX9+vpHzzVrvOQe1eTPs5T7V64cMHRo0cPc3oUfa/1//Vdd91lDsG3ej8QU4wOg/fRf4o6hAEAAFiJGiAAAGA7BCAAAGA7BCAAAGA7BCAAAGA7BCAAAGA7BCAAAGA7BCAAAGA7BCAAAGA7BCAAltKLkOplEUJCQqRcuXLmmlAfffRRUa8WAJvhWmAALKXXDpo9e7a5arVeKHTTpk3y2GOPmWsO6U8AsAI9QAAs9corr0jbtm3l2muvNT1Bffv2lU6dOsn69evN/dWqVTMX/M2uf//+OS6euWLFCmnTpo2ULVvWBKq//e1v5kr2TvPmzTP3ZdeuXTtzYU4nvUq2XhBW10MvQNuyZUtzAdz8nuPgwYPi4+NjLgapdHm9ffbs2RxXyNZ5y5cvz5qnIU97uq655hpzn05Nmza9qtcRwNUhAAEoMnopwu3bt8vGjRvljjvuKPDjkpOTZcSIEbJt2zb5/vvvxdfXV3r06GGu9l5Qw4YNM8Fk4cKF8tNPP5mrZus67Nu37wq3Rsy2fP7555fM1yuER0dHmyvBHz9+XJ566qkrbgOAexCAAFhOe0e0NyQwMFBuuOEGefTRR01PUEH16tVLevbsaYbRtCdlzpw58vPPP8uvv/5q7i9VqpRcvHgxz8cfOnRI5s6dK4sXL5ZbbrlFatasaXqDtFdJ518pDWUjR47MMe/UqVNy7Ngx0/ukPV6RkZFm2wEULQIQAMvdfvvtZhhp69atMmPGDJk6darMnDkz6/5Ro0aZkOCcPv744xyP116a3r17S40aNSQsLMwMmzmDjWrQoIEZ4lqyZInL9jUsZWRkSO3atXO0s27duhxDaYmJiTnu1+fNL9T98ccfl/TulC9fXsqUKSOffPKJpKenX+ErBsDdKIIGYDmtudHeG6U9OPHx8fL6669nFUFrL4rW/WQPRBpYnLp16yYxMTEya9YsiYqKMkNfDRs2lLS0NHO//q6P0WGt4OBgM0R24cKFrLobLbj28/MzQ1b6M7vsvTOhoaGyY8eOrNtHjx41tUS5abD55z//KS+++KLpfcrO399fPvzwQxkyZIhMmzbNrI+uZ/369a/6dQRw5QhAALyiFih7/U54eHhWQHIGEWeh8enTp2Xv3r0m/OjwldqwYcMlz/nyyy/LmDFjzBCUeuCBB7Lua9asmQlUep/zOVzR4JR9PTTMuKK9WBqctADaFQ1sGoI0KL322mvy9ttvZxV9AygaBCAAlklKSjLnAXrkkUekTp06plfmhx9+MKFg7NixBXoOPXeQHvn13nvvSeXKlc2w1zPPPONyWR0e00ll75nRoS8NRFp39MYbb5hApL1QWlDduHFj6dq1a6G269VXX5UvvvjCHN3lypQpU7KG/HQ4TIfFABQtAhAAy+jwj4YXrZPRQ8p1+KlRo0by/vvvm+GqgtBeGT1y6/HHHzdDXRqktEfF1dBUfrTY+YUXXjDrokNb2ut00003mUPqC6t9+/ZmckUD3sSJE00vlYYfAN7Bx6F9zwAAADbCUWAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2CEAAAMB2/hfzsZL7I1wb0gAAAABJRU5ErkJggg==",
|
||
"text/plain": [
|
||
"<Figure size 640x480 with 1 Axes>"
|
||
]
|
||
},
|
||
"metadata": {},
|
||
"output_type": "display_data"
|
||
}
|
||
],
|
||
"source": [
|
||
"plt.hist(data, bins=np.arange(min(data)-0.5, max(data)+1.5, 1), edgecolor='black', alpha=0.7)\n",
|
||
"plt.title(\"Гистограмма частот\")\n",
|
||
"plt.xlabel(\"Значения\")\n",
|
||
"plt.ylabel(\"Частота\")\n",
|
||
"plt.xticks(np.arange(min(data), max(data)+1))\n",
|
||
"plt.grid(axis='y')\n",
|
||
"plt.show()"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "44f7e836",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Пункт b)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "c32cd292",
|
||
"metadata": {},
|
||
"source": [
|
||
"### (i) Выборочное среднее (математическое ожидание)\n",
|
||
"Выборочное среднее — оценка теоретического математического ожидания.\n",
|
||
"$$\n",
|
||
"\\bar{X} = \\frac{1}{n} \\sum_{i=1}^{n} X_i.\n",
|
||
"$$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 6,
|
||
"id": "ead66cb6",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Выборочное среднее: 1.96\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"import numpy as np\n",
|
||
"mean = np.mean(data)\n",
|
||
"print(f\"Выборочное среднее: {mean:.2f}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "83c9665b",
|
||
"metadata": {},
|
||
"source": [
|
||
"### (ii) Выборочная дисперсия\n",
|
||
"Несмещённая оценка дисперсии:\n",
|
||
"$$\n",
|
||
"s^2 = \\frac{1}{n-1} \\sum_{i=1}^{n}(X_i-\\bar{X})^2.\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Смещенная оценка дисперсии:\n",
|
||
"$$\n",
|
||
"s^2_{\\text{смещенная}} = \\frac{1}{n} \\sum_{i=1}^{n}(X_i - \\bar{X})^2\n",
|
||
"$$\n",
|
||
"\n",
|
||
"где:\n",
|
||
"- $ n $ — общее количество наблюдений,\n",
|
||
"- $X_i$ — каждое отдельное наблюдение,\n",
|
||
"- $\\bar{X}$ — среднее значение выборки."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 4,
|
||
"id": "a24ea7eb",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Несмещённая оценка дисперсии: 7.67\n",
|
||
"Смещённая оценка дисперсии: 7.52\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"variance = np.var(data, ddof=1)\n",
|
||
"print(f\"Несмещённая оценка дисперсии: {variance:.2f}\")\n",
|
||
"print(f\"Смещённая оценка дисперсии: {(np.var(data, ddof=0)):.2f}\")\n",
|
||
"# print(sum((x - mean) ** 2 for x in data) / (n - 1))"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "bd8ee128",
|
||
"metadata": {},
|
||
"source": [
|
||
"### (iii) Медиана\n",
|
||
"Значение, разделяющее выборку на две равные части."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 8,
|
||
"id": "e8490052",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Медиана: 1.0\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"median = np.median(data)\n",
|
||
"print(f\"Медиана: {median}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "34384b8f",
|
||
"metadata": {},
|
||
"source": [
|
||
"### (iv) Ассиметрия\n",
|
||
"$$\n",
|
||
"Skewness = \\frac{\\frac{1}{n}\\sum_{i=1}^{n}(X_i-\\bar{X})^3}{s^3}.\n",
|
||
"$$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 9,
|
||
"id": "cc21a5b6",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Асимметрия: 2.25\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"from scipy.stats import skew\n",
|
||
"skewness = skew(data)\n",
|
||
"print(f\"Асимметрия: {skewness:.2f}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "ddd4b8a7",
|
||
"metadata": {},
|
||
"source": [
|
||
"### (v) Эксцесс\n",
|
||
"$$\n",
|
||
"Kurtosis = \\frac{\\frac{1}{n}\\sum_{i=1}^{n}(X_i-\\bar{X})^4}{s^4} - 3.\n",
|
||
"$$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 10,
|
||
"id": "118d475e",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Эксцесс: 5.92\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"from scipy.stats import kurtosis\n",
|
||
"excess_kurtosis = kurtosis(data)\n",
|
||
"print(f\"Эксцесс: {excess_kurtosis:.2f}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "93fd7cc5",
|
||
"metadata": {},
|
||
"source": [
|
||
"### (vi) Вероятность $P(X \\in [0.00, 2.49])$\n",
|
||
"$$\n",
|
||
"P(X \\in [a, b]) = \\frac{\\text{число элементов выборки} \\in [a, b]}{n}.\n",
|
||
"$$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 11,
|
||
"id": "08ea631c",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"P(X ∈ [0.0, 2.49]): 0.74\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"count = np.sum((data >= a) & (data <= b))\n",
|
||
"probability = count / len(data)\n",
|
||
"print(f\"P(X ∈ [{a}, {b}]): {probability:.2f}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "26424ded",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Пункт c)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "f6b509ff",
|
||
"metadata": {},
|
||
"source": [
|
||
"### 1. Оценка максимального правдоподобия (ОМП)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "c40e8461",
|
||
"metadata": {},
|
||
"source": [
|
||
"Функция правдоподобия для Пуассона:\n",
|
||
"$$\n",
|
||
"L(λ) = \\prod_{i=1}^{n}\\frac{λ^{X_i}e^{-λ}}{X_i!}.\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Логарифмируя, получаем:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\ln L(\\lambda) = \\sum_{i=1}^{n} \\left( X_i \\ln \\lambda - \\lambda - \\ln X_i! \\right).\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Дифференцируя по $\\lambda$, приравнивая к нулю:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\frac{d}{d\\lambda} \\ln L(\\lambda) = \\frac{1}{\\lambda} \\sum_{i=1}^{n} X_i - n = 0 \n",
|
||
"\\Longrightarrow \\hat{\\lambda}_{\\text{ОМП}} = \\frac{1}{n} \\sum_{i=1}^{n} X_i = \\bar{X}.\n",
|
||
"$$\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 12,
|
||
"id": "7fa556a6",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"ОМП для λ: 1.96\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"lambda_ml = np.mean(data)\n",
|
||
"print(f\"ОМП для λ: {lambda_ml:.2f}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "5a5e2f27",
|
||
"metadata": {},
|
||
"source": [
|
||
"\n",
|
||
"**Смещение ОМП:** \n",
|
||
"В случае распределения Пуассона оценка максимального правдоподобия (ОМП) параметра $\\lambda$ совпадает с выборочным средним:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\hat{\\lambda}_{\\text{ОМП}} = \\bar{x} = \\frac{1}{n} \\sum_{i=1}^{n} x_i.\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Найдём математическое ожидание этой оценки:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\mathbb{E}[\\hat{\\lambda}_{\\text{ОМП}}] = \\mathbb{E} \\left[ \\frac{1}{n} \\sum_{i=1}^{n} x_i \\right] = \\frac{1}{n} \\sum_{i=1}^{n} \\mathbb{E}[x_i].\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Так как для распределения Пуассона $\\mathbb{E}[x_i] = \\lambda$, то:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\mathbb{E}[\\hat{\\lambda}_{\\text{ОМП}}] = \\frac{1}{n} \\cdot n \\lambda = \\lambda.\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Отсюда следует:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\text{Смещение}(\\hat{\\lambda}_{\\text{ОМП}}) = \\lambda - \\lambda = 0.\n",
|
||
"$$\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "545f29e7",
|
||
"metadata": {},
|
||
"source": [
|
||
"### 2. Оценка по методу моментов (ОММ)\n",
|
||
"Приравниваем теоретическое математическое ожидание к выборочному:\n",
|
||
"$$\n",
|
||
"E[X]=λ \\Longrightarrow \\hat{λ}_{\\text{MM}} = \\bar{X}. \\\n",
|
||
"$$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 13,
|
||
"id": "96484e1c",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"ОММ для λ: 1.96\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"lambda_mm = np.mean(data)\n",
|
||
"print(f\"ОММ для λ: {lambda_mm:.2f}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "525cff2b",
|
||
"metadata": {},
|
||
"source": [
|
||
"\n",
|
||
"**Смещение ОММ:** \n",
|
||
"Метод моментов приводит к той же оценке:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\hat{\\lambda}_{\\text{ММ}} = \\bar{x}.\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Математическое ожидание:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\mathbb{E}[\\hat{\\lambda}_{\\text{ММ}}] = \\lambda \\\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Смещение этой оценки:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\text{Смещение}(\\hat{\\lambda}_{\\text{ММ}}) = \\lambda - \\lambda = 0.\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Таким образом, обе оценки ($\\hat{\\lambda}_{\\text{ОМП}}$ и $\\hat{\\lambda}_{\\text{ММ}}$) являются несмещёнными.\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "289e0726",
|
||
"metadata": {},
|
||
"source": [
|
||
"# d) Aсимптотический доверительный интервал уровня значимости α1=0.02 для параметра λ на базе оценки максимального правдоподобия\n",
|
||
"\n",
|
||
"## Шаги построения\n",
|
||
"\n",
|
||
"### 1. Оценка $\\hat{\\lambda}$\n",
|
||
"ОМП параметра $\\lambda$ равна выборочному среднему:\n",
|
||
"$$ \\hat{\\lambda} = \\bar{x} $$\n",
|
||
"\n",
|
||
"### 2. Стандартная ошибка\n",
|
||
"Для распределения Пуассона дисперсия равна $\\lambda$:\n",
|
||
"$$ SE = \\sqrt{\\frac{\\hat{\\lambda}}{n}} $$\n",
|
||
"\n",
|
||
"### 3. Квантиль нормального распределения\n",
|
||
"Для уровня значимости $\\alpha_{1} = 0.02$:\n",
|
||
"$$ z_{1-\\alpha/2} = z_{0.99} $$\n",
|
||
"\n",
|
||
"### 4. Границы интервала\n",
|
||
"$$ \\hat{\\lambda} \\pm z_{0.99} \\cdot SE $$\n",
|
||
"\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 16,
|
||
"id": "7f3db200",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"z = 2.326\n",
|
||
"se = 0.198\n",
|
||
"Доверительный интервал (98%): (1.499, 2.421)\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"import numpy as np\n",
|
||
"from scipy.stats import norm\n",
|
||
"\n",
|
||
"z = norm.ppf(1 - alpha/2)\n",
|
||
"se = np.sqrt(lambda_ml / len(data))\n",
|
||
"lower = lambda_ml - z * se\n",
|
||
"upper = lambda_ml + z * se\n",
|
||
"\n",
|
||
"print(f\"z = {z:.3f}\")\n",
|
||
"print(f\"se = {se:.3f}\")\n",
|
||
"print(f\"Доверительный интервал (98%): ({lower:.3f}, {upper:.3f})\")\n",
|
||
"\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "4604ecf9",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Пункт e) Критерий $\\chi^2$ для проверки гипотезы согласия с распределением Пуассона ($λ0 = 2.00$)\n",
|
||
"Критерий $\\chi^2$ проверяет, насколько эмпирические частоты $O_i$ соответствуют теоретическим частотам $E_i$ при заданном распределении.\n",
|
||
"\n",
|
||
"1. **Расчёт наблюдаемых и теоретических частот:** \n",
|
||
" $O_i$ - наблюдаемые частоты для каждого интервала,\n",
|
||
"\n",
|
||
" $$\n",
|
||
" E_i = n \\cdot P(X = k\\ |\\ λ = λ_0),\n",
|
||
" $$\n",
|
||
" где $P(X=k)$ — вероятность по распределению Пуассона.\n",
|
||
"\n",
|
||
"2. **Группировка данных:** Объединить значения так, чтобы $E_i \\geq 5$.\n",
|
||
"\n",
|
||
"3. **Статистика $\\chi^2$:**\n",
|
||
" $$\n",
|
||
" \\chi^2 = \\sum_{i=1}^{k}\\frac{(O_i - E_i)^2}{E_i}.\n",
|
||
" $$\n",
|
||
"4. **Степени свободы:**\n",
|
||
" $$\n",
|
||
" df = k - 1 - m,\n",
|
||
" $$\n",
|
||
" где $k$ — число категорий, $m=0$. \n",
|
||
"\n",
|
||
"**Критическое значение:** Сравнение с $χ_{\\text{крит}}^2(df, α)$. \n",
|
||
"**p-значение:** Вероятность $P(χ^2 \\geq χ_{\\text{набл}}^2$)."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 17,
|
||
"id": "d881725f",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Таблица до группировки категорий:\n",
|
||
" Значение k Наблюдаемая частота (O_i) Теоретическая вероятность (P(X=k)) \\\n",
|
||
"0 0 19 0.1353 \n",
|
||
"1 1 11 0.2707 \n",
|
||
"2 2 7 0.2707 \n",
|
||
"3 3 4 0.1804 \n",
|
||
"4 4 3 0.0902 \n",
|
||
"5 6 2 0.0120 \n",
|
||
"6 7 1 0.0034 \n",
|
||
"7 8 2 0.0009 \n",
|
||
"8 14 1 0.0000 \n",
|
||
"\n",
|
||
" Теоретическая частота (E_i) \n",
|
||
"0 6.767 \n",
|
||
"1 13.534 \n",
|
||
"2 13.534 \n",
|
||
"3 9.022 \n",
|
||
"4 4.511 \n",
|
||
"5 0.601 \n",
|
||
"6 0.172 \n",
|
||
"7 0.043 \n",
|
||
"8 0.000 \n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"import pandas as pd\n",
|
||
"from scipy.stats import poisson, chi2\n",
|
||
"\n",
|
||
"# Теоретические вероятности для каждого k\n",
|
||
"probs_individual = [poisson.pmf(k, lambda0) for k in unique_values]\n",
|
||
"\n",
|
||
"# Теоретические частоты\n",
|
||
"expected_individual = np.array(probs_individual) * n\n",
|
||
"\n",
|
||
"df_individual = pd.DataFrame({\n",
|
||
" \"Значение k\": unique_values,\n",
|
||
" \"Наблюдаемая частота (O_i)\": counts,\n",
|
||
" \"Теоретическая вероятность (P(X=k))\": np.round(probs_individual, 4),\n",
|
||
" \"Теоретическая частота (E_i)\": np.round(expected_individual, 3)\n",
|
||
"})\n",
|
||
"\n",
|
||
"print(\"Таблица до группировки категорий:\")\n",
|
||
"print(df_individual)\n",
|
||
"\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "c669819b",
|
||
"metadata": {},
|
||
"source": [
|
||
"### Интерпретация\n",
|
||
"- **Наблюдаемые частоты** $O_i$ — количество раз, когда значение $k$ встречается в выборке.\n",
|
||
"- **Теоретическая вероятность** $P(X=k)$ — вероятность по распределению Пуассона с $λ=2.0$.\n",
|
||
"- **Теоретическая частота** $E_i$ — ожидаемое количество значений $k$ при условии, что данные следуют распределению Пуассона ($E_i = n \\cdot P(X = k)$).\n",
|
||
"\n",
|
||
"После группировки категорий (чтобы $E_i ≥ 5$) таблица принимает вид:"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 33,
|
||
"id": "74f3d6a5",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"\n",
|
||
"Таблица после группировки категорий:\n",
|
||
" Группа Наблюдаемая частота (O_i) Теоретическая вероятность \\\n",
|
||
"0 0 19 0.1353 \n",
|
||
"1 1 11 0.2707 \n",
|
||
"2 2 7 0.2707 \n",
|
||
"3 3 4 0.1804 \n",
|
||
"4 4,5,6,7,8 9 0.1429 \n",
|
||
"\n",
|
||
" Теоретическая частота (E_i) \n",
|
||
"0 6.767 \n",
|
||
"1 13.534 \n",
|
||
"2 13.534 \n",
|
||
"3 9.022 \n",
|
||
"4 7.144 \n",
|
||
"\n",
|
||
"χ² наблюдаемое: 29.022\n",
|
||
"Критическое значение (α=0.02): 11.668\n",
|
||
"p-значение: 0.0000077\n",
|
||
"Отвергаем гипотезу на уровне 0.02\n",
|
||
"Наибольший уровень значимости, на котором ещё нет оснований отвергнуть гипотезу: 0.0000077\n",
|
||
"Это означает, что гипотеза отвергается на любом уровне значимости α ≥ 0.0000077\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"from scipy.stats import chi2\n",
|
||
"# Группировка категорий (для E_i ≥ 5)\n",
|
||
"groups = [\n",
|
||
" [0], # Группа 1: k=0\n",
|
||
" [1], # Группа 2: k=1\n",
|
||
" [2], # Группа 3: k=2\n",
|
||
" [3], # Группа 4: k=3\n",
|
||
" [4, 5, 6, 7, 8] # Группа 5: k=4,5,6,7,8\n",
|
||
"]\n",
|
||
"\n",
|
||
"# Расчёт наблюдаемых частот по группам\n",
|
||
"# observed_grouped = np.array([19, 11, 7, 4, 3+2+1+2+1]) # ??\n",
|
||
"observed_grouped = np.array([np.sum(data==k) for k in [0,1,2,3]] + [np.sum(data>=4)])\n",
|
||
"\n",
|
||
"# Расчёт теоретических вероятностей по группам\n",
|
||
"probs_grouped = [\n",
|
||
" poisson.pmf(0, lambda0),\n",
|
||
" poisson.pmf(1, lambda0),\n",
|
||
" poisson.pmf(2, lambda0),\n",
|
||
" poisson.pmf(3, lambda0),\n",
|
||
" 1 - poisson.cdf(3, lambda0) # sum(poisson.pmf(k, lambda0) for k in groups[4])\n",
|
||
"]\n",
|
||
"\n",
|
||
"# Теоретические частоты\n",
|
||
"expected_grouped = np.array(probs_grouped) * n\n",
|
||
"\n",
|
||
"# Создание таблицы после группировки\n",
|
||
"df_grouped = pd.DataFrame({\n",
|
||
" \"Группа\": [\"0\", \"1\", \"2\", \"3\", \"4,5,6,7,8\"],\n",
|
||
" \"Наблюдаемая частота (O_i)\": observed_grouped,\n",
|
||
" \"Теоретическая вероятность\": np.round(probs_grouped, 4),\n",
|
||
" \"Теоретическая частота (E_i)\": np.round(expected_grouped, 3)\n",
|
||
"})\n",
|
||
"\n",
|
||
"print(\"\\nТаблица после группировки категорий:\")\n",
|
||
"print(df_grouped)\n",
|
||
"\n",
|
||
"# Статистика χ²\n",
|
||
"chi2_stat = np.sum((observed_grouped - expected_grouped)**2 / expected_grouped)\n",
|
||
"\n",
|
||
"# Степени свободы\n",
|
||
"df = 5 - 1 - 0 # 4\n",
|
||
"\n",
|
||
"# Критическое значение и p-значение\n",
|
||
"chi2_crit = chi2.ppf(1 - alpha, df)\n",
|
||
"p_value = 1 - chi2.cdf(chi2_stat, df)\n",
|
||
"\n",
|
||
"print(f\"\\nχ² наблюдаемое: {chi2_stat:.3f}\")\n",
|
||
"print(f\"Критическое значение (α=0.02): {chi2_crit:.3f}\")\n",
|
||
"print(f\"p-значение: {p_value:.7f}\")\n",
|
||
"\n",
|
||
"if chi2_stat > chi2_crit:\n",
|
||
" print(\"Отвергаем гипотезу на уровне 0.02\")\n",
|
||
"else:\n",
|
||
" print(\"Нет оснований отвергнуть гипотезу на уровне 0.02\")\n",
|
||
"print(f\"\"\"Наибольший уровень значимости, на котором ещё нет оснований отвергнуть гипотезу: {p_value:.7f}\n",
|
||
"Это означает, что гипотеза отвергается на любом уровне значимости α ≥ {p_value:.7f}\"\"\")\n",
|
||
"\n",
|
||
"# observed = np.array([np.sum(data==k) for k in [0,1,2,3]] + [np.sum(data>=4)])\n",
|
||
"# expected = np.array([poisson.pmf(k,2)*n for k in [0,1,2,3]] + [n*(1 - poisson.cdf(3,2))])\n",
|
||
"# chi2_stat = np.sum((observed - expected)**2 / expected)\n",
|
||
"# df = 4\n",
|
||
"# crit = chi2.ppf(1-0.02, df)\n",
|
||
"# p_val = 1 - chi2.cdf(chi2_stat, df)\n",
|
||
"# print(f\"\\ne) χ²: {chi2_stat:.2f}, крит: {crit:.2f}, p-value: {p_val:.4f}\")\n",
|
||
"\n",
|
||
"\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 37,
|
||
"id": "b4be90df",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
" Group Lower Upper O_i P_i E_i O_i - E_i \\\n",
|
||
"0 0 -inf 0.0 19 0.1353 6.767 12.233 \n",
|
||
"1 1 1.0 1.0 11 0.2707 13.534 -2.534 \n",
|
||
"2 2 2.0 2.0 7 0.2707 13.534 -6.534 \n",
|
||
"3 3 3.0 3.0 4 0.1804 9.022 -5.022 \n",
|
||
"4 4, 5, 6, 7, 8 4.0 inf 9 0.1429 7.144 1.856 \n",
|
||
"\n",
|
||
" (O_i - E_i)^2 / E_i \n",
|
||
"0 22.1157 \n",
|
||
"1 0.4743 \n",
|
||
"2 3.1542 \n",
|
||
"3 2.7957 \n",
|
||
"4 0.4823 \n",
|
||
"\n",
|
||
"χ² наблюдаемое: 29.022\n",
|
||
"Критическое значение (α=0.02): 11.668\n",
|
||
"p-значение: 0.0000077\n",
|
||
"Отвергаем H₀\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"import numpy as np\n",
|
||
"import pandas as pd\n",
|
||
"from scipy.stats import poisson, chi2\n",
|
||
"\n",
|
||
"# Пусть заданы:\n",
|
||
"# counts — частоты\n",
|
||
"# unique_values — уникальные значения\n",
|
||
"# n — общее число наблюдений\n",
|
||
"# lambda0 — параметр Пуассона\n",
|
||
"# alpha — уровень значимости - вычислить его\n",
|
||
"\n",
|
||
"# Группы значений\n",
|
||
"groups = [\n",
|
||
" [0], # Группа 1\n",
|
||
" [1], # Группа 2\n",
|
||
" [2], # Группа 3\n",
|
||
" [3], # Группа 4\n",
|
||
" [4, 5, 6, 7, 8] # Группа 5\n",
|
||
"]\n",
|
||
"\n",
|
||
"observed_grouped = np.array([np.sum(data==k) for k in [0,1,2,3]] + [np.sum(data>=4)])\n",
|
||
"probs_grouped = [\n",
|
||
" poisson.pmf(0, lambda0),\n",
|
||
" poisson.pmf(1, lambda0),\n",
|
||
" poisson.pmf(2, lambda0),\n",
|
||
" poisson.pmf(3, lambda0),\n",
|
||
" 1 - poisson.cdf(3, lambda0) # sum(poisson.pmf(k, lambda0) for k in groups[4])\n",
|
||
"]\n",
|
||
"\n",
|
||
"# Теоретические частоты\n",
|
||
"expected_grouped = np.array(probs_grouped) * n\n",
|
||
"lower_bounds = []\n",
|
||
"upper_bounds = []\n",
|
||
"\n",
|
||
"for i, group in enumerate(groups):\n",
|
||
" # obs = sum(counts[np.where(unique_values == k)[0][0]] for k in group if k in unique_values)\n",
|
||
" # prob = sum(poisson.pmf(k, lambda0) for k in group)\n",
|
||
" # exp = prob * n\n",
|
||
" # observed_grouped.append(obs)\n",
|
||
" # expected_grouped.append(exp)\n",
|
||
"\n",
|
||
" # Нижняя и верхняя границы\n",
|
||
" lower = -np.inf if i == 0 else min(group)\n",
|
||
" upper = np.inf if i == len(groups) - 1 else max(group)\n",
|
||
" lower_bounds.append(lower)\n",
|
||
" upper_bounds.append(upper)\n",
|
||
"\n",
|
||
"# Разности и вклад в статистику\n",
|
||
"diff = np.array(observed_grouped) - np.array(expected_grouped)\n",
|
||
"chi2_terms = diff**2 / expected_grouped\n",
|
||
"\n",
|
||
"# Таблица\n",
|
||
"df_final = pd.DataFrame({\n",
|
||
" \"Group\": [\", \".join(map(str, g)) for g in groups],\n",
|
||
" \"Lower\": lower_bounds,\n",
|
||
" \"Upper\": upper_bounds,\n",
|
||
" \"O_i\": observed_grouped,\n",
|
||
" \"P_i\": np.round(np.array(expected_grouped) / n, 4),\n",
|
||
" \"E_i\": np.round(expected_grouped, 3),\n",
|
||
" \"O_i - E_i\": np.round(diff, 3),\n",
|
||
" \"(O_i - E_i)^2 / E_i\": np.round(chi2_terms, 4)\n",
|
||
"})\n",
|
||
"\n",
|
||
"print(df_final)\n",
|
||
"\n",
|
||
"# Хи-квадрат статистика и p-value\n",
|
||
"chi2_stat = np.sum(chi2_terms)\n",
|
||
"df = len(groups) - 1 # без оценки параметров — простая гипотеза\n",
|
||
"chi2_crit = chi2.ppf(1 - alpha, df)\n",
|
||
"p_value = 1 - chi2.cdf(chi2_stat, df)\n",
|
||
"\n",
|
||
"print(f\"\\nχ² наблюдаемое: {chi2_stat:.3f}\")\n",
|
||
"print(f\"Критическое значение (α={alpha:.2f}): {chi2_crit:.3f}\")\n",
|
||
"print(f\"p-значение: {p_value:.7f}\")\n",
|
||
"print(\"Отвергаем H₀\" if chi2_stat > chi2_crit else \"Нет оснований отвергнуть H₀\")\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "f9ef2691",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Пункт f) Критерий $χ^2$ для проверки сложной гипотезы согласия с распределением Пуассона\n",
|
||
"\n",
|
||
"**Оценка параметра $\\lambda$** \n",
|
||
"Если параметр $\\lambda$ неизвестен, его оценивают по выборке (например, через выборочное среднее): \n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\hat{\\lambda} = \\frac{1}{n} \\sum_{i=1}^n x_i,\n",
|
||
"$$\n",
|
||
"\n",
|
||
"где $x_i$ — значения выборки, $n$ — объем выборки.\n",
|
||
"\n",
|
||
"**Степени свободы** \n",
|
||
"Число степеней свободы для критерия хи-квадрат: \n",
|
||
"\n",
|
||
"$$\n",
|
||
"df = k - 1 - m,\n",
|
||
"$$\n",
|
||
"\n",
|
||
"где: \n",
|
||
"- \\( k \\) — количество интервалов, \n",
|
||
"- \\( m \\) — количество оцененных параметров (в данном случае \\( m = 1 \\), так как оценивается $\\lambda$).\n",
|
||
"\n",
|
||
"**Критическое значение:** Сравнение с $χ_{\\text{крит}}^2(df, α)$. \n",
|
||
"**p-значение:** Вероятность $P(χ^2 \\geq χ_{\\text{набл}}^2$).\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 62,
|
||
"id": "4383629c",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Таблица до группировки категорий:\n",
|
||
" Значение k Наблюдаемая частота (O_i) Теоретическая вероятность (P(X=k)) \\\n",
|
||
"0 0 19 0.1409 \n",
|
||
"1 1 11 0.2761 \n",
|
||
"2 2 7 0.2706 \n",
|
||
"3 3 4 0.1768 \n",
|
||
"4 4 3 0.0866 \n",
|
||
"5 6 2 0.0111 \n",
|
||
"6 7 1 0.0031 \n",
|
||
"7 8 2 0.0008 \n",
|
||
"8 14 1 0.0000 \n",
|
||
"\n",
|
||
" Теоретическая частота (E_i) \n",
|
||
"0 7.043 \n",
|
||
"1 13.804 \n",
|
||
"2 13.528 \n",
|
||
"3 8.838 \n",
|
||
"4 4.331 \n",
|
||
"5 0.555 \n",
|
||
"6 0.155 \n",
|
||
"7 0.038 \n",
|
||
"8 0.000 \n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"import pandas as pd\n",
|
||
"from scipy.stats import poisson\n",
|
||
"\n",
|
||
"# Теоретические вероятности для каждого k\n",
|
||
"probs_individual = [poisson.pmf(k, mean) for k in unique_values]\n",
|
||
"\n",
|
||
"# Теоретические частоты\n",
|
||
"expected_individual = np.array(probs_individual) * n\n",
|
||
"\n",
|
||
"df_individual = pd.DataFrame({\n",
|
||
" \"Значение k\": unique_values,\n",
|
||
" \"Наблюдаемая частота (O_i)\": counts,\n",
|
||
" \"Теоретическая вероятность (P(X=k))\": np.round(probs_individual, 4),\n",
|
||
" \"Теоретическая частота (E_i)\": np.round(expected_individual, 3)\n",
|
||
"})\n",
|
||
"\n",
|
||
"print(\"Таблица до группировки категорий:\")\n",
|
||
"print(df_individual)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "7c23d34d",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# from scipy.stats import chi2\n",
|
||
"# # Группировка категорий (для E_i ≥ 5)\n",
|
||
"# groups = [\n",
|
||
"# [0], # Группа 1: k=0,1,2\n",
|
||
"# [1], # Группа 2: k=3\n",
|
||
"# [2], # Группа 3: k=4\n",
|
||
"# [3], # Группа 4: k=5\n",
|
||
"# [4, 5, 6, 7, 8] # Группа 5: k=6,8,9\n",
|
||
"# ]\n",
|
||
"\n",
|
||
"# # Расчёт наблюдаемых частот по группам\n",
|
||
"# # observed_grouped = np.array([np.sum(data==k) for k in [0,1,2,3]] + [np.sum(data>=4)])\n",
|
||
"# probs_grouped = [\n",
|
||
"# poisson.pmf(0, lambda0),\n",
|
||
"# poisson.pmf(1, lambda0),\n",
|
||
"# poisson.pmf(2, lambda0),\n",
|
||
"# poisson.pmf(3, lambda0),\n",
|
||
"# 1 - poisson.cdf(3, lambda0) # sum(poisson.pmf(k, lambda0) for k in groups[4])\n",
|
||
"# ]\n",
|
||
"\n",
|
||
"# # Расчёт теоретических вероятностей по группам\n",
|
||
"# probs_grouped = [\n",
|
||
"# poisson.pmf(3, mean),\n",
|
||
"# poisson.pmf(3, mean),\n",
|
||
"# poisson.pmf(4, mean),\n",
|
||
"# poisson.pmf(5, mean),\n",
|
||
"# sum(poisson.pmf(k, mean) for k in groups[4])\n",
|
||
"# ]\n",
|
||
"\n",
|
||
"# # Теоретические частоты\n",
|
||
"# expected_grouped = np.array(probs_grouped) * n\n",
|
||
"\n",
|
||
"# # Создание таблицы после группировки\n",
|
||
"# df_grouped = pd.DataFrame({\n",
|
||
"# \"Группа\": [\"0,1,2\", \"3\", \"4\", \"5\", \"6,8,9\"],\n",
|
||
"# \"Наблюдаемая частота (O_i)\": observed_grouped,\n",
|
||
"# \"Теоретическая вероятность\": np.round(probs_grouped, 4),\n",
|
||
"# \"Теоретическая частота (E_i)\": np.round(expected_grouped, 3)\n",
|
||
"# })\n",
|
||
"\n",
|
||
"# print(\"\\nТаблица после группировки категорий:\")\n",
|
||
"# print(df_grouped)\n",
|
||
"\n",
|
||
"# # Статистика χ²\n",
|
||
"# chi2_stat = np.sum((observed_grouped - expected_grouped)**2 / expected_grouped)\n",
|
||
"\n",
|
||
"# # Степени свободы\n",
|
||
"# df = 5 - 1 - 1 # 3\n",
|
||
"\n",
|
||
"# # Критическое значение и p-значение\n",
|
||
"# chi2_crit = chi2.ppf(1 - alpha, df)\n",
|
||
"# p_value = 1 - chi2.cdf(chi2_stat, df)\n",
|
||
"\n",
|
||
"# print(f\"\\nχ² наблюдаемое: {chi2_stat:.3f}\")\n",
|
||
"# print(f\"Критическое значение (α=0.10): {chi2_crit:.3f}\")\n",
|
||
"# print(f\"p-значение: {p_value:.3f}\")\n",
|
||
"\n",
|
||
"# if chi2_stat > chi2_crit:\n",
|
||
"# print(\"Отвергаем гипотезу на уровне 0.10\")\n",
|
||
"# else:\n",
|
||
"# print(\"Нет оснований отвергнуть гипотезу на уровне 0.10\")\n",
|
||
"\n",
|
||
"# lambda_hat = np.mean(data)\n",
|
||
"# expected = np.array([poisson.pmf(k,lambda_hat)*n for k in [0,1,2,3]] + [n*(1 - poisson.cdf(3,lambda_hat))])\n",
|
||
"# chi2_stat = np.sum((observed - expected)**2 / expected)\n",
|
||
"# df = 3\n",
|
||
"# crit = chi2.ppf(1-0.02, df)\n",
|
||
"# p_val = 1 - chi2.cdf(chi2_stat, df)\n",
|
||
"# print(f\"\\nf) χ²: {chi2_stat:.2f}, крит: {crit:.2f}, p-value: {p_val:.4f}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 75,
|
||
"id": "9d39fc45",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# import numpy as np\n",
|
||
"# from scipy.stats import poisson, chi2\n",
|
||
"# groups = [\n",
|
||
"# [0], # Группа 1: k=0,1,2\n",
|
||
"# [1], # Группа 2: k=3\n",
|
||
"# [2], # Группа 3: k=4\n",
|
||
"# [3], # Группа 4: k=5\n",
|
||
"# [4, 5, 6, 7, 8] # Группа 5: k=6,8,9\n",
|
||
"# ]\n",
|
||
"\n",
|
||
"# # Оценка параметра λ\n",
|
||
"# lambda_hat = np.mean(data)\n",
|
||
"# print(f\"Оценка λ: {lambda_hat:.4f}\")\n",
|
||
"\n",
|
||
"# # Разбиение на интервалы (пример)\n",
|
||
"# intervals = [\n",
|
||
"# (-np.inf, 0),\n",
|
||
"# (1, 1),\n",
|
||
"# (2, 2),\n",
|
||
"# (3, 3),\n",
|
||
"# (4, np.inf)\n",
|
||
"# ]\n",
|
||
"\n",
|
||
"# # Наблюдаемые частоты\n",
|
||
"# observed = [19, 11, 7, 4, 9] # Пример из таблицы\n",
|
||
"\n",
|
||
"# # Ожидаемые частоты для λ0\n",
|
||
"# expected = []\n",
|
||
"# n = len(data)\n",
|
||
"# for interval in intervals:\n",
|
||
"# if interval[0] == -np.inf:\n",
|
||
"# prob = poisson.cdf(0, lambda0)\n",
|
||
"# elif interval[1] == np.inf:\n",
|
||
"# prob = 1 - poisson.cdf(interval[0] - 1, lambda0)\n",
|
||
"# else:\n",
|
||
"# prob = poisson.pmf(interval[0], lambda0)\n",
|
||
"# expected.append(n * prob)\n",
|
||
"\n",
|
||
"# # Статистика хи-квадрат\n",
|
||
"# chi2_stat = sum((o - e)**2 / e for o, e in zip(observed, expected))\n",
|
||
"# print(f\"Наблюдаемое χ²: {chi2_stat:.4f}\")\n",
|
||
"\n",
|
||
"# # Степени свободы\n",
|
||
"# k = len(intervals)\n",
|
||
"# m = 1 # Оценен один параметр\n",
|
||
"# df = k - 1 - m\n",
|
||
"# print(f\"Степени свободы: {df}\")\n",
|
||
"\n",
|
||
"# # Критическое значение и p-значение\n",
|
||
"# chi2_crit = chi2.ppf(1 - alpha, df)\n",
|
||
"# p_value = 1 - chi2.cdf(chi2_stat, df)\n",
|
||
"# print(f\"Критическое значение: {chi2_crit:.4f}\")\n",
|
||
"# print(f\"p-значение: {p_value:.4f}\")\n",
|
||
"\n",
|
||
"# # Вывод решения\n",
|
||
"# if chi2_stat > chi2_crit:\n",
|
||
"# print(\"Отвергаем H₀\")\n",
|
||
"# else:\n",
|
||
"# print(\"Не отвергаем H₀\")\n",
|
||
"\n",
|
||
"# # Обновим expected_grouped, чтобы быть уверенными, что это numpy-массив\n",
|
||
"# expected_grouped = np.array(expected_grouped)\n",
|
||
"\n",
|
||
"# # Вычислим границы для групп\n",
|
||
"# lower_bounds = [float('-inf')] + [min(g) for g in groups[1:]]\n",
|
||
"# upper_bounds = [max(g) for g in groups[:-1]] + [float('inf')]\n",
|
||
"\n",
|
||
"# # Вычислим разности и хи-квадрат члены\n",
|
||
"# diff = np.array(observed_grouped) - expected_grouped\n",
|
||
"# chi2_terms = diff**2 / expected_grouped\n",
|
||
"\n",
|
||
"# # Построим финальную таблицу\n",
|
||
"# df_final = pd.DataFrame({\n",
|
||
"# \"Group\": [\", \".join(map(str, g)) for g in groups],\n",
|
||
"# \"Lower\": lower_bounds,\n",
|
||
"# \"Upper\": upper_bounds,\n",
|
||
"# \"O_i\": observed_grouped,\n",
|
||
"# \"P_i\": np.round(expected_grouped / n, 4),\n",
|
||
"# \"E_i\": np.round(expected_grouped, 3),\n",
|
||
"# \"O_i - E_i\": np.round(diff, 3),\n",
|
||
"# \"(O_i - E_i)^2 / E_i\": np.round(chi2_terms, 4)\n",
|
||
"# })\n",
|
||
"\n",
|
||
"# print(\"\\nПодробная таблица для χ² при сложной гипотезе:\")\n",
|
||
"# print(df_final)\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 74,
|
||
"id": "d5e3c152",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# import numpy as np\n",
|
||
"# from scipy.stats import poisson, chi2\n",
|
||
"\n",
|
||
"\n",
|
||
"# # Оценка параметра λ\n",
|
||
"# lambda_hat = np.mean(data)\n",
|
||
"# print(f\"Оценка λ: {lambda_hat:.4f}\")\n",
|
||
"\n",
|
||
"# # Разбиение на интервалы (пример)\n",
|
||
"# intervals = [\n",
|
||
"# (-np.inf, 0),\n",
|
||
"# (1, 1),\n",
|
||
"# (2, 2),\n",
|
||
"# (3, 3),\n",
|
||
"# (4, np.inf)\n",
|
||
"# ]\n",
|
||
"\n",
|
||
"# # Наблюдаемые частоты\n",
|
||
"# observed = [19, 11, 7, 4, 9] # Пример из таблицы\n",
|
||
"\n",
|
||
"# # Ожидаемые частоты для λ0\n",
|
||
"# expected = []\n",
|
||
"# n = len(data)\n",
|
||
"# for interval in intervals:\n",
|
||
"# if interval[0] == -np.inf:\n",
|
||
"# prob = poisson.cdf(0, lambda0)\n",
|
||
"# elif interval[1] == np.inf:\n",
|
||
"# prob = 1 - poisson.cdf(interval[0] - 1, lambda0)\n",
|
||
"# else:\n",
|
||
"# prob = poisson.pmf(interval[0], lambda0)\n",
|
||
"# expected.append(n * prob)\n",
|
||
"\n",
|
||
"# # Статистика хи-квадрат\n",
|
||
"# chi2_stat = sum((o - e)**2 / e for o, e in zip(observed, expected))\n",
|
||
"# print(f\"Наблюдаемое χ²: {chi2_stat:.4f}\")\n",
|
||
"\n",
|
||
"# # Степени свободы\n",
|
||
"# k = len(intervals)\n",
|
||
"# m = 1 # Оценен один параметр\n",
|
||
"# df = k - 1 - m\n",
|
||
"# print(f\"Степени свободы: {df}\")\n",
|
||
"\n",
|
||
"# # Критическое значение и p-значение\n",
|
||
"# chi2_crit = chi2.ppf(1 - alpha, df)\n",
|
||
"# p_value = 1 - chi2.cdf(chi2_stat, df)\n",
|
||
"# print(f\"Критическое значение: {chi2_crit:.4f}\")\n",
|
||
"# print(f\"p-значение: {p_value:.4f}\")\n",
|
||
"\n",
|
||
"# # Вывод решения\n",
|
||
"# if chi2_stat > chi2_crit:\n",
|
||
"# print(\"Отвергаем H₀\")\n",
|
||
"# else:\n",
|
||
"# print(\"Не отвергаем H₀\")\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 8,
|
||
"id": "a937cbce",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Оценка λ: 1.9600\n",
|
||
"\n",
|
||
"Таблица до объединения категорий:\n",
|
||
" Значение k Наблюдаемая частота (Oᵢ) Теоретическая вероятность P(X=k) \\\n",
|
||
"0 0 19 0.1409 \n",
|
||
"1 1 11 0.2761 \n",
|
||
"2 2 7 0.2706 \n",
|
||
"3 3 4 0.1768 \n",
|
||
"4 4 3 0.0866 \n",
|
||
"5 6 2 0.0111 \n",
|
||
"6 7 1 0.0031 \n",
|
||
"7 8 2 0.0008 \n",
|
||
"8 14 1 0.0000 \n",
|
||
"\n",
|
||
" Теоретическая частота (Eᵢ) \n",
|
||
"0 7.043 \n",
|
||
"1 13.804 \n",
|
||
"2 13.528 \n",
|
||
"3 8.838 \n",
|
||
"4 4.331 \n",
|
||
"5 0.555 \n",
|
||
"6 0.155 \n",
|
||
"7 0.038 \n",
|
||
"8 0.000 \n",
|
||
"\n",
|
||
"Таблица после объединения категорий:\n",
|
||
" Группа Наблюдаемая частота (Oᵢ) Теоретическая вероятность \\\n",
|
||
"0 0 19 0.1409 \n",
|
||
"1 1 11 0.2761 \n",
|
||
"2 2 7 0.2706 \n",
|
||
"3 3 4 0.1768 \n",
|
||
"4 4,6,7,8,14 9 0.1016 \n",
|
||
"\n",
|
||
" Теоретическая частота (Eᵢ) \n",
|
||
"0 7.043 \n",
|
||
"1 13.804 \n",
|
||
"2 13.528 \n",
|
||
"3 8.838 \n",
|
||
"4 6.787 \n",
|
||
"[ 7.04292105 13.80412525 13.52804275 8.83832126 6.7865897 ]\n",
|
||
"\n",
|
||
"Подробная таблица для χ²:\n",
|
||
" Группа Oᵢ Pᵢ Eᵢ Oᵢ - Eᵢ (Oᵢ - Eᵢ)² / Eᵢ\n",
|
||
"0 0 19 0.1409 7.043 11.957 20.3001\n",
|
||
"1 1 11 0.2761 13.804 -2.804 0.5696\n",
|
||
"2 2 7 0.2706 13.528 -6.528 3.1501\n",
|
||
"3 3 4 0.1768 8.838 -4.838 2.6486\n",
|
||
"4 4,6,7,8,14 9 0.1357 6.787 2.213 0.7219\n",
|
||
"\n",
|
||
"Хи-квадрат статистика: 27.3903\n",
|
||
"Критическое значение (α=0.02): 9.8374\n",
|
||
"p-value: 0.000005\n",
|
||
"Вывод: Отвергаем нулевую гипотезу\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"import numpy as np\n",
|
||
"import pandas as pd\n",
|
||
"from scipy.stats import poisson, chi2\n",
|
||
"\n",
|
||
"# Данные: частоты по значениям\n",
|
||
"freq_table = {\n",
|
||
" 0: 19,\n",
|
||
" 1: 11,\n",
|
||
" 2: 7,\n",
|
||
" 3: 4,\n",
|
||
" 4: 3,\n",
|
||
" 6: 2,\n",
|
||
" 7: 1,\n",
|
||
" 8: 2,\n",
|
||
" 14: 1\n",
|
||
"}\n",
|
||
"\n",
|
||
"# Общее количество наблюдений\n",
|
||
"n = sum(freq_table.values())\n",
|
||
"\n",
|
||
"# Оценка параметра λ\n",
|
||
"lambda_hat = sum(k * v for k, v in freq_table.items()) / n\n",
|
||
"print(f\"Оценка λ: {lambda_hat:.4f}\")\n",
|
||
"\n",
|
||
"# Теоретические вероятности и частоты до объединения\n",
|
||
"unique_values = list(freq_table.keys())\n",
|
||
"probs = [poisson.pmf(k, lambda_hat) for k in unique_values]\n",
|
||
"expected = np.array(probs) * n\n",
|
||
"\n",
|
||
"# Таблица до объединения\n",
|
||
"df_individual = pd.DataFrame({\n",
|
||
" 'Значение k': unique_values,\n",
|
||
" 'Наблюдаемая частота (Oᵢ)': [freq_table[k] for k in unique_values],\n",
|
||
" 'Теоретическая вероятность P(X=k)': np.round(probs, 4),\n",
|
||
" 'Теоретическая частота (Eᵢ)': np.round(expected, 3)\n",
|
||
"})\n",
|
||
"\n",
|
||
"print(\"\\nТаблица до объединения категорий:\")\n",
|
||
"print(df_individual)\n",
|
||
"\n",
|
||
"# Группировка категорий:\n",
|
||
"groups = {\n",
|
||
" '0': [0],\n",
|
||
" '1': [1],\n",
|
||
" '2': [2],\n",
|
||
" '3': [3],\n",
|
||
" '4,6,7,8,14': [4, 6, 7, 8, 14]\n",
|
||
"}\n",
|
||
"\n",
|
||
"# Наблюдаемые и ожидаемые частоты по группам\n",
|
||
"observed_grouped = []\n",
|
||
"expected_grouped = []\n",
|
||
"expected_grouped = [\n",
|
||
" poisson.pmf(0, lambda_hat)*n,\n",
|
||
" poisson.pmf(1, lambda_hat)*n,\n",
|
||
" poisson.pmf(2, lambda_hat)*n,\n",
|
||
" poisson.pmf(3, lambda_hat)*n,\n",
|
||
" (1 - poisson.cdf(3, lambda_hat))*n\n",
|
||
"]\n",
|
||
"\n",
|
||
"for group, values in groups.items():\n",
|
||
" O_i = sum(freq_table.get(k, 0) for k in values)\n",
|
||
" p_i = sum(poisson.pmf(k, lambda_hat) for k in values)\n",
|
||
" \n",
|
||
" observed_grouped.append(O_i)\n",
|
||
"\n",
|
||
"# Таблица после объединения\n",
|
||
"df_grouped = pd.DataFrame({\n",
|
||
" 'Группа': list(groups.keys()),\n",
|
||
" 'Наблюдаемая частота (Oᵢ)': observed_grouped,\n",
|
||
" 'Теоретическая вероятность': np.round(probs_grouped, 4),\n",
|
||
" 'Теоретическая частота (Eᵢ)': np.round(expected_grouped, 3)\n",
|
||
"})\n",
|
||
"\n",
|
||
"print(\"\\nТаблица после объединения категорий:\")\n",
|
||
"print(df_grouped)\n",
|
||
"\n",
|
||
"# Расчёт χ²\n",
|
||
"observed_grouped = np.array(observed_grouped)\n",
|
||
"expected_grouped = np.array(expected_grouped)\n",
|
||
"diff = observed_grouped - expected_grouped\n",
|
||
"chi2_terms = diff**2 / expected_grouped\n",
|
||
"chi2_stat = np.sum(chi2_terms)\n",
|
||
"print(expected_grouped)\n",
|
||
"\n",
|
||
"# Степени свободы\n",
|
||
"k = len(groups)\n",
|
||
"m = 1 # число оцененных параметров\n",
|
||
"df_chi2 = k - 1 - m\n",
|
||
"\n",
|
||
"# Критическое значение и p-value\n",
|
||
"alpha = 0.02\n",
|
||
"chi2_crit = chi2.ppf(1 - alpha, df_chi2)\n",
|
||
"p_value = 1 - chi2.cdf(chi2_stat, df_chi2)\n",
|
||
"\n",
|
||
"# Подробная таблица расчётов\n",
|
||
"df_final = pd.DataFrame({\n",
|
||
" 'Группа': list(groups.keys()),\n",
|
||
" 'Oᵢ': observed_grouped,\n",
|
||
" 'Pᵢ': np.round(expected_grouped / n, 4),\n",
|
||
" 'Eᵢ': np.round(expected_grouped, 3),\n",
|
||
" 'Oᵢ - Eᵢ': np.round(diff, 3),\n",
|
||
" '(Oᵢ - Eᵢ)² / Eᵢ': np.round(chi2_terms, 4)\n",
|
||
"})\n",
|
||
"\n",
|
||
"print(\"\\nПодробная таблица для χ²:\")\n",
|
||
"print(df_final)\n",
|
||
"\n",
|
||
"# Вывод результатов\n",
|
||
"print(f\"\\nХи-квадрат статистика: {chi2_stat:.4f}\")\n",
|
||
"print(f\"Критическое значение (α={alpha}): {chi2_crit:.4f}\")\n",
|
||
"print(f\"p-value: {p_value:.6f}\")\n",
|
||
"\n",
|
||
"if chi2_stat > chi2_crit:\n",
|
||
" print(\"Вывод: Отвергаем нулевую гипотезу\")\n",
|
||
"else:\n",
|
||
" print(\"Вывод: Нет оснований отвергнуть нулевую гипотезу H0\")\n"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "07c231d4",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Пункт g) Наиболее мощный критерий проверки гипотезы $H_0 : λ = λ_0 = 2.0$ против $H_1 : λ = λ_1 = 4.0$\n",
|
||
"### Пункт g) Наиболее мощный критерий проверки гипотезы\n",
|
||
"\n",
|
||
"**Логарифм отношения правдоподобия**\n",
|
||
"\n",
|
||
"Функция правдоподобия для распределения Пуассона:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"L(\\lambda) = \\prod_{i=1}^n \\frac{\\lambda^{X_i} e^{-\\lambda}}{X_i!}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Логарифм отношения правдоподобия:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\ln \\left( \\frac{L(\\lambda_1)}{L(\\lambda_0)} \\right) = \\sum_{i=1}^n \\left( X_i \\ln \\left( \\frac{\\lambda_1}{\\lambda_0} \\right) - (\\lambda_1 - \\lambda_0) \\right).\n",
|
||
"$$\n",
|
||
"\n",
|
||
"**Критерий отношения правдоподобия**\n",
|
||
"\n",
|
||
"Для проверки $H_0$ против $H_1$ используется сумма наблюдений $T = \\sum_{i=1}^n X_i$. Критерий принимает $H_1$, если:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"T > k,\n",
|
||
"$$\n",
|
||
"\n",
|
||
"где $k$ определяется как:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"k = \\text{qpois}(1 - \\alpha, n\\lambda_0).\n",
|
||
"$$\n",
|
||
"\n",
|
||
"**Смена гипотез**\n",
|
||
"\n",
|
||
"Если поменять местами гипотезы, новая нулевая гипотеза $H_0 : \\lambda = \\lambda_1$, а альтернатива $H_1 : \\lambda = \\lambda_0$. В этом случае критерий принимает $H_0$, если:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"T < k',\n",
|
||
"$$\n",
|
||
"\n",
|
||
"где $k'$ определяется как:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"k' = \\text{qpois}(\\alpha, n\\lambda_1).\n",
|
||
"$$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 47,
|
||
"id": "635fbf6b",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Сумма наблюдений: T_obs = 98\n",
|
||
"Порог для H0:λ=2.00: k = 121\n",
|
||
"Порог для H0:λ=4.00: k' = 172\n",
|
||
"\n",
|
||
"Проверка H0:λ=2.00 vs H1:λ=4.00:\n",
|
||
"Не отклоняем H0: T_obs = 98 ≤ 121\n",
|
||
"\n",
|
||
"Проверка H0:λ=4.00 vs H1:λ=2.00:\n",
|
||
"Отклоняем H0: T_obs = 98 < 172\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"import numpy as np\n",
|
||
"from scipy.stats import norm\n",
|
||
"\n",
|
||
"# sum_X = np.sum(data)\n",
|
||
"\n",
|
||
"# # Критическое значение для H0: λ=2.0\n",
|
||
"# mu_H0 = n * lambda0\n",
|
||
"# sigma_H0 = np.sqrt(mu_H0)\n",
|
||
"# z_crit = norm.ppf(1 - alpha)\n",
|
||
"# C = mu_H0 + z_crit * sigma_H0\n",
|
||
"\n",
|
||
"# print(f\"Сумма наблюдений: {sum_X}\")\n",
|
||
"# print(f\"Критическое значение C: {C:.1f}\")\n",
|
||
"\n",
|
||
"# if sum_X > C:\n",
|
||
"# print(\"Отвергаем H0: λ=2.0\")\n",
|
||
"# else:\n",
|
||
"# print(\"Нет оснований отвергнуть H0: λ=2.0\")\n",
|
||
"\n",
|
||
"# print(\"\\nПоменяли гипотезы местами\")\n",
|
||
"# # Смена гипотез местами\n",
|
||
"# mu_H0_swapped = n * lambda1\n",
|
||
"# sigma_H0_swapped = np.sqrt(mu_H0_swapped)\n",
|
||
"# z_crit_swapped = norm.ppf(alpha)\n",
|
||
"# C_swapped = mu_H0_swapped + z_crit_swapped * sigma_H0_swapped\n",
|
||
"\n",
|
||
"# print(f\"\\nКритическое значение C' (при H0: λ=4.0): {C_swapped:.1f}\")\n",
|
||
"\n",
|
||
"# if sum_X < C_swapped:\n",
|
||
"# print(\"Отвергаем H0: λ=4.0\")\n",
|
||
"# else:\n",
|
||
"# print(\"Нет оснований отвергнуть H0: λ=4.0\")\n",
|
||
"\n",
|
||
"# print(\"---\")\n",
|
||
"\n",
|
||
"# sum_data = np.sum(data)\n",
|
||
"# mu0 = 2*n\n",
|
||
"# mu1 = 4*n\n",
|
||
"# c = norm.ppf(1-0.02, mu0, np.sqrt(mu0))\n",
|
||
"# print(f\"\\ng) Критич. значение: {c:.1f}, сумма: {sum_data}\")\n",
|
||
"# print(\"Отвергаем H0\" if sum_data > c else \"Не отвергаем H0\")\n",
|
||
"\n",
|
||
"\n",
|
||
"# # Сумма наблюдений\n",
|
||
"# T_obs = np.sum(data)\n",
|
||
"\n",
|
||
"# # Критерий для H0: λ = λ0 vs H1: λ = λ1\n",
|
||
"# k = poisson.ppf(1 - alpha, n * lambda0) # Квантиль для порога k\n",
|
||
"# decision_H0 = \"Отклоняем H0 в пользу H1\" if T_obs > k else \"Не отклоняем H0\"\n",
|
||
"\n",
|
||
"# # Критерий для H0: λ = λ1 vs H1: λ = λ0\n",
|
||
"# k_prime = poisson.ppf(alpha, n * lambda1) # Квантиль для порога k'\n",
|
||
"# decision_H1 = \"Отклоняем H0 в пользу H1\" if T_obs < k_prime else \"Не отклоняем H0\"\n",
|
||
"\n",
|
||
"# # Вывод результатов\n",
|
||
"# print(f\"Сумма наблюдений T_obs: {T_obs}\")\n",
|
||
"# print(f\"Порог k для H0: λ = {lambda0}: {k}\")\n",
|
||
"# print(f\"Решение для H0: {decision_H0}\")\n",
|
||
"# print(f\"Порог k' для H0: λ = {lambda1}: {k_prime}\")\n",
|
||
"# print(f\"Решение для H1: {decision_H1}\")\n",
|
||
"# print(\"---\")\n",
|
||
"from scipy.stats import poisson\n",
|
||
"\n",
|
||
"# Данные наблюдений\n",
|
||
"data = list(map(int, \"0 1 2 0 0 7 1 0 2 1 0 1 2 2 0 0 1 8 0 0 14 4 3 0 0 3 0 6 2 2 1 0 0 2 0 4 0 0 3 3 1 1 0 0 6 8 1 4 1 1\".split()))\n",
|
||
"n = len(data)\n",
|
||
"T_obs = sum(data)\n",
|
||
"\n",
|
||
"# Параметры\n",
|
||
"alpha = 0.02\n",
|
||
"lambda0 = 2.00\n",
|
||
"lambda1 = 4.00\n",
|
||
"\n",
|
||
"# Вычисление порогов\n",
|
||
"k = poisson.ppf(1 - alpha, n * lambda0)\n",
|
||
"k_prime = poisson.ppf(alpha, n * lambda1)\n",
|
||
"\n",
|
||
"# Результаты\n",
|
||
"# print(f\"Количество наблюдений: n = {n}\")\n",
|
||
"print(f\"Сумма наблюдений: T_obs = {T_obs}\")\n",
|
||
"print(f\"Порог для H0:λ=2.00: k = {int(k)}\")\n",
|
||
"print(f\"Порог для H0:λ=4.00: k' = {int(k_prime)}\")\n",
|
||
"\n",
|
||
"# Проверка исходных гипотез\n",
|
||
"print(\"\\nПроверка H0:λ=2.00 vs H1:λ=4.00:\")\n",
|
||
"if T_obs > k:\n",
|
||
" print(f\"Отклоняем H0: T_obs = {T_obs} > {int(k)}\")\n",
|
||
"else:\n",
|
||
" print(f\"Не отклоняем H0: T_obs = {T_obs} ≤ {int(k)}\")\n",
|
||
"\n",
|
||
"# Проверка инвертированных гипотез\n",
|
||
"print(\"\\nПроверка H0:λ=4.00 vs H1:λ=2.00:\")\n",
|
||
"if T_obs < k_prime:\n",
|
||
" print(f\"Отклоняем H0: T_obs = {T_obs} < {int(k_prime)}\")\n",
|
||
"else:\n",
|
||
" print(f\"Не отклоняем H0: T_obs = {T_obs} ≥ {int(k_prime)}\")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "863bdc6c",
|
||
"metadata": {},
|
||
"source": []
|
||
}
|
||
],
|
||
"metadata": {
|
||
"kernelspec": {
|
||
"display_name": "venv",
|
||
"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.13.2"
|
||
}
|
||
},
|
||
"nbformat": 4,
|
||
"nbformat_minor": 5
|
||
}
|