Skip to content

Instantly share code, notes, and snippets.

@imadmali
Created June 30, 2021 00:01
Show Gist options
  • Select an option

  • Save imadmali/9c78b8fd230e8e91b71a9935298a13aa to your computer and use it in GitHub Desktop.

Select an option

Save imadmali/9c78b8fd230e8e91b71a9935298a13aa to your computer and use it in GitHub Desktop.
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# P-values and Big Data\n",
"Author: Imad Ali \n",
"Date: 2021-06-23\n",
"\n",
"AB testing, experimentation, randomized control trials, or whatever else you wanna call it, are useful for measuring the effect of an intervention (e.g. marketing campaigns, medical treatments). This is done by collecting data on a control group that didn't receive the intervention and the variant groups that did. Typically a test statistic and associated p-value are calculated using the collected data to determine whether there is a \"significant\" difference between the groups. Some p-value calculations have behavior under large sample sizes that isn't commonly known. Specifically, p-values can asymptotically tend towards a statistically significant result. For a fixed set of parameters this means that a statistically significant result can be achieved simply by having a larger sample size. This is an important result to be aware of given the large scale experiments that currently take place in industry. This doesn't mean that p-values are useless, but rather that they should be interpreted within the context of the effect size (the difference in what you're measuring between the two groups).\n",
"\n",
"The relationship between p-values and large samples is shown in two scenarios: Welch's t-test and the chi-squared test. Welch's t-test is useful in situations where your measurement is normally distributed (e.g. engagement on a platform in hours). The chi-squared test is useful when your measurement refers to counts (e.g. number of users that engaged with an advertisement). "
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import scipy\n",
"from scipy import stats\n",
"import matplotlib.pyplot as plt"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Welch's t-test"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The graphics below illustrate the relationship between p-values computed from Welch's t-test and sample size for different standard deviation values and a trivially small effect size ($2.5-2=0.5$). The _stochastic_ graphic shows the relationship via simulation and the _deterministic_ graphic shows the result using the analytical solution.\n",
"\n",
"As a recap, the formula for Welch's t-test is provdided below. This assumes normally distributed samples $x_a$ and $x_b$ that do not have equal variance. The null hypothesis asserts that the means of each group are identical. The test statistic is calculated as follows,\n",
"\n",
"$$\n",
"\\text{stat} = \\frac{\\bar{x}_a-\\bar{x}_b}{\\sqrt{s_a^2 + s_b^2}}\n",
"$$\n",
"\n",
"where $\\bar{x}_i$ is the mean and $s_i$ is the standard error of sample $i$. The standard error is calculated as $\\sigma_i^2 / N_i$, where $\\sigma_i^2$ is the variance and $n_i$ is the sample size of sample $i$.\n",
"\n",
"The test statistic follows the t-distribution with the following degrees of freedom calculation,\n",
"\n",
"$$\n",
"\\eta_i = \\frac{s_i^2}{n_i}\\\\\n",
"\\text{dof} = \\frac{(\\eta_a + \\eta_b )^2}\n",
" {\\frac{\\eta_a^2}{n_a-1} + \\frac{\\eta_b^2}{n_b-1}}\n",
"$$\n",
"\n",
"Using this information the (two-sided) p-value is calculated using the t-distribution CDF $F$.\n",
"\n",
"$$\n",
"\\text{pvalue} = 2 \\cdot F(-|\\text{stat}|, \\text{dof})\n",
"$$\n",
"\n",
"`scipy.stats.ttest_ind` can be used to implement this test. A manual Python implementation of the calculation is provided in the `welch_ttest` class."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"mean_a = 2\n",
"mean_b = 2.5\n",
"N = np.arange(100,20_001,100)\n",
"S = np.arange(1,31,1)\n",
"sim = np.zeros((len(S),len(N)))\n",
"for n in range(0,len(N)):\n",
" for s in range(0,len(S)):\n",
" rep=[]\n",
" for i in range(0,100):\n",
" a = np.random.normal(loc=mean_a, scale=S[s], size=N[n])\n",
" b = np.random.normal(loc=mean_b, scale=S[s], size=N[n])\n",
" test = stats.ttest_ind(a, b, equal_var=False)\n",
" rep.append(test.pvalue)\n",
" sim[s,n] = np.mean(rep)"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAtsAAAFrCAYAAAAAUkVAAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAB9UUlEQVR4nO2dd5gkV3W3f6fT9OSZzVmrHBESkkBgMhiTg0nGfDYYm/AZMBiwjcE2svkcsY2xMWBhMBgTJLLIIkmAEUI559VKm3dy7Ol4vj+qpu+5d7t7ekPNzI5+7/Pss7e7bt17bqjq2zVVb4uqghBCCCGEEHLsSS11AIQQQgghhKxUuNgmhBBCCCEkIbjYJoQQQgghJCG42CaEEEIIISQhuNgmhBBCCCEkIbjYJoQQQgghJCG42CaErAhE5LUi8rOljuNIOR7iF5GdIvLMY1zmHSLy1GNZZlzu50Xkxce63Bb1XSIi/7MI9TxJRO5pI98LROSypOMhhCwMF9uErDDaWRCJyFUi8nvHoK6nisjuoy2HPHJR1bNV9apjWaaInAvg0QC+3kZeFZFTjmX9x5IwPlX9qaqevtB+qvoNAGfHfUEIWUK42CaEELLSeCOAzyp/te3zAN6w1EEQ8kiHi21CVhAi8hkA2wB8Q0SmReSPG+T5awBPAvDhOM+H4/fPEJHvi8ioiNwjIq8w+zxXRO4UkSkR2SMi7xKRbgDfAbApLmdaRDY1qO9TIvKxuOwpEblaRE5oEv9HReQfg/e+LiLviNPvFpEH4nLuFJGXNClne3xFMGPe867mi8jrROQuERkTke/NxyQRHxSRgyIyKSK3icg5Ter5nbiMKRHZISJvNNueKiK7ReSdcVn7ROR3zPbVInJFXMcvAZzcqI44b15E/kdERkRkXESuE5H1hxHDH5sYXhyP573xWL/H5L9ERL4kIpfF5d0oIo9uElPKjMeIiFwuIqua5F0jIt+MYx8VkZ+KSCreVv9LTLx9fi7NxGO4Pd72fBG5Oc7z8wWu2D4HwNWm/lPieTchIsMS314hIj+Js9wS1/nK+P3Xi8j9caxX2HktImeLO04O2P4DkBOR/4777g4RudDs13TuHk58Evw1SUS2ishXRGQoHocPm3iuAvC8Fv1ECFkMVJX/+I//VtA/ADsBPHOBPFcB+D3zuhvALgC/AyAD4HwAwwDOirfvA/CkOD0I4DFx+qkAdi9Q16cATAF4MoAOAB8C8LMmeZ8cxyGmrgKATfHrlwPYhOhCwSsBzADYGG977Xy5ALYDUACZRm0G8CIA9wM4M27vnwH4ebzt1wDcAGAAgMR5NjaJ93mIFskC4CkAZoO+qQD4KwBZAM+Ntw/G278A4PK4788BsKdFv7wRwDcAdAFIA7gAQN9hxPAXcQyvBzAE4HMAegGcHffviXH+SwCUAbwszv8uAA8CyIZzC8DbAPwCwJZ4XP8DwOebxP+3AD4Wl5lF9GVPwjKDff4GwE/i/OcDOAjgcXH7XxPv19Fgv+547Nea9z4P4L2I5k0ewBPNNgVwinn9dERz/zFxu/4NwE/ibb2IjoV3xuX0Anic6bu5eJzTcZt/YcptNXcPJ76nIj7m4npuAfDBuN3hvqvi/fuW+rzEf/z3SP7HK9uEEAB4PoCdqvpfqlpR1ZsAfBnRAgGIFmBniUifqo6p6o2HWf63VPUnqlpEtKh4vIhsbZDvp4gWB0+KX78MwDWquhcAVPWLqrpXVWuqehmA+wA89jBjAYA3AfhbVb1LVSuIFnbnSXR1u4xoEXUGogXhXaq6r1EhqvotVX1AI64GcKWJHXFZf6WqZVX9NoBpAKeLSBrASwH8harOqOrtAD7dIt4ygNWIFl1VVb1BVScPI4a/VtUyogX+GgAfUtUpVb0DwJ2I7m+e5wZV/VKc/58RLeAubtKH71XV3fG4XgLgZWL+mhDEsBHACXFf/FRVm97iEV9h/k0AL43jeAOA/1DVa+P2fxpAsUlcA/H/U0H9JyD60janqq0eRH01gE+q6o1xu/4U0Xzdjug42a+q/xSXM6Wq15p9f6aq31bVKoDPwPTrAnP3cOKzPBbRAv6P4nkU7jvfBwOH7EkIWTS42CZkhSPRLRzzf5p/T5NsJwB4XPwn+nERGUe06NgQb38poit2D8V/7n78YYaxaz6hqtMARhHdfvIeE9vH4gXYFwC8Ks7+mwA+a9ry2+ZWgnFEV4TXHGYsQNTeD5lyRhFdGd6sqj8C8GEA/w7goIhcKiJ9jQoRkeeIyC/iWwrGEfWRjWckXszPMwugB8BaRFfUd5ltD7WI9zMAvgfgCyKyV0T+QUSyhxFDNU4X4v8PmO2FOKZ57FjVAOxGtKALOQHAV00f3gWgCmB9g7wfQPSXhCslutXl3c0aKiLnI+r/l6jqkKnrncH83NokrvH4/17z3h8jGt9fxrd3vK5Z/XGZ9bGI5+sIgM1xnQ+02He/Sc8CyM9/+Vhg7h5OfJatAB4K5phlvg/G2yyPEJIAXGwTsvLwrhiq6ptUtSf+9zeN8iBaYF2tqgPmX4+q/t+4jOtU9UUA1gH4GqLbHxqV04z6VWwR6UH05+29qvo3JrY3xVk+j+gK6QmIbhv4crzfCQA+DuAtAFar6gCA2xEtUkJm4v+7zHsbTHoXgDcG7e1U1Z/H7f1XVb0AwFkATgPwR2EFItIRx/aPANbH8Xy7STwhQ4hu77BX97c1yxxfDf5LVT0LwBMQXWH97aOMoRl2rFKIbhPZ2yDfLgDPCfowr6p7GsQ/parvVNWTALwQwDtE5BlhPhGZn19vjv+6Yuv666CuLlX9fIO6ZhAtiE8z7+1X1der6iZEt+R8RJobSPYiWtzPx9SN6K8Ke+I4TmqyX1MWmruHGZ9lF4BtTf6aAES3QO2c/ysIIWRp4GKbkJXHASy8IAjzfBPAaSLyWyKSjf9dJCJnikhORF4tIv3xn/QnAdRMOatFpH+B+p4rIk8UkRyA9yO6l3VXo4zxImsYwH8C+J6qjseb5u/FHQKiBwMRXR1sVMYQosXR/xGRdHyl0D6A+DEAfyoiZ8dl9YvIy+P0RSLyuPjK8Qyi+3BrOJQcont6hwBUROQ5AJ61QD/Mx1cF8BUAl4hIl4icheg+5IaIyNNE5FHx7SeTiG47qB1NDC24QER+PV7AvR3R7Rq/aJDvYwD+WtyDpWtF5EVN4n9+/BCgAJhAdAW8FuTJAPgSgP9R1cuDIj4O4E3xuIiIdIvI80SkF435NqL71+fLfrmIbIlfjiGaR3YO22Ph8wB+R0TOi7/M/A2Aa1V1J6LjZKOIvF1EOkSkV0Qe1yQGS8u5e5jxWX6J6B7yv4v7JC8iv2K2PwXRQ8yEkCWEi21CVh5/C+DP4j9Xv6tJng8huno8JiL/qqpTiBZpv4Hoyt5+AH+PaCEHAL8FYKeITCK6V/fVAKCqdyNanOyI62v0Z30geiDvfYhu17gAwP9ZoA2fA/DM+H/Edd0J4J8AXINoAfIoAP/boozXI7oiPYLoQcCfm7K+GrfvC3GbbkdksACAPkSLuzFEtxOMILoNwiPusz9AdJV/DNEtL1cs0C7LWxDdvrEf0UOk/9Ui7wZEC9FJRLdrXA3gM8cghkZ8HdEDfGOIxv3X4y9ZIR+K67pSRKYQLcibLTxPBfADRPesXwPgI6r64yDPFkT3mr/d3Fo0LSLbVPV6ROP54Tiu+xE9ENuMSwG8Ol7cA8BFAK4Vkek45rep6o542yUAPh3P31eo6g8A/DmivxjsQ/Ql7TeA+pj/KoAXIBq3+wA8rUUciPdbaO62HV9QbjWO5RQADyO65eeVJsurED24SghZQuafBieEkEQQkU8hsif82VLHQlojIpcgeghzoS9Dyx4R+RyAy1X1a0sdy1IgIi8A8Fuq+ooFMxNCEqXZfV6EEELIcYuq/uZSx7CUaPQLkt9Y6jgIIbyNhBBCCCGEkMTgbSSEEEIIIYQkBK9sE0IIIYQQkhBcbBNCCCGEEJIQiT0gKSJ5AD9BpA7LAPiSqr5PRE5E9AtxqwHcgOhp6VKrsnoHs7pmc2Qgy4r/Q1kdUq2np2sd9bQGv+nQkXLmqpK6ZmdQ9fLZ/cT8XsesKTssLx0oeOc0V0/3pQr1dEGzXj4bYYe48gq1nJfPxlTSNJqRFhfHXM3VlU/5fZYxfVasub5Ii39L0UQ5X0+nzLaM+O3tMOVLUEZV3fe5zpQb5uFSj5fP1p1Pu74I6ypUXbv6Mq5vh0q+brfTlJE17R2Z6/by9XXMNcw3Wcl7+XozLt+k6ZeujG9Es+2tBXdo1cx323LVjWMm5c/BbMqMYzVj3vfzFStuW2e23PB9AKiYulKm7HTK79t82o2jnRcS/G6NnQuFOTdXO/P+YWzn40zJP34sWjPHXMqvq9eMz+RspwnCz5fNuL4Rc2CVypkgn2tjueS2dXT441g0+3XmzLFZ8NshGdOHpl+06l/HSJt8YX9WiiZG065ch3/cVkyZqq6RWvXPdems6ws7VnYeAEDO9oXZViv7sUvGtMuGrsHv5thj35QRnBKQ7nL1VsrB+cxOSVt8UIZUTPttvvDyka081eR9AFIyfZs224Impsw42n6Ssp9RO0wZJtbgdAa10/OQNppN5qND/NNA0z7TbHgCMnGYMjRcBdj2m3aF42j73bZLM0HfVmxQjWMFAE0132Y+OlDLu0JSxaDfbRmmrhYfm37sHf4ApWfM+dz/+PbLsP3ZYj567TAf8+GhZNsfLHm8cbDbDonPdruJIx2stuw2W161q71xDOdFOMfru4SnC9tndnwO4+e4ms/BIGPj0/Qh2PYXhncPq+ra9qNpTZI2kiKAp6vqdPzjED8Tke8AeAeAD6rqF0TkYwB+F8BHWxW0ZnMH3veVRwEANmfHvG0nZSbq6Z8WttfT5aC3T+1wv6K7q7y6nl6dnvby2cWsXSzcNLvdy3dKh/u14950wdt2z5xTDf9azx319G1FX0FsF3en5g7W0zfPbfHy2bY8XHKxp4Kzc39mtp6+a2ZjPX161wEv3/qs67N759yP6vUH7fjOgbPr6U6zqFzdMePlO7lrqJ62X0IAYKLifsDv3K6H6+lP7HqSl6835xZVp/W4vliT9cfnjmnXh7866Pr240F5Zw/sq6c35NyPp33qTl8D/OxT7mqY74cHT/fyPXnt/fX09/edUU+ft9r/sbzxslsQlmr+GX624s6ueybcb8Cs6fH7c33nVD19/7j71e0N3VNevh1jq+rpc9bub7gPAAyPuC8ivX1ujHvyRS/faQNuHB+eHqynw0V+V8adrW+6a3s9/agz/N+oGci5+fiLnSeiGeWC+5To6PFjevqJ99XT37n5UfV0ptufZ+tWubHLm0Xkg3v9vti4brye3vuwO5ZOO3mfl2/HAbffmZtd395++wlevswqN28zZpE7N+N/Ye4fcH2RSfufRiM7XF9rlytj+wlDfr4ZdyzNzbk+K0/4XwAGN7q+6OlwY3Vg3P9CeuLakXp636T7NfrJ/X6+7IBrY7Xi5rRWwkW5+UKxz30hDT/cB8937Tqwf8DfWLKf/Ob8VvM/gfMH3Dmx0uny1XL+ObGWNzGZRX74pS69y8VbHjTzPeuPVdegO35mD7gv7p37/M+bwkluHqdH3FilC347yoOm/OBLU8eo64viOhdTZtLv98yM269qFqKljf4xIrNu7LLj5ov/qmB11GsuJO13cytc9NXMF4r0rLkgtNo/X+QPmi/udjHX4Y9BtcU49ux0sU+d7uLr3uGvMM3HDexHUak/WDial3ZMyif7n4G917rzeWGD+dIZLKKzk9JwW7hg7d7t0tPm92LDLwM180XJzgPAb1fHqLkguNGfP7aMqjlF9Ozy81XNdaXcuNtn9Dx/HDuG3DjaGIKPfJhrYN4CuNLp57Oxl/pdTK2+NIRfSM01IWTdKRZza4Lj24yxjSn84pofdvvd9J/vfAjHkMRuI9GI+ZVSNv6nAJ6O6McZAODTAF6cVAyEEEIIIYQsJYnesx3/TPLNAA4C+D6ABwCMq+r899vdADY32fcNInK9iFw/Ndbox8sIIYQQQghZ3iS62FbVqqqeh+hneB8L4IzWe3j7XqqqF6rqhb2DLW6WIoQQQgghZJmyKDYSVR0H8GMAjwcwICLzd9psAbCn2X6EEEIIIYQczyT2ozYishZAWVXHRaQTwJUA/h7AawB82TwgeauqfqRVWeeem9Urvh09sHRV8KCiNTyMVpzhwj6YCAAz5iHD7eYBr9Up/679G0rujvmfzbiH5NZkJr18F+Xdw3654PHb/zUPau42DzQ+NLfKy/fO9T+opy+fuKCefn7vLV6+78+cVU9fM3pSPV0LntJ49IB7+qLDPInyYMF/SKw74x7eqZnHeXvS/sNp68wDg60euEyZ9u8K2viYHveMwY/GzqynR4qBFSTnnqqomAcLTzUPSwLA3rmBenpfwT3UVak1/9747PVuLuwpDnrbbhzdWk9v7HLt7cvOefmuO+CeZlnVOYtm2Acfn7z1AW/bDx44rZ7O5dw82zIw7uWzpom1nf4DopbZivuLT8aYRW7Z49+Z1d3pxtVaJ05cNerl29Ll4hgvueNix8RqL9+BXa4PO1e7cXv8lp1evl/scQ8TbhmY8LbZB2JvesiNQWeX/zTdeRvcd/F9s268R83DggAwPeOe8jl1o5sz9+xe7+WzVpDyrOu/zZv9vujJuT7bOezmdC2wjNgyela5eaHBo/ezky6+LRv9uoYm3Xmr+oBL958z4uWbLbq65h52DzHqYPAgnHkAs888EDs+4huA7MNGqSl3fuze7o/V1D7zwKR9sDB4oC834uaWffDPPvQJAKlxV1etK9RzmHxFYwgJHmi0GgLtM+ahvP8UX3nKPagqxh7S/aD/RFpxtbG29LiYOob9fBnzLPPMdtOuoBkdpi+K5gFE7QzMQ8Pm4cGgiaVN7lhID7t2BKdplPsbx5sKHkwtDZoH5jpNwKExY871bWa6uVXGPuBo0/lhv0D7cKJ9gLNjJHjA1oQUPiBpx7u0yjwsOu2XkTZ2Etv+THDKnjrZldG515ih/NM+Zja5OHofav4Q35yZP50HGj8sCQAdEy7fjHmgMZCdoWuvybfZP85sW+wDfmFMWTNX03PGILU2MLiYKW4feu3zP768OOw+qWA+2gcuOw82fh8AstPmmMubh1T957O9dnXv8eeFbUt+xG2bWx08UD1q56rbli765XWYB0R/+bl33aCqF+IYkaSNZCOAT4tIGtHhfLmqflNE7gTwBRH5fwBuAvCJBGMghBBCCCFkyUhssa2qtwI4v8H7OxDdv00IIYQQQsiKhr8gSQghhBBCSEJwsU0IIYQQQkhCcLFNCCGEEEJIQnCxTQghhBBCSEIkpv47lmw/p1ff++XzAADpQLNXMkq/iarTlZ2UG/LyZY0K74bZE+vpp/Xc5eW7t7Shnp42rpqetO8E2jnndHqnde73tlkFYUfKaalCVd99hXX19PMHbq6nvzl+npfv1wdvqKe/M3luPX1+10Nevp9NOrXcpo7xero/8B7dOr2lnn7O4G319PUzJ3r5ysbvc/Oo22d9l69BnKs6N8/Dk75ar1By2/70zO/W05cfaG7U2TnuVGtb+n0N2bn9TgU3XXW+pBuHt3r5UsZT9aYTrq6nP3Dvr3r5JiadgvAJJznX0e1DG718q7tdH6aMn6wGXzF0Zr+bCx0pX0N2zUHXv+u6plwMJV8/aZktu/7b2O33+y7T1zVzGOezfr0n9w/X0zfuc/2Uy/j5LtrgdJY/ut9pLzevGffyTc65frfKxWza15qNDjuH06aNY962EdPvxdHm7c/0lRq+Lyn/PGBVimmzrTDn+7DK006h1jXotHhzszkvnx1WHXHtldW+5yrf6eIr7HfHfXZNwctXKbljqVbydXKpSXcOS61355mebv+cMzHhdIdiVHg6EcTeY1SA06b9vb4i0GoQq6aNCD4SspNujEsbXRndg34bCw8ZHWHGFSKVQBE4Ycob8Mcxv9UdF96YjPhutJrRCVqlX2bCP8d27XN1zxqNW3bSj6lmqir3GVWfP1Serq5klHvZQEFXMYq7qtEbSs2vt2u37Qu/460mr7bWzbP8/b5DzVomi2tcv3Tv9oOveYo3l7btAPzxth8dYV/Yj7PCZncuGbzFzzi31qrW3PuhMs6OT/Ax742Xp74LdYSmayrGDir+qQ724zxjpvGcbzk9pM3zdPr2W79eczoL29E5bPR0q1ybxD91onuf27GW8efMhPuYx8DdLp2qBPPHKO5mN7h054FAn7fOlG+SHaN+vlrWaBXLblvwMYep7SYmo2Ls2R2W59LZWbdtamughDS7dYz5ZWTMfjObjKYy6E9L55A5rgJFSMrsd/2n33lM1X+8sk0IIYQQQkhCcLFNCCGEEEJIQnCxTQghhBBCSEJwsU0IIYQQQkhCcLFNCCGEEEJIQnCxTQghhBBCSEIcF+q/c87N6Ze+Fan2bi5u8rZNGd3fQNo5gX4xfbKX73n9t9TTD5edWm6q5mvH1macXq3buIl2ltZ6+TZnncrsnjlfE3ew5BRYW/Oj5v0+L9+z+m5vGNPdBb+NWeOjGSk5vdj5Pb7678Gii/GUvHMTXTlytpdvY97p9HYVrD4u0FJlnG7q7J599fTlO8/38r3xlJ/V098d8us6ucdp524d2+xiHVrl5fuV7Q/W02ty0y7fjO9i2pB3ajCrVbx2aLuX78xB1/5VuZl6ulD1VXAPz7g4hgtOR3eS0eUBwI4Jp3q0R8yJfaNevtsOurkwO+vryk5a78osVp1z6OBkj5dvXZ9rfz7j2rhz2O+zPqOGG51wsdfK/nfop51+bz09UXaOqlt3b/byWcdSR4fzORWCdtSMym3LBncchOq/B+9xfdG1cdrbZrWD4w8NuBAGfNVfT69r4+SoayPK/lyVTld3ZreLt/tsXzk4vscdg13r3byYK/j6vJxpv5h+KYwFmkIzGVKzzhOWG/XHYG6TKy877vvEMqe4OZ0ySr+Z4S4vH7JOB5Yy2r5DtHj9bs7IlJtn2u2PT8psy0y5/iz3BXqt9e68mrnDzdVQm2U1fpkZozULPmJSJbet0hnowDa78c487OZqqGCr5tx+aVNeLVS1bXZeN93vyus8GOjFjL6s1G+0ZsE8q3Y21oalt814+cpDbp5kx4zeb1OgXxxzhXQEc8b2TaXHtHfOj6k84Ma1795gUAz248fqDXt3+PXadpXtR1YwjmUTk42olvEzdu535VvFYs63unrx5ab8bVWzn1XrZf3Tiqfa69nn+qWwym/j1EmmDKMVtFo4wO+L4kBj/SAAdIy7/YpG6VfLBvmMuq7Y7/J17/cdgVb3V+r1x9seC4UN5pw94ufr2WPKNM1Kl/02jp7hCrQaxLAvJk5x6fyQ0QAGysGsORSk1ljNB/h9kZ2xx5XfjnTRbcvM+f1UWOVir3Tb/vTPdVNbXb7stNFyZv26rBLz9n9+B9V/hBBCCCGEHA9wsU0IIYQQQkhCcLFNCCGEEEJIQnCxTQghhBBCSEJwsU0IIYQQQkhCcLFNCCGEEEJIQjT3BC0jpmod+NHMaQCA7TlfyXav0e6d2rG/nl6Xm0QzHiyuq6cf1/2At21nySneZlPOA3Nxp5/vv4afVE+f3rXf2/bywevq6S+MPa6ergReqm+Mn1dPP6p7dz1tVX8AsDE7Xk9vy43U07tLvgpu56zT5F03ckI9fcGqh7189884ReDTV91dT+8rD3j5piuu/deMOlfSmi7fe/S1/a4dA7mCt+2hWRfjYN7td/GZD3r5HjSx/2jPqfX0het3efkOFp16bLbifFCTc77+zCruiqbfd00Pevl6c07v2GPSYb7tRvFn9X6VXv/76syMq3ftKn8OdhqN38Fp147CqK+TG8u48R809rdy0T9c071GBZdy6YtP3enlu2N0Qz1dqri+SKUD3VTNaZBmdzoP17Zz9nn55ioujt373PhmjC4PANSo6iTwv40NOz1mepXT/VWnfFdWKe/KTI27emuBxk6rLvbqNqMLnPT7NlVy4zVr1Hodg3NevrkRt196yvWZdPl9ZttY63GxzvmWT+QfdnO11O+XUT7o4hDTDqT9PkvNur5Jr3fHWTbQ2OGAOxbmzjLH44SvN7SKwHLaKAKDeq3ur9zrtgWmUI+c0anNbvHHSoyqT7N+XV13uH6vukMJtZqfT9JGLWiK1w4/X+4e17dzW1x7K9OBOtLMfaunC07tsJK7ojn96p2+vjOTN/pAM6bd9/ljMLvJ7xtLftjoGI3uT/zpg+49ZuzMVJjZ4vdFx6gro+d2l54+wcuG7l1Gw9ZhFHR9oeLN6N+M0TAz6/dt38NGIzpojqWgHRAzL4IpbU8fPbutTs6vKz/ito2das51vnERXXtdOlVy+3QO+0GVzPm9c8Rtq+T9emfXu9dd+5vrlMvGXpofNWrHTr+8udXudbrobfI0iKtvdWWUAytpucvENFQ17/udaxWMVosXqj2tWjAz5/Llx4LzWadRPZrjdPAe//PBKv2sznBu0F8nFVa78uxxCgAdk437sNQTaAZNX2dMvaVuLxu6hsJJeezglW1CCCGEEEISgottQgghhBBCEoKLbUIIIYQQQhKCi21CCCGEEEISgottQgghhBBCEoKLbUIIIYQQQhKCi21CCCGEEEISIjHPtohsBfDfANYDUACXquqHROQSAK8HMBRnfY+qfrtVWSlR9KYjD+6q9LS37WDJ+XovzDlf721zJS/fnDpHbc24Uq1XGwAu6nT+57KRP3514gIv30tWXV9PXzNzqrdtf7W/nu7POM/tfdPrvHwX9j9UT2/NOn/2aMV3tl47aRzXOdf+e6bWe/m6Mq7Nq/Mz9bT1agPApk4n1vzJ2Gn1dCpwIT846USyT1l/fz290Yo5AVwz7uI7q9d3Mt84vrWe3tY9Vk/PVn3f7GjRCS/PW+skqD/e4fdtT5eTjnblnDz1lFW+fz2fdtuG51x/7h7y/dn5TtNn3c4D3tfhe5fvHHZ9/bSt99XTV+8+xcu3Za1ro/V2A8DDEwP19JoeNz4TY77sc3LYvT7hVFfexrV+v+/Z68bnotPcvL1x7xYv37o+N2cO7HHtl1zgPzZ+5XVnDNXTOx/25485fCAFd4x09PuO9U7Tt+Hcwpz7np8fcPsVAs1p0Tioe09y7Z8K/Nk65+LQSXesd4z41xPmTjAxmfaX9wTC1T7nhK2ucfnSw/68rQ64gLM9ruzytJ9vboMrz3NpA8hMG/ewUdGmyqEL2tRbcP1S7vH7tmZU5eldTlYdeo3VxmiKyI379ZYG3EZbRmXAnz+du4wH3cTQMex7cyvWQT3bvI0w6fyYn6+w3pRhvN1dE36+ipkmPfe6oKq+lt/rw+yUK8PWA/iuYetuFl8hjIxxUFeMKz/0PQ/c5ean9YqH+/U/6Dpj8gR/TpuPQGTMaat7t98Xts0zW9y23JiXDeUet63zoHFQH/DziRqn9RnufTsegO9dtmVbVzMA9Ox182lymz9n0gXjk+423vK9fhnpsvFE51y+jnF/8hcHXEy5KbfP9OZgrpr5k7cfMYFj3rq1vRj84jw3uf0ItH7r8HVhjT/eXaO2fDvP/KDSxh9eaTIGANCzz5zfzJgUVvvBWwd1dtalQ293x6Qrr5Y1Xvp+v7yU6Sc7F9L+0g3dB1x54ceIdXDbMbFtD8kPu4MwM+fHVMu0+PGAoyTJH7WpAHinqt4oIr0AbhCR78fbPqiq/5hg3YQQQgghhCw5iS22VXUfgH1xekpE7gKwOan6CCGEEEIIWW4syj3bIrIdwPkAro3feouI3CoinxSRwSb7vEFErheR66dGy42yEEIIIYQQsqxJfLEtIj0Avgzg7ao6CeCjAE4GcB6iK9//1Gg/Vb1UVS9U1Qt7V2UbZSGEEEIIIWRZk+hiW0SyiBban1XVrwCAqh5Q1aqq1gB8HMBjk4yBEEIIIYSQpSKxxbaICIBPALhLVf/ZvL/RZHsJgNuTioEQQgghhJClRFSbK1KOqmCRJwL4KYDb4CRO7wHwKkS3kCiAnQDeGD9M2ZRt5/TpO78YXQB/6+BD3rarCu77wunZyXp6NHDu5MXpY74xfU49vTU76uVbm3FlDFX66ulvj57rx9Tp9ttb7Pe2ndnduDnFmn87zOO6HqinfzbjFHzbcr7G7oGi086VTbusmhAAtnW4mPaVXEyh+m9kzmnOTulzireOlO+vmqk4V9TNQ5vq6a19voIun3H31O+cWOVt6zcKvd6cS89V/b7IGOfXUMHFd/7qPV6+m0bcM7ZD407pl836GrKtg+P19P37Xfu7u3wdX7Xm5s8569243bjL1+edtdF5ryrq9rnz4Y1evmyH68PuTr+uiqmrXHbPJnfl/XwzBdfvxXHnA0tP+ePdfbIbB9uOUtF/7rk84/p61Xo3v8fHfd2dTjgXlWasns3/Tq45o2xa47R9+nCgzzOnlvSJvrKzOGO8V0YfKGW/Lqthq5ldKj3+eKPT6KGmTPs1UDkZd5TmjQ6sFijojL4sN+5iqnT558vcmBnTXqOPm/HLq5mQiuv84yxt+tfWNbfWb2Nmxm3LTLvyS4O+1qxj1MTbacYxUAl273PbSr1u28w2vzw1SkirdUsfou1zrzXV/HNFTV+k5/wy8s6AioI5bdkYACBj6s6aqTWzKchnyjcfDwhOnZ7izV6CSvsGUO+1VfNJMB3zw2aemfbOrg/aaz5+NLj0lZkx6rYuaZrPtiUza8Z0IJj7pmu6DrgXuSl/vEfPcAV2GfVfqIyzdXnKtKBa+1Hcvd/VZecc4Gsfc4EKz+oOO4z6Lj/hxz6zwVWWNUq//Lg/QKUeU54pY26VPzHs6WN2g3uRH/Xj6xxxZVgdXaj0s1rAzIyJPeiz6U1u0oQqvPyYUeuZfi/1+hOjc8ScZ0xVpUDBZ/vdjlXXAf9ZubnVLiar6ptZ75fXu8vtp0ZNaNPRa5e2ysHpjf7nV+8eV97MBn/d0L3X6Fb73H52fAEgVbHnwcbKQcDvz59/5Y9uUNULcYxI0kbyMxwyhQAALZ3ahBBCCCGErBT4C5KEEEIIIYQkBBfbhBBCCCGEJAQX24QQQgghhCQEF9uEEEIIIYQkBBfbhBBCCCGEJERiNpJjSX+qgOf33AEA+OzUCd623SWnmst3311P7yyv8fLdNefUdc/ovaOevq+4wcv344kz6+k/Xv/Devq8jXu9fN+ePrue3pgb97YdLDtl4N7iQD398Iz/y/RV43Cyur+7C5u8fJPGSzVacr6pVblZL9/NU1vr6e6M08mlxNfbdGedLueX+1x/vvzEm7x8d5edcnBt90w9fWC2x8uXNuXnM77WrL/DqeHuG3Eur839vj5w/0xvPV0oOb3P9UNbvXwHdpk+NHq6tZt9heOqDtc3523dXU/fM7zOy7dlYLye3j01UE+Xx/JevvFVbgwe3ufmXG9/wcs3PeX26+z31Un77nPt79/m2j96oM/Lt2aDc5QVgziaMXPQaPcCTZpkndtpdJ9TQq7dPO7lm9zh4isNOr2U1QCG1IxyULb4fVEpuNNLdX+Xty03ZZR5/S4+7fIVXWrKsHq1UBGYmXQeqapRE4Y+JKvPK2VdPqn4Ga0mr7TaxZeeCq5PmN0yBfeisMk/Dnp2uHZoxldldR5w+02c7eZM525fc5UyCrBqi2lR6jNatwmj6gsUYnOr3Tar/+rd4bfRavGswlB865r3WqqNlY1RTI3LC19nzHQq9/r5slNmmzkd5Sb9cUwbq6Y1rwanRG+O91i7bDB/Sn1WCWneH/DzVfMun9WOZfxTNvKjrtMqeb+yuTWmDHMqyY37wZd7G+fr2+kP0NygUdcZFV6xz5+P6250hRTWuAHJTfrlWXWfGNVjNee3I2cUfJmiK6OW8edZhym/3B1sM2227a3lwnyujKxRE86u9dvYNWT0eSbe/Kh//rFjlzPnjontgYZ1v9V+mnTB77NK3pVh650bDFTFY6YdM34ZmYI7t8yud5pYqTY/T2cKrl2VoG/Lne51307ntiz3+Aen1RZ2DLsDS6r+AV7pcm1Jl+x4h9pLd0Ka2eDaMXiv/zkyt8aV33+ffwDNbnYnQqv069njn+xSFRfH7DpXXm7S/4yudCa3JOaVbUIIIYQQQhKCi21CCCGEEEISgottQgghhBBCEoKLbUIIIYQQQhKCi21CCCGEEEISgottQgghhBBCEkJUm+tilguPOjerX/tWpPK7r+zr856Yd0q6/548sZ62SkAAeHKP0wJ+ZM/T6+mnrrnHy/ez0VPq6bP79tXTF3Q96OW7p7ixnr5m9CRv24bOKTTihPyI93rXnIux03i5bh7b4uV73Oqd9fRQyXmutubHvHwPzDp1W1WdZme63OHl68o43c14yXm9erJFL98dB50WcatR5BWrvh7naevuraevH/PVjJmUUw71mvJvOejrDaemXRwnbxiqp+95cKOXb/s2t80qDO/e4yscU2lXb0+Xq3d8vNvL19vnNEOT+5xfrGvdjJevWm38vbQ4ETjYjEKud6M/D6bHnf5Ojbou2+v3u1idXMa1o7Db959pp9uWHXFjUumpNc2XnnD5UiVfxVQzKrxah9HiVf18mSmjkysafdXawAVndGC5UV9tVdzk5mDXDudkmz3N74v8Q27ulk278sP+eGSnXdocVij3+uc3NVrEtGl/Ne/n69rjyrfqu2qXn8/vC/N+oHiz6rriqqCu/a6MjNGGzW4IvHOmyT0PuXyVLj/fzFajWpsxWrNhv95Kpxm7NY3HHgDyQ65iqx/sCBV03a68ojlN54e9bMhNu/1KPcHcKrhtVk2Y8g1dXt2zG12+nl3+HJxdb44z07dhebPrGysSu/aH6jajhBywGsBAt2l2s8o0TQXHnDmVdkwEdZnxmd7s2tEx5tdl29I1bJR2gWrNmGZR6nEvuoZ8TeXMBqPsNFa3bDCn7XFhNY3pYjC/D7oOnd7sChTfsofOIaMcXOdrL/MjLsbigPVP+mVY7WBuyvWn1ILjdsb0U9b1xex6/7PN9k3VaAarHX7FPQ+7DigNmE4LlldWOSlGR6eZ5tc9reoPAGppE4fR9tWyfkzZqaCDm2Db7ykrS37wVlVoFYb5UT++dNHVmxlz/VLt9T8r1cxP26ZUOfgcMWvUSrc/Pvl97sSvWfMZk/L7U8pGfdjjxicsL2tUgD/8+Z/foKoX4hjBK9uEEEIIIYQkBBfbhBBCCCGEJMSCP5cjImsBvB7AdptfVV+XXFiEEEIIIYQc/7Tz25RfB/BTAD8A0N5NQIQQQgghhJC2FttdqvoniUdCCCGEEELICqOde7a/KSLPTTwSQgghhBBCVhgLqv9EZApAN4ASgHkviqpqX8Kx1TnxnB5931ceBQBIix+vVdw9qXNnPf1QxQ/vhMxkPX1dcXM9fe+cr5Zbk3G6tqmaU9XU1P9ekjJup4cKa5puGy87P9Jc1dcZ3T28rp4+a+2Bero34+vP7ptwSr/f2vqLevqWma1evqLxSN0x6tq1oXvSyzdWdAq6lOnPg1M9Xr6c0c5t7XOawVARuCrnnFA7p1d72+aMJnCi6PqzWvP7czBvFHwlp3sbHvN1d6m069uT1zmn2P0H/DE4zegDx+fcGOwf6UczakbHl+nwdUblCaNPNEo7yft3VvX3u77oyPpljE66fs/nnWJoZocfkxr1mna5MiQdKL8OuJhSpqpKoKfTHrcxPermoGaaa/Eys2Z8AhNTedC0OWs2Bhqujt1OsZT1TYqYPt3pwGTWKZtSRb+QbqPgK5luqgZ6uuyk20/N3+sKm/0x6L/TbZzd5MroGPHrteq24qBRGE6Eei23zerAKl1eNqixUlm1HOAr2azCMDjlQExTbF+ob1VEtrF5FDPb/IEcuMvFO2WMnd17/P2sUsyq+Wr+6cxT0FltXarS/DMm7Z9KUOpz+3UOuXhnNoWqR6NwnHPvl3v98SkOuHT3HjOOM0FM5qUdx1AZl511rwuDRitY8PNZ/d3cgNGkjQdaTtOscnfwGWP6rWz0jvkxvwz7kTizzk2Gnv3+ualmzWhmU7Ev6FvTxt4H3YE7tyYf5HMTstLpjqvCGn9C2r6wmrhyj58vO+2CquaDvii6Nlc7gwlv6yq4MuwxHJ6bSn0u3sys2aca6DFNjJWO5tcmc5OujPywm5CF9Z1+vnF38Fe7mt9YUO529XbtK3jbakYTaMe+NOAfkNlp19elXrctU/DnRcd+d9KZ3e7WTR1Dc16+apcrIzvuttVy/nhUu12+tNEWpmb9E19xndPw2pisihDwFZaZmUCDmHd1ZybdyaQ06M/VtNEJSsmlUxX/WLJz/Kff/pNjqv5b8DYSVe1dKA8hhBBCCCHkUNq5Zxsi8kIAT45fXqWq30wuJEIIIYQQQlYGC96zLSJ/B+BtAO6M/71NRP426cAIIYQQQgg53mnnyvZzAZynqjUAEJFPA7gJwJ8mGRghhBBCCCHHO+3+guSASTd/wowQQgghhBBSp50r238L4CYR+TGiZ3qfDODdiUZFCCGEEELICmBB9R8AiMhGABfFL3+pqvvb2GcrgP8GsB6RWOlSVf2QiKwCcBmin3/fCeAVqjrWrBwAOO/ROb3y25HarRzE+49DT66nt+ZH6+mzOnx/1fWzJ9XTF3btqKdXp30n2c1z2+rpuwtOn9cV+Lq6UkbhEzi6+jNGhTfnlHSTFV9HUzKqvhSaj8OOKafTu3D1w/X0vjn/jwwVE0fGeI+u3+0rAlf3uTb3ZF07xuZ8TdHouFMBPuGkB+rp63af4OVLGx1fKuWrdFZ1OW3Rw3tdO7r6fK2Q3a8w65R2Ya8M9rm+HRk1qsIgY3a3K6N6ootBa4Hizby0UyvMl8kZzdVB10+pkp+vOmBUfbO+Eim/yfV7YdLNhVQ2UHSV3X6dDzp9Xqnf71tbt9VclYN86YKbF7WMVab5sVt9YNVMVQnUf9WcK8Oq+abOKHv5MiNufudH/boqdqpZ7Vo+1BG6dM/DRq3X7WVDdtKo1ta4fKUBP/jMrO0zl+7d6ddbMGXkTNlWCxcVYmI3dsjOIb+8mU2mvHG/iDlnAEXPQ2Z8yn4Zo+e4Mlbd7ral/G7HzEY3Jr27m6jQEGjnulIm7eezmjir6guVfpk593p2rSvP9l8UiEtWc35/Wj2d7euOiSB4g6bMOO7yXYJjpxs9pumnztHgx5BNiGI2TW31j+HOEaOxG3YFVrr9fFaXmB9x59iRM/3PgPyYGe+S309W+Ta31u1n9X4AkDd907XHnVdLgzm0RTg8NTu3jH6wz683Y1R95V53rFfy/ph2jBulX6ebF537/c+AUr+LNxUo+DJTrq9rHS6O3MFpL19xg5OnpaquX0p9vhYvN+7K07Q5JwT1eupH8wERKiErVotnFHS1bDAvTOy23lTJn99Wf6eZ4LNo2sWeKprPm0BjV+0y/Wm2FTb4n/OlHldX34Pu8zU17a95xLRfTbtC9Z/dlpozsRb8E1V5tVH/Tbq5IHN+vuqgOyFJ2T9uU1NuP80brW3Ov45stYOVftf+zPgs/IyuL753+18fU/Vf09tIROSM+P/HANgIYHf8b1P83kJUALxTVc8CcDGAN4vIWYiuiv9QVU8F8EPwKjkhhBBCCFmhtLqN5B0A3gDgnxpsUwBPb1Wwqu4DsC9OT4nIXQA2A3gRgKfG2T4N4CoA/Dl4QgghhBCy4mi62FbVN8TJ56iq9/ceEck32KUpIrIdwPkArgWwPl6IA8B+RLeZEEIIIYQQsuJox0by8zbfa4iI9AD4MoC3q6r3u+Ea3TDe8GZlEXmDiFwvItePjDS/X48QQgghhJDlStMr2yKyAdFtH50icj7cYy19ALqa7ReUkUW00P6sqn4lfvuAiGxU1X3xg5cHG+2rqpcCuBSIHpBspz5CCCGEEEKWE63u2f41AK8FsAXAP5v3pwC8Z6GCRUQAfALAXapq978CwGsA/F38/9cPL2RCCCGEEEKODxZU/4nIS1X1y4ddsMgTAfwUwG0A5u8DeQ+i+7YvB7ANwEOI1H+jDQuJ2XD2Kv2tzz0DAPDY7h3etq3ZEZfPeKk2pn29zbv3X1RP/2r/HfX0F4cv8vLNVJ0u55zevfX0ttywl++WGacI3DM34G0bmXNKmxduuLWe7g80gz8eP7Oe3jvjNH7P2XC7l+9no6fU0zV1GqC90776b1Wn09iMFtwfH3o7fB3WdKmxEmrEqP4AYKBvtmG+kNEJ195q0dcArV8/UU/PFF29MxP++Oic2y/T5yuHLKsHnOqpWHbfFWfn/DaVh43eZ7VTaFltHwDIoKtLR41yMB+qwVy/d61141i+t8/LVt1qVETjQT83O9SCm7k05TJmptIN3weAVMUqq9z7NX8IUBlwG7sfdH02tyYIyJRvtYDlnlBJZhSTZkpXgyc5yr2mvEJzZZ5V0FWCurz9zKaOcb84z75p8lWCv8Hl3HRE2QxdqM9Lm6dUMgVXYMekf0tbYbVRdJlQOwLdXanHbew+6M+tuX6j0TKXPzJFv4zsrHs9s97VWxzw+7Znt9HTGe3a7JpgYhis9jGs1+rQskb3NrvO16mZUyc6R2smnz/BO8ZNfGMVb9vMeldm56jbVur1Y7eaPKsPTAexW71cacAFOLndv860+nZz3BolWzXv12v1appx7ap0BW0cbawFTBUDPVunKT88HMsub6XL5esY9yermHxza905rHu3f/6ePsEdDFZbmCoG6lGjp5vZ6Pqs/+4pL1+l3+hVO4zqccz/vPHUcFZ3F9Q7u9Gdm8M2qnG05g66k47m/XGcPsF9FnXvcef9UDtX63TzzNP9Beuhqinf9kt2MtDimfKrvUY3OefPb6m5sap1mLI7/XZkxlzsCJWdVq1X8su3lAddf2ZH3Vyo9vkn6swBd1LUzo7m+Q6au4DLrt7aql4vX8XsZ9WE6Slf9QijTyycOFhP53f780zmjEpxwHe+psbdXJAZ12fVDav9fHNmvIzeTwr+XLUa0e/d94Fjqv5b8EdtVPXLIvI8AGcDyJv3/2qB/X4G7yPV4xmHEyQhhBBCCCHHIws+ICkiHwPwSgBvRbR4fjmAE1ruRAghhBBCCGnLRvIEVf1tAGOq+pcAHg/gtGTDIoQQQggh5PinncX2/I0wsyKyCUAZ0S9KEkIIIYQQQlqw4D3bAL4pIgMAPgDgRkSPcfxnkkERQgghhBCyEmjnAcn3x8kvi8g3AeRVdaLVPoQQQgghhJDWP2rzdFX9kYj8eoNtMD9SkzgbMtN495roRyu/PrPd32Z0fzcWN9TTWfGVOJmU8+esTTu1zOvW/cTL9x8HnlpPb8k5I+G3R8718k2XnSLn4lUPett+MHNGPT1acaqabx04p2lMq/NOYXPlwbO8fOs7Xbxl43V7+sZ7vXwHin0N97nl4Ca/3rSrt1wxajn15TFWGThbdqqksSnfp+btNeersg7sckqfVLcbk/Xrx718Y79cX0/X+l295bEOL9/IrrWuXmOOqvQ1/5VRq/vTTKDP2+80RbrO1Zvf4WuPSoOu/OIO18+do36fFZ2wB+VV/hzsecD14cw2F3zHQb/PNO1izE248qe3BbEbm5XNVwuO6lTZqKIq9n0/X6bQWGOXn/PbOLvR9UXnAbdPdtrLhpwxRdV8SxwqxsDY95Arb2qLf2db1Qx/xtiwchN+X9jyrL7JavAAIDfj6ho3Srbuvf78qWVdGbbsYl/zO++sMq4SaBBz0y6OUrdfRrZgt7kyrOoPAGbXGL2aUQv27vIVamrGv2zqKq7yx3HNLW6+z61xA9S1z9dhFVe5bZlZV1e24M/bTNGVb9Vt6UAHmp1028p9/sTomHTlW51cR6AILPW5MnsechOjsMHveFt+2mj31l3vT9ZKj1HBGZXe7Fr//DN4u5vUMyc6VWp2MlC8GZ1cbsxpx6qB4i1dcO3NjfhqtNJaN/Fs7FLz50Vm2pVvW1/u89WjfXeNoxFWgwcAKaM3HLjDKeM07Y9jquRiz8y4Ma1l/fmdGXftqnWZ8Rj31YT2UyU97OvftNvo5GZdedVOXzvXf+N+t22VGx+r2QOA9KTp65KLXbv88a71udcdu801xkyg0TTKwIyJPVTV2c9Yq8LLHPDHXvOuXqu+i2I0o1w1561s0EYzJlbj57UdAIyO0PZFetLPZnV/tv2pCV9pnDGaRikY5V4tdBi6Puu6y42bjSHMJ8H4YMrUnTVza9QPXs2c0fVOCyjFQOHY0ViLfCxodWX7KQB+BOAFDbYpgEVbbBNCCCGEEHI80nSxrarvi5O/p6rVZvkIIYQQQgghjWnHRvKgiFwqIs+If4KdEEIIIYQQ0gbtLLbPAPADAG9GtPD+cPxT7IQQQgghhJAWLLjYVtVZVb1cVX8dwPkA+gBcnXhkhBBCCCGEHOe0c2UbIvIUEfkIgBsQPfD8ikSjIoQQQgghZAUgqto6g8hOADcBuBzAFao603KHBDjlUV36D187HQBwsNLnbduWHamnb5vbUk8/oes+L9+NhRPr6Q7jPLt9ZrOX7+GZVfV0Spyqpifr63dW5Zy26Mbhrd62F2++pZ6+Z9Yp7e4ZX+/le8LaHfX0nZPuRzmHCr4uqFpz34lGJ9y2tYO+Hmmy4PQ+J65y2sJixX8O9qERp+M7d9PeevrWvb4isLfLtbloFIGzs75+pzrryl+13lfujD3s6oJR/0kmUK1NOW1PdszVVV4VPJtrpqt0Gm3Ww76yJ2U0ZGYYIcF0r9rdzFfPUr8fX/6g22jLKA74BdqHGirdfhl997p2FU23ZHwDlle+VbyVu/xHJjqM/m74fJfu3O9/hy6ZGPvvd++n/SkNY5VEpdNo7Hr8fBljTbMasuKgH1/3PrctXfL7qdTj8qoJ1+r9wjK7DhhlZY9fV6nXKLVMu2z/AUCq4l4XjEqv66A/VnYMrAYx1K5ZJZ1VnqXKfnmlAXeMWDUh4GvibL2VvJ/P9mGm0EJ12eniKKxy6VX3BGq5XnPMzZhjM4i9uMYdJPmDrnPDvqgYrVs1b5SQU74Wz6oEU8G8yMwZndyU03LVcs21c1YNVsv557rSoIs9YzR7qVJzDZk9x6Tm/NirRqcnZi6VAoVhbtzFblV6tUxwfcsMcaro15UyajhPSRco1FJT7qApb3InFquWA4BKv1MJ5vaMubKzzaVkVo1W6/OVr7DnVRtT2W8HUqbNaaMXDcZKpq3bM3SFmvE2j45poPSzqrnqoPusTM36OjmrIMwEmkGLWsWfaaMU/fJ8HZ+JNeWPt0yak6dRKWp3p5+vYE5iab8MNXo6r7y8/7mseZvPLdu8WAFgzCkNpcesPULloFXwmXbpjP8BJv1mjWbnRbWFZyOTaSufBjGJbbOJT2eDD5IOo1K07Q81g6Zvv/vgP9+gqhc2D/rwaOcXJM9V1dC4SAghhBBCCFmAdm4j2SAiPxSR2wFARM4VkT9LOC5CCCGEEEKOe9pZbH8cwJ8CKAOAqt4K4DeSDIoQQgghhJCVQDuL7S5V/WXwXqVhTkIIIYQQQkiddhbbwyJyMuLHRkTkZQD2JRoVIYQQQgghK4B2HpB8M4BLAZwhInsAPAjg1YlGRQghhBBCyApgwcW2qu4A8EwR6QaQUtXmnpyEKNRyuL2wpeG2C/I76+nrxrc3LeOmiW319LNX31ZPr8v5zTm9a389fdeMU+HdOupr8U7eMFxPv2bbNd627wyfU0/vme6vp7f2jnv5bhh1MXUbteD0nK/w6cw5PU025+7g6c35GpzhCedomyk7hU0q8N0Vh51m6Lopp0Ts7PdVUZMzTpFT3u+0T9oRqKd6XHyj+/u9behwGp/0sNMt1ToCZV7NqOZ6Xfmbtw97+fYeGHD7jLk2VnxzEqprXD91DBntWtqvN11w9ZZ7tOH7UXwubVV12Rk/39x6197seKDgM0YkNUapdAlN81kFX+dQMI79btu669y2ye1+eeuud8FPnOQq7jrgl2fblZsyfRHYkaxazyr8enb75ZWNRarrgK9zkqqLY87o/Yqr/P7s3WUVby45l/X7dvWdLshyj9tW7vLzZYquEKv7yxT9OT3X7+JbdafTa02e7Gs5Lflhd/zMbgpUXrb5gW7VOzxNumPc77PMrHvt6QP9LkN2wvVF927XrjCmrt1O2aVZ195Kp6/Zyx9wA54qG33eTKDhsgrDmjvW07P+BJcBoxyc9e9ItDq94mp3/skf8PViavRvqZIrI1T1pQvlxvsE6rZap4vJtgMVf16kp01bjP6s6/5A2GWUcbUuox8c98+xVmlXG2g+t9QoA1MzfuzVNe6Ekd3rlH6V9f65OPewfy6dRwKVoKe7s/nmghNBCxVeW7RQ0GHNKm+TGEWbjS81FvSn0RhmDpryAgUfTH9qp9HCTfpmYymY8tP2pB2UZ/vQzDNMBMulLncM6qTZVgw8rB3+GsCLadocC610d6Z8TbnYJdAq1ozeUadd+6XbVz1K1u2nZTcenuoPvp5P51z/pXr8+V0zdVndn3T7+XTK9FPan5tacrGLVVgGOksx41U7MOTe7wwWDuUjmMdt0nKxLSKnA3gDop9sB4C7RORSVb03sYgIIYQQQghZITS9Z1tEHg/gKgBTiG4j+TiAGQBXicjFixIdIYQQQgghxzGtrmz/BYBXqepV5r2viciPALwPwHOSDIwQQgghhJDjnVY2kpODhTYAQFWvBnBSYhERQgghhBCyQmi12G71IORMi22EEEIIIYQQtL6NZKuI/GuD9wXA5oTiIYQQQgghZMXQarH9Ry22XX+sAyGEEEIIIWSlIRo4X49ZwSKfBPB8AAdV9Zz4vUsAvB7AvOjwPar67YXKOv3cvP7HFVsBALcUTvC2PaHrvnr6gfLaevpDO57h5fv9E6+upx8qrqmnD5Z7vXz9GeeqHCk5b/U9k+u8fMWK+56yrXfM23bPqMs70OnKe3D/Gi/flrVuv4mCc45W1Rfnzky7bY87cWc9fc3tp3j5kHOuz7TxW2vNL2+g390FNDrkHJnpMf+7V2qz83lWDjofZWat7zatDLv4uh/2PZizm62g2iVr3YFDeNTVnTKq3JqvBEXGeK2tn7o4EHq7XTo75fapBfpSMXV1jLt0NefnK6x35XfvMf0ZHD7WkZ3ym+h5rFNG55kOFKv9D7qgZta7/ix3++Noy+iYdA2udASOcKsNNt7u3t2+k9g6qW2sncN+vrlVmYb5uvf7jtLJba4TB+/1HbBjp7v5lC65APMjfqdV866C7JQfRzOsk7iWDSTUpi9yxkdd7vPnfnbS1VXudduyM34MVeP7tv1c7fDv0MsPNXdw5yatJ9octxm/DOu4tp7tSpcfe9o4wzNTxnkblFfNu/2sj7raGXh4c24OZq3XOB3MR+PT1lzz6zhScPnKa3v8baYPM+NmzpSDsTfe3Lmtzifdeb/vkta8m4MaupZtvbb8EXdelrzvgi6f4D5jsjsPurJ7fCex51K3aQmOTePZ1tANbNzS1sethcCnPGBOOunGjmwAEOM/Rs34macDt7Qtz8ZeCU9oTdoVerpt7MbpjFrg3s+ZsQrqkm7jpzbua+t+BgCdcZ9ZYr3OlWD+mNhrU86jL53+eEvGjIltYzCX1I5V2RwHVd9hLsazbccg7AsYLzZSwTnM0swDHmDnTOjF9sowc9D2MxDEbr3n4e8GBB7vOjU/n5ox8fo54x8HtUnjsK/6/WTd3d4YBMeZ3Rb6vj1M+d+b+OQNqnph88yHRzs/136kfArAsxu8/0FVPS/+t+BCmxBCCCGEkOOVxBbbqvoTAKNJlU8IIYQQQshyJ8kr2814i4jcKiKfFJHBZplE5A0icr2IXD8R/GmZEEIIIYSQ44GmN9aJyL/hkDtSHar6B0dQ30cBvD8u9/0A/gnA65qUfymiX67E6efmk7mxnBBCCCGEkARpdWX7egA3AMgDeAyA++J/5wHINd+tOap6QFWrqlpD9PPvjz2ScgghhBBCCDkeaHplW1U/DQAi8n8BPFFVK/HrjwH46ZFUJiIbVXVf/PIlAG4/knIIIYQQQgg5Hmjl2Z5nEEAf3MOOPfF7LRGRzwN4KoA1IrIbwPsAPFVEzkN0G8lOAG9sJ8iJahe+NXEeAODxPfd72z498iv19NP776qnz1m1z8t3V2FTPV1Wp8hJw1fz3DW1oZ7e1Dnh9qn6Wp2hSaesqgWqvmLF5Z0uuT8CWNUfADx8l6nrdKeRmin6fzioTTmVzoGCUxVu2OY/fzo+7fRT6bRr19wOX284WnTxpSat/stvR3WPK69jyijOhnx1TnXQ1TXzKF8XlBpybcmNuzJSB/w/qsxudffld93r4isNeNlQ7nd1iSlv4N4gn9HkTZ7qyu7c549jz57GdygV1vp9sfo2l0+sNiv425AaJVnvbv9Zg8ltru7Be532qJYLVI9G95edcXVZ1SEAZApm7hrVUS0wQHVMGRXcXpcOtXj5MRdvLeO2acbP1zHu8uUmXFCVHl/51DlqlF+BvarTPIeRKrs2pqr+eOR3Oi1XtduVXw10dzbeVMm1sWuvrzWr9BuF5aSbq7Wcr27LTBRMPtMXgdJOjBYvPeHKK633jxExar1U2deLddzvjv3yltX1dG7XiJevusYpu7TDxZHf6//Yb63D9ZPV/UnZn4+5YafU0g53nKbmfE2azDRW8NXWDnj5vL4dMTq1QuC2NKqwbKDysseWpyQLtX01N8ad4679oYJP5sxBM2zOl6Fmr8vs123SgYIuM2Tq6nQeUZn1z3s6YXRlRicngdbM6s8OUaMZhZyndQvK8MbH9JOnQgM8JZ+WzbZALVc76PSJYvupI/Cm2vExujcJ8tXG3eeo1fuFGje1uruiP2e8flJzDpv1NYhi5rEdAw3nme1P02dWAwgEmjg7PuqvGyCm380x4rUXgWbRzulaUJ6NN9Ab2m22Xamg37Xk5o90GnWinZtBjGr0kJ7qD0Hf2Hq7/GPOllEzcySVD+aPxcyf2thY02zh3LJqxWbtBQIdoRm7UM0YzsljSTuL7b8DcJOI/BjRr0c+GcAlC+2kqq9q8PYnDis6QgghhBBCjmNaLrZFJAXgHgCPi/8BwJ+o6v6kAyOEEEIIIeR4p+ViW1VrIvLvqno+gK8vUkyEEEIIIYSsCNrxbP9QRF4qIi1+M5QQQgghhBAS0s5i+40AvgigKCKTIjIlIpML7UQIIYQQQsgjnQUfkFTV3oXyEEIIIYQQQg5FVBf+ccb4Z9VPRfQDNwAAVf1JgnF5nHFuh378ii0AgLuLm7xtDxbX1tPndj1cT+8urfby/WTk1Hr6/IFd9fQtE5u9fKWa+/5hdX/dWV9FtHtqwOWr+K6109c4lde1955YTw+u9rVC07NGHWVu0hHxx6Q4bfIVjLZvrvmdPdlJ90eLSrdfXi1jNXbu/dyE/4eOcnfjfLWcX151lVMd5R/yVUe27sy0izcz68dbMfYgMfYdCcw8OWM5m13vyu59yM83t8bVZcvITvqxV7pcvp49LmM10PFZJZndli4Hui4zkN37fFdfpduN3exalx64359btZwbh6nNTllkFX4AkJkxGkQTRi1Q9Vm1oNXs5Q/42qy59U6XJEbBV+7250Xvg05fVTPqO4R3mplzS6gFTBfdIJe73TGXnfK1c3Nr3HzqNhrA9Lh/LJU3DtTTqZKdQH5M6YPj7oUd0zX9fr5R88c7o4eqBrq7VMlovib8mCzabVRUwTlXJs1+gW7Lw+jB1FNZBedwoxSzWjirqgPgqeAwPO72SQd/8LTqMauMC8fbKOk8DVegqvMUdHOBFtCowsQo+DRQ8InViBktXJjPU7QZPVuou7NzQWx8YXlWx1dr8dlplHk64052tk2H0OJOTauMC8vw+tDMhdQq39Cr026eWT3dIWFYXV0rTZrFKNMOGW/b10bdFirorN7QU/MBqBkVoHdHa6iEtNtMX3jqN/hKOhuvBspB2xc2Xy3MZ9ts6g2Vg948s3MkHA/T755+EEDVnC/E9q2EHlozdk3qDWM6JF5Dyur57HgEY2Vp1WdePlOGtjiuwrpsvLYvWik2W7XR9uEPqpfdoKoXNs98eCx4ZVtEfg/A2wBsAXAzgIsBXAPg6ccqCEIIIYQQQlYi7dyz/TYAFwF4SFWfBuB8AONJBkUIIYQQQshKoJ3F9pyqzgGAiHSo6t0ATk82LEIIIYQQQo5/2vkFyd0iMgDgawC+LyJjAB5quQchhBBCCCGkLRvJS+LkJfFPtvcD+G6iURFCCCGEELICaLrYFpFVDd6+Lf6/B8BoIhERQgghhBCyQmh1ZfsGAApAAGwDMBanBwA8DODEpnseY6qawngt0h1NVH3t0bP6bqunvzVxXj29b67Py3dyz1A9fcfUxnp657j/naIv7/Q0ByacYvyiLf6dMyNjPfV0R77sbbtp95Z6Otvpto3t82OSDqfm6eh2SqlyKRiWqtH7GPtSfsi/5b6wweiCjKEqNPOkjCInO+XSNd+OhNykUcYZM1Et7SuqUvudHqnvYV8PNbvOxZgfddumt/ix9+w2GjZj9coGNrWUUdKliy4ODdrY/4DT+2RnXb3lbj9jddqosipGizfij6nV9lk1IQJbl1XXza3xO9Sq+gbvdkq2at6PKTPl6u7f6faZXe9rFdMlt63S5coI9XmZKTen1WjdNOOPQff9Y/V0rdsNQi4XdK4Zg8yEU2gdoqAzGq7MQV/zVet3x3Fu57Areo1/jPTd6LYdMpENqaJRO2VdvvTknJ/RKqZ6XBtTc76mEaWyyediTQ+N+/maqFNrgUowdXDMvPD7XfvNTxmMm34aCM4XVvE26uKQzkAXaPV8VqcWtFELRv1otWF5XxGo4xMNy7NKOyBQ9Rm1Wqgz1EmjBQx1d515l8/oAw+pK9QTzuebDfIZdZsWXXulI2hjqcUJ02L0Z1q0GkD/vNdM4lebbK6HbNYmAN6cCfvC6nttf1aHhr18Vofm9UspmPtWGVdqoXpsprgL9WwFM2/tvAgUbLaulno6ba4ZbKYnrIWKyWYE+jzbLi/dQp1ot0momLTHT4s5p0bRWg3mTDMFH6SFMs/EEY63p9PT5nrH2pwZu1aayiZ9lgrOU7WCr55tFoOnBaz4n8t2vJrVG73RpG8OmdMt9JZHSdOjW1VPVNWTAPwAwAtUdY2qrgbwfABXJhYRIYQQQgghK4R2bCQXq+q351+o6ncAPCG5kAghhBBCCFkZtGMj2Ssifwbgf+LXrwawN7mQCCGEEEIIWRm0c2X7VQDWAvhq/G9d/B4hhBBCCCGkBe2o/0YR/YokIYQQQggh5DBYcLEtIqcBeBeA7Ta/qj49ubAIIYQQQgg5/hFtpkSZzyByC4CPIVIB1n0qqnpDsqE5Np09oL/7hacCAM7t3OVtGzcqwB+OnllPv3jtjV6+B4vr6ukv7Lignl7V7WuUamoUOSZdrPpqnslZp6gKVX3lWafmSU25bTrga2sy+4x+yRSfG/N1NKV+N0aVAae0SU/5MfU87ParGENix7g/xpMnmX0eaq7Pswq+rgOujOJgoP4zmqJMYPMpDri8PXudVkdqfkxq7YYZo4wr+vlKPe7OJ6v0s9o+AMjYbWWr/vPHSozGrus+p8qa2+4rIaudrnNse0OlUGbWqI7KgX7IkCoYzVUl0A0ZdVRljVNMZsb9zrWKO093NzLl58ubeZYz7a+2UEWVyk23wSrKrGorE0ygkXFXXqB/s2o9T4UX9Kd2ukkoRiGm5SC+1QMuPWQ0e61UTlZ5FbbX7mf0UlZvd0i8LfRntRH3swSpNav9IuYCPeH8+7P+eItRfsGo68J8ajRnqZ5uF0OgxUsNOD2hLeMQNZjRxFnF3SF6NqvhKrZQrRndllYqwSazzc6tUAdm22+VcTlft2nL9zRxoXbO1GUVfOHnoy3DqzfU9qUaK8kOodZEnxfgqduCc6ftm5Tpl1qrY1hte4Pj1o5jxfatrx719HdNjpeWhMem3a9dBVu7dbWq+0jLaKfslvnMOIYKulYxNSu/1VrOlt9uvjCGduttpgVcYK15RBzrulLuWPhB9bIbVPXCIyvoUNp5QLKiqh89VhUSQgghhBDySKGdr3TfEJHfF5GNIrJq/l/ikRFCCCGEEHKc086V7dfE//+ReU8BnNQgLyGEEEIIISSmHRvJov0sOyGEEEIIISuJdq5sQ0TOAXAWgPpTgar630kFRQghhBBCyEqgHfXf+wA8FdFi+9sAngPgZwC42CaEEEIIIaQF7VzZfhmARwO4SVV/R0TWw/10e1NE5JMAng/goKqeE7+3CsBliJzdOwG8QlXHmpUxT1eqhPO7HgIAnJU74G37wP5n1dOPH9hRT3/4QV8Dfubg/nq6N++0VAenerx8z9h2bz19xW2PrqdP3DLk5dt/0Om70tOBOmmNU2fJOqf1Su3y9We5caetKWywWp1AZ2MeY+2/ww3Z3Fpfb1M0RrG+B9y2qW1+eXZbqc+9n532y+scdq/LXa6MgQd8pVSxz7U/W/D1QJk5F3zHqNuvlvWfza3lXPkdY04ppWk/9vxBN3ZizUmBZq/W4fqplnN15Yd9zVrm4KSry6jMOoYD7dqM2c/qybo6vHxqtkmgH0pNzDQuo9MvQ4xSKzPkNH4yNePnyzrNmVdG2depWSVfatoo3gqBp9HEK/k8mlEbnzD5jJrPaOaieo0mbXzS32ZVaUY9Fmrs4KrydH/S3+fn2+eOT08nVwt0VU2UUC31bAYJtHgeRvcmoZLKbKse8M8lXh9adV0Qu1X62THWUPFmNXF2jIM2Wh1hSy2eqdfrp9lAx5cx87HSXDvnKeTCuprsE+rudMYdC6luozcsBBpFoyvzlIaB1izV6eb7IWXY4pqU0aq9njoy0CV6/dlCGaftTU9/7ofltT33G1d2iM6xmU4uDLaZnq2VWq6VTq5VXUdCqzLaVeYddQyhzvEYtOto42g1jp6mMciXZD+Fc+lY11VLrt/bsZEUVLUGoCIifQAOAtjaxn6fAvDs4L13A/ihqp4K4Ifxa0IIIYQQQlYk7Sy2rxeRAQAfR/TDNjcCuGahnVT1JwBGg7dfBODTcfrTAF7cbqCEEEIIIYQcb7RjI/n9OPkxEfkugD5VvfUI61uvqvvi9H4A65tlFJE3AHgDAKzZlG2WjRBCCCGEkGXLgle2ReSH82lV3amqt9r3jhSNfge36Q03qnqpql6oqhf2r2pLmkIIIYQQQsiyoukqVkTyALoArBGRQbin9voAbD7C+g6IyEZV3SciGxHd/00IIYQQQsiKpNWV7Tciukf7jPj/+X9fB/DhI6zvCrhfpHxNXBYhhBBCCCErEtEF1Cki8lZV/bfDLljk84j83GsAHADwPgBfA3A5gG0AHkKk/gsfojyE7ef06F98JdLwXTlytrdtIOfUVjumnPtufM7X7G3rc4bBB8ZWoxlT026/SskpyfoHZr184/ucekzKvo4mN+r2s3q6jhG/ro4Jt7GwxpVR9U1wyE2gIZk5f+zEWGvSZbdNg69UHeONFVOZWV97M3mCU3T173Tap3K3/weR3FRz7ZVV4WXHjIYsUPhI0ZVR7XUaLikHsWZcY9LDTidXHfAVjqlZo6kymjmZCdRyVilmVGva0+Xns1q8otV/BSqiilV5BceWVWwZ1ZhOTqEpVuMXaMOQMspJo+U65Ji2yrze3ob7AIHWzOj4pMvvC52aapgvbK/YNoZ6unKTOZP1n8+wMYmtK+x327dWJRjUK2mra3N966nvgm2e5irVRGMGQM1YebHiMNRwVttXC3Vg5ljw1FstFG/NdF2N9pvPFsbuHSOmn5sp3cIYDqmgzZhsGalAr9pM0dWuGqxVvla6t2OtglsstVy7MSxlHCsR9u1xyw/0Szeo6oXHqrxWt5FcBGDX/EJbRH4bwEsRLZIvWWiRrKqvarLpGUcYKyGEEEIIIccVrW4j+Q8AJQAQkScD+DtEvxo5AeDS5EMjhBBCCCHk+KaV5iNtrl6/EsClqvplAF8WkZsTj4wQQgghhJDjnFZXttMiMr8YfwaAH5ltdPERQgghhBCyAK0WzZ8HcLWIDAMoAPgpAIjIKYhuJSGEEEIIIYS0oOliW1X/Ov7xmo0ArlSnOEgBeOtiBEcIIYQQQsjxTMvbQVT1Fw3euze5cBpT1gz2lgYBAJXAY5dLOd3WaMEpyob293v5hkad8qxWdBqpXE/Jy6cPuTJ0ndN1Td816OWTLqeoGrjLj2l2o9P75Iec+qdjzFdtVbrctu59bpsEdqBSt8tnmouefb5ObHatG86eXb7WzZIqO21Wqd/p/Tr2+wq67g7Xh5nxYsN09IZRl5V9JVctZxRyBaNxKwb93uU0cekJp+eTaV+5qH3d7oXRuqWHxv18c3OuDKOT00oTZRjg68Smpv1tHc7HqOWSSVe8bKkB12c2hmijUdINO5mPp9wLEKuOqvh1WV1drWDam/bnY830k4ybP0oF6jdbnlXmacn/Q5Yt39Z7SOy2b4K6PJ2c2SZBG60arlYqNHz/kPLDbTabUfBJxs2LWqgm9OJtMWe8fYxus4lWL8wX127ia1cN1qL8pnW1pwjUcAzaKvswsPtpm33bTPXXquwjzXek25rRUiVoFY5ttvFYQx1dcrBvScyCP9dOCCGEEEIIOTK42CaEEEIIISQhuNgmhBBCCCEkIbjYJoQQQgghJCG42CaEEEIIISQhuNgmhBBCCCEkIY6LX4Is1jJ4sLAGAJBL+XqkHdNr6unZotPYITDu6IhTt60/dbienrh2nZcvbYrPzLnyMjN+efl73PeUjilfqZUpuG3ZWaM1q/lBWcVf55BTj2VmfPVWZsIqz5xGqrihx8s3eMt4Pa2dTmsmBV9rljI6vfSoyRfoz7rvMUq6UI3mBeiUdghVeCZdM9o5GfDVjDg40rCMQ8RJM24g1GizQn2epF1Mtelg8JognZ1un2KgNzSvre5PUuJlq42NN8wXolU30cIyvPLMnDmkrlmjRfT0Ymk0w9Z7yLawza5iP1+5uarQy2dftKluO0Q710ybFmrSbL5262qzHUfE4Si/jkQPdqz3oaIsOVr1bbtKQ0LIcQ2vbBNCCCGEEJIQXGwTQgghhBCSEFxsE0IIIYQQkhBcbBNCCCGEEJIQXGwTQgghhBCSEFxsE0IIIYQQkhDHhfovLTV0ZyIt2e2jG7xt+YxThc2MO3VbZiTr5UvPOTVY4Uqn+6ut9rVMnWMuX+ew0/aVO33tWmbO7Zeb9PVN+RGnySv1uy7OD/mqsdSs0eml7PuBkswqA9MuY/6O3X4+o7vDtBnaOV/pplWjIyy6unRm1suHrCnDxhDo4zydXCrQxNm6jbqueuAg2iFUwUmmyZRN+7o7u5+n4FNf0+hh9YahkitlyjdlaC1o72zQh5YmGrswJK+NJo5WoXvlhfq8oyXU7LXLsdCatauko7qOEELIMoVXtgkhhBBCCEkILrYJIYQQQghJCC62CSGEEEIISQgutgkhhBBCCEkILrYJIYQQQghJCC62CSGEEEIISYglUf+JyE4AUwCqACqqeuFSxEEIIYQQQkiSLKVn+2mqOtxOxslSHj/adSoAIJv23b1Dd66tp63iuWPU92Kr2daz1wmLs1N+vp59zrUsVefu7Z0qe/nSk84fbV3VACAVF2POuKplasaPyTqojav6ELd0l/OHS9b5w2uTU355pgwRv11ePuugNv5sCVzVKNidjFs68GxDUg3ztaRdL3LQjqYO6WPhlm7lk27mjD4cB3WbbT7mnmxCCCGELBm8jYQQQgghhJCEWKrFtgK4UkRuEJE3LFEMhBBCCCGEJMpS3UbyRFXdIyLrAHxfRO5W1Z/YDPEi/A0AkF3bvxQxEkIIIYQQclQsyZVtVd0T/38QwFcBPLZBnktV9UJVvTDd17XYIRJCCCGEEHLULPpiW0S6RaR3Pg3gWQBuX+w4CCGEEEIISZqluI1kPYCvxraMDIDPqep3lyAOQgghhBBCEmXRF9uqugPAow9nn77cHH516z0AgK/cer63rXPcXZzvGHPvd+/zlWw9DzntnmbcPv1js35lVmtntX2TvrbPKuk0UPBVi0U0JFDmedq9lCkvVOsVnINPMk79pxVfR2i31cpNYgB8nZ7R0Wkzvd1CHI7+7rDLblMRSAghhBCyDKH6jxBCCCGEkITgYpsQQgghhJCE4GKbEEIIIYSQhOBimxBCCCGEkITgYpsQQgghhJCE4GKbEEIIIYSQhFiqn2s/LMYLXfjqbZHyb+CGDm9b/4NOf9cxNFdPp2Z9LV5qYrqe1lmn+9PCnJdPOvP1dG1i0uWr+Qo6T9UXbMMRKPS01ma+cunwt1nVH0CdHiGEEELIIsEr24QQQgghhCQEF9uEEEIIIYQkBBfbhBBCCCGEJAQX24QQQgghhCQEF9uEEEIIIYQkBBfbhBBCCCGEJMRxof5Lzwj6r4uUf2tuLXjbsnc+VE9bjV+t4Oertau7M1rAVrSr6lsWUPVHCCGEELIk8Mo2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJcVyo/zJDM1j30WujF7Wqt63aID8hhBBCCCHLAV7ZJoQQQgghJCG42CaEEEIIISQhuNgmhBBCCCEkIbjYJoQQQgghJCG42CaEEEIIISQhuNgmhBBCCCEkIURVlzqGBRGRIQAPLWEIawAML2H9yyUGYHnEsRxiAJZHHMshBmB5xLEcYgCWRxzLIQZgecSxHGIAlkccjMGxHOJYDjEAyyOO5RAD4OI4QVXXHqtCj4vF9lIjIter6oWP9BiWSxzLIYblEsdyiGG5xLEcYlgucSyHGJZLHMshhuUSB2NYXnEshxiWSxzLIYYk4+BtJIQQQgghhCQEF9uEEEIIIYQkBBfb7XHpUgeA5REDsDziWA4xAMsjjuUQA7A84lgOMQDLI47lEAOwPOJYDjEAyyMOxuBYDnEshxiA5RHHcogBSCgO3rNNCCGEEEJIQvDKNiGEEEIIIUmhqo+4fwA+CeAggNvNe6sAfB/AffH/g/H7AuBfAdwP4FYAjzH7vCbOfx+A1xxmDFsB/BjAnQDuAPC2JYojD+CXAG6J4/jL+P0TAVwb13cZgFz8fkf8+v54+3ZT1p/G798D4NeOYFzSAG4C8M0ljGEngNsA3Azg+iUakwEAXwJwN4C7ADx+CWI4Pe6D+X+TAN6+BHH8YTwvbwfw+Xi+LsW8eFscwx0A3r5Y8wIJn6sAXIBovt8f7yttxvDyuC9qAC4M8jfsawDPjt+7H8C7zfsNx7ONGD6A6Bi5FcBXAQwkGUOLON4fx3AzgCsBbFrs8TDb3glAAaxJMoYWfXEJgD1w543nLva8iN9/azw37gDwD0s0Ly4z/bATwM1LcIycB+AXcQzXA3jsEs2LRwO4Jt7/GwD6Eu6LxNdX7fZHPX+rjSv1H4AnA3hMMBn+YX5AAbwbwN/H6ecC+E48GBcDuNYM2o74/8E4PXgYMWycH1AAvQDuBXDWEsQhAHridDaexBcDuBzAb8TvfwzA/43Tvw/gY3H6NwBcFqfPQrRg74gPhgcApA9zXN4B4HNwi+2liGEn4g+qJZwbnwbwe3E6h2jxvagxBPGkAewHcMJixgFgM4AHAXSa+fDaxZ4XAM5BtNDuApAB8AMApyxGXyDhcxWiL9oXx/t8B8Bz2ozhTERfyK6CWWw36+v43wMATkI0p28BcFar47yNGJ4FIBOn/970QyIxtIjDLhz+AG4OLtp4xO9vBfA9RL9JsSbJGFr0xSUA3tUg72LOi6chOkY74tfrlmJeBNv/CcBfLEFfXDk/fvFcuGqJ5sV1AJ4Sp18H4P0J90Xi66t2+6MeU6uNK/kfgO3BZLgHwEYzUPfE6f8A8KowH4BXAfgP876X7wji+TqAX13KOBAtJm4E8DhEUvf5D7LHA/henP4egMfH6UycTxB9O/1TU1Y9X5t1bwHwQwBPB/DNuMxFjSHeZycOXWwv2pgA6Ee0wJSliqFBTM8C8L9L0BebAexCdKLLxPPi15Zgbr4cwCfM6z8H8MeL1RdI6FwVb7vbvO/laxWDef8q+Ivthn1tx8nmQ4vjvN0Y4m0vAfDZpGNoI44/BfDRpRgPRH8NezTMOSzJGJrMzUvQeLG9aPMC0WLsmYsZwwLHiCA6j526BH3xPQCvNGP5uSWaFxNwzwhuBXDnYoyJ2f+Yrq8Otz9UlfdsG9ar6r44vR/A+jg9/4E/z+74vWbvHzYish3A+YiuKi96HCKSFpGbEf3p5/uIvlGOq2qlQZn1+uLtEwBWH4M4/gXRAqYWv169BDEA0Z9grxSRG0TkDfF7izkmJwIYAvBfInKTiPyniHQvcgwhv4HoFg4sZhyqugfAPwJ4GMA+RON8AxZ/XtwO4EkislpEuhBdBdmKpRuTY1Xv5jh9tPFYDjeGVsf54fA6RFeXliQGEflrEdkF4NUA/uII4zji8RCRFwHYo6q3BJuWYk68RURuFZFPisjgEcZxNGNyGqLj9VoRuVpELlqCGCxPAnBAVe9bgjjeDuAD8dz8R0QL1iOJ4WjnxR0AXhSnX47o/HkkcRx2XyS0vjrs/uBiuwEafVXRxahLRHoAfBnRfaCTSxGHqlZV9TxEV5cfC+CMpOu0iMjzARxU1RsWs94mPFFVHwPgOQDeLCJPthsXYUwyiP4E91FVPR/ADKI/dy1mDHVEJAfghQC+GG5LOo74g/pFiL6AbALQjeg+vkVFVe9CdJvClQC+i+j+x2qQZ9HGZDnUu5wQkfcCqAD47FLFoKrvVdWtcQxvWcy64y+A74Fb5C8lHwVwMqJ7hfchun1isckg+mvYxQD+CMDlIiJLEMc8r4K7WLHY/F8AfxjPzT8E8IkliuN1AH5fRG5AdFtHaTEqXQ7rq3m42HYcEJGNABD/fzB+fw/ctzAgWpDuafF+24hIFtFE+KyqfmWp4phHVccRPVTweAADIpJpUGa9vnh7P4CRo4zjVwC8UER2AvgColtJPrTIMQCoX02Fqh5E9NDVY7G4Y7IbwG5VvTZ+/SVEi++lmhfPAXCjqh6IXy9mHM8E8KCqDqlqGcBXEM2VpZgXn1DVC1T1yQDGEN0DuFRjcqzq3ROnjzYey+HGMILm47kgIvJaAM8H8Or4w3PRYwj4LICXHmEcRzoeJyP6QnpLfA7dAuBGEdmwiDEAAFT1QHzxpgbg44jOnziCOI5mTHYD+IpG/BLRX0vXLHIMAOrnoV9H9CDfPIsZx2sQnTeB6ILJkY7H0c6Lu1X1Wap6AaIvHg8cYRxt90XC66vD74+F7nVZqf9w6D1FH4B/4/w/xOnnwb9x/pfx+6sQ3Vs7GP97EMCqw6hfAPw3gH8J3l/sONYifoofQCeAnyL68Poi/IcQfj9Ovxn+Q2iXx+mz4T/osAOH+XBiXM5T4R6QXNQYEF057TXpnyO6krrYY/JTAKfH6Uvi+hc1BhPLFwD8zlLMT0TPDtyB6FkCQfTg6FuXYm7CPWS1DZHlYGCx+gIJnqtw6EM+z20nBvP+VfDv2W7Y14iuNu6I35t/4OnsVsd5G/3wbES2gbVBvsRiaBLHqSb9VgBfWqrxiLfthLtnO7EYmvTFRpP+QwBfWIJ58SYAfxWnT0N0G4As9rwwc/TqxZqfDfriLgBPjdPPAHDDEs2L+fNnCtG653VJ9gUWYX11OP2hqo/MxTaib1b7AJQRfQv+XUT3Av0Qkd7lB6ZDBcC/I/omdhv8D5bXIdK+3A+zIGkzhici+hPGvDLqZkT3gi52HOci0u3diuje1Pknpk+KJ9P98eSef7I7H7++P95+kinrvXF892CBJ3NbxPNUuMX2osYQ13cLnAbxvfH7iz0m5yHSNN0K4GuIDvJFjSHevxvRlYR+895i98VfIlrc3g7gM4hOyos+NxF9AboznhvPWKy+QMLnKgAXxn37AIAPo7FqrlEML4nTRQAH4D/M1LCvEZ3f7o23vTc47g4ZzzZiuB/RQurm+N/HkoyhRRxfjvvwVkRas82LPR7B9p3w1X/HPIYWffGZuJ5bAVwBf/G9WPMiB+B/4jbcCODpSzEv4vc/BeBNDfIvVl88EdFzLrcgumf5giWaF2+L23UvgL+z+ybUF4mvr9rtj/l//AVJQgghhBBCEoL3bBNCCCGEEJIQXGwTQgghhBCSEFxsE0IIIYQQkhBcbBNCCCGEEJIQXGwTQgghhBCSEFxsE0JIwojIe0XkjvhnrG8WkcclXN9VInLhYeS/OP5565tF5C4RuSR+/4Ui8u4FdieEENKCzMJZCCGEHCki8nhEPxT1GFUtisgaRC7g5cSnAbxCVW8RkTSA0wFAVa9A5EsmhBByhPDKNiGEJMtGAMOqWgQAVR1W1b0AICJ/ISLXicjtInKpiEj8/lUi8kERuT6+0nyRiHxFRO4Tkf8X59kuIneLyGfjPF8Ska6wchF5lohcIyI3isgXRaSnQYzrEP0QBTT6ue07431fKyIfjtM3m38FEXmKiHSLyCdF5JcicpOIvCiB/iOEkOMaLrYJISRZrgSwVUTuFZGPiMhTzLYPq+pFqnoOgE5EV8DnKanqhYh+kvjriH6O/hwArxWR1XGe0wF8RFXPBDAJ4PdtxfFV9D8D8ExVfQyiXyZ9R4MYPwjgHhH5qoi8UUTyYQZVPU9VzwPw53E5P0f0628/UtXHAngagA+ISHf7XUMIISsfLrYJISRBVHUawAUA3gBgCMBlIvLaePPT4nulbwPwdABnm13nb9+4DcAdqrovvjq+A8DWeNsuVf3fOP0/iH6m2HIxgLMA/K+I3AzgNQBOaBDjXyH6+eErAfwmgO82aouInArgA4huOSkDeBaAd8dlXwUgD2Bbi+4ghJBHHLxnmxBCEkZVq4gWo1fFC+vXiMgXAHwEwIWquit+KNFeUS7G/9dMev71/Llbw6qC1wLg+6r6qjZifADAR0Xk4wCGzNXzqKDo9pPLAbxeVfeZ8l+qqvcsVD4hhDxS4ZVtQghJEBE5Pb4iPM95AB6CW1gPxwvZlx1B8dviBzCB6Ir0z4LtvwDwKyJyShxLt4ic1iDG583fLw7gVABVAONBtk8C+C9V/al573sA3mruNT//CNpACCErGl7ZJoSQZOkB8G8iMgCgAuB+AG9Q1fH4KvLtAPYDuO4Iyr4HwJtF5JMA7gTwUbtRVYfiW1Y+LyId8dt/BuDeoJzfAvBBEZmNY3y1qlbn198icgKiLwOnicjr4n1+D8D7AfwLgFtFJAXgQfj3nRNCyCMeUQ3/6kgIIWS5IyLbAXwzfriSEELIMoW3kRBCCCGEEJIQvLJNCCGEEEJIQvDKNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCZFZ6gDaYY1s0BJKjTeKuGTD7U1ftCzr6MppM2/TTdJWNV6GBfM2ync47fDRI6n7MKput216WGUeSRw2rxxa51GX2XxT03oSaGfDupLuz0Z1L1KdjfY5puN6jPbVxeqPI9pXj7y+I42x4X4LjFy7p9MF6pGFZ0iUb8GKXDkLfWY1LurQOBaqs1HsC30ULLhPo/flcPZptx2t+6vdcpqNn8u6QD0Ltu0I4mywf7OxnN+/3bL9vIvTtgXnx+HEVj9uG5faukz73qHvSoNXYa4bbi1+T1Wf3bDyI+C4WGyXUMLj5BmAuAvxkpofiQbv2fe99+JOTaUOeS/Km2rwXryPqadRmQ3LabZ9Pt2onAVib1iOzWveU297o3qOfJ+G+ZptlwbbG+5z6Ht2/0Z1NqunXpbXNhyaN9Vge6NymsbZPF/TMr22H149C8bZNPb59xbqr4ViO/w4FtwHDfIeQduOOo4jKRMLbD+CPrQLliMp08XeoBw0yOfVow23N+qPet6mfXjodmm3HjTeLg3LbFGPqcsuKNypr3U9/mGph273TsGN6jm0zlSD7Skcmq9ZmQ33t9vRep9Uo30axNHuPnY/f3utaT6/zFr9vXTDOGqH7JNeaLtX56FxpNFoH/OeKT9db4eJEw1ibxBTw3pg23tonf57h/aNX+ah9TRqe9qLs9a0jV7bmsSZbjDW82Wmm8yf+ff97Tj0PTjS9e3SZLvE77nt83m998Tun2qwPeVtq7+/8b41OIbwNhJCCCGEEEISgottQgghhBBCEoKLbUIIIYQQQhKCi21CCCGEEEISgottQgghhBBCEoKLbUIIIYQQQhKCi21CCCGEEEISgottQgghhBBCEoKLbUIIIYQQQhKCi21CCCGEEEISgottQgghhBBCEkJUdeFcS4yI3A5gbqnjIImyBsDwUgdBEoVjvPLhGK98OMYrH44xMKyqzz5WhWWOVUEJM6eqFy51ECQ5ROR6jvHKhmO88uEYr3w4xisfjvGxh7eREEIIIYQQkhBcbBNCCCGEEJIQx8ti+9KlDoAkDsd45cMxXvlwjFc+HOOVD8f4GHNcPCBJCCGEEELI8cjxcmWbEEIIIYSQ4w4utgkhhBBCCEmIZbXYFpFni8g9InK/iLy7wfYOEbks3n6tiGxfgjDJUdDGGD9ZRG4UkYqIvGwpYiRHRxtj/A4RuVNEbhWRH4rICUsRJzly2hjjN4nIbSJys4j8TETOWoo4yZGz0BibfC8VERURquKOM9o4jl8rIkPxcXyziPzeUsS5Elg2i20RSQP4dwDPAXAWgFc1OEH/LoAxVT0FwAcB/P3iRkmOhjbH+GEArwXwucWNjhwL2hzjmwBcqKrnAvgSgH9Y3CjJ0dDmGH9OVR+lquchGt9/XtwoydHQ5hhDRHoBvA3AtYsbITla2h1jAJep6nnxv/9c1CBXEMtmsQ3gsQDuV9UdqloC8AUALwryvAjAp+P0lwA8Q0RkEWMkR8eCY6yqO1X1VgC1pQiQHDXtjPGPVXU2fvkLAFsWOUZydLQzxpPmZTcAPol/fNHO5zEAvB/RRS/+wvPxR7tjTI4By2mxvRnALvN6d/xewzyqWgEwAWD1okRHjgXtjDE5vjncMf5dAN9JNCJyrGlrjEXkzSLyAKIr23+wSLGRY8OCYywijwGwVVW/tZiBkWNGu+fql8a3/H1JRLYuTmgrj+W02CaEPIIQkf8D4EIAH1jqWMixR1X/XVVPBvAnAP5sqeMhxw4RSSG6NeidSx0LSZRvANge3/L3fbg7C8hhspwW23sA2G9NW+L3GuYRkQyAfgAjixIdORa0M8bk+KatMRaRZwJ4L4AXqmpxkWIjx4bDPY6/AODFSQZEjjkLjXEvgHMAXCUiOwFcDOAKPiR5XLHgcayqI+b8/J8ALlik2FYcy2mxfR2AU0XkRBHJAfgNAFcEea4A8Jo4/TIAP1L+Ks/xRDtjTI5vFhxjETkfwH8gWmgfXIIYydHRzhifal4+D8B9ixgfOXpajrGqTqjqGlXdrqrbET178UJVvX5pwiVHQDvH8Ubz8oUA7lrE+FYUmaUOYB5VrYjIWwB8D0AawCdV9Q4R+SsA16vqFQA+AeAzInI/gFFEk4McJ7QzxiJyEYCvAhgE8AIR+UtVPXsJwyaHQZvH8QcA9AD4Yvx888Oq+sIlC5ocFm2O8Vviv16UAYzBXSQhxwFtjjE5jmlzjP9ARF4IoIJozfXaJQv4OIc/104IIYQQQkhCLKfbSAghhBBCCFlRcLFNCCGEEEJIQnCxTQghhBBCSEJwsU0IIYQQQkhCcLFNCCGEEEJIQnCxTQghjwBE5LUi8uGljoMQQh5pcLFNCCGEEEJIQnCxTQghywwR2S4id4vIZ0XkLhH5koh0me0pEdkpIgPmvftEZL2IvEBErhWRm0TkByKyvkH5nxKRl5nX0yb9RyJynYjcKiJ/mWAzCSHkEQEX24QQsjw5HcBHVPVMAJMAfn9+g6rWAHwdwEsAQEQeB+AhVT0A4GcALlbV8wF8AcAft1uhiDwLwKkAHgvgPAAXiMiTj0lrCCHkEQoX24QQsjzZpar/G6f/B8ATg+2XAXhlnP6N+DUAbAHwPRG5DcAfATj7MOp8VvzvJgA3AjgD0eKbEELIEcLFNiGELE80eN0vIjfH/14I4BoAp4jIWgAvBvCVON+/Afiwqj4KwBsB5BuUXUF8/heRFIBc/L4A+FtVPS/+d4qqfuKYtooQQh5hcLFNCCHLk20i8vg4/ZsAvmkWwVeoqgL4KoB/BnCXqo7EefsB7InTr2lS9k4AF8TpFwLIxunvAXidiPQAgIhsFpF1x6xFhBDyCISLbUIIWZ7cA+DNInIXgEEAH22Q5zIA/wfuFhIAuATAF0XkBgDDTcr+OICniMgtAB4PYAYAVPVKAJ8DcE18G8qXAPQefVMIIeSRi0QXRwghhCwXRGQ7oivZ5yx1LIQQQo4OXtkmhBBCCCEkIXhlmxBCCCGEkITglW1CyIpDRDpF5GoRSS91LEeDiJwhIteISFFE3hVse7aI3CMi94vIu837J8Y/anO/iFwmIrlDSz7quLaLyO1x+kIR+dcjLOftwY/1fNv+UM/RIiLPF5G/OlblEULIkcDFNiFkJfI6AF9R1epSB3I4iEgmeGsUwB8A+McgXxrAvwN4DoCzALxKRM6KN/89gA+q6ikAxgD87mHWeVio6vWq+gdHuPvbAdQX26r6XFUdP5p4Ar4F4AV2QU8IIYsNF9uEkJXIqxH9wiJEpEdEfigiN4rIbSLyovj9vxORN8/vICKXiMi74p9C/0j8c+nfj6+2viysQETOE5FfxD9r/lURGYyvRP/S5NkeWz0gIhfEV9tvEJHvicjG+P2rRORfROR6AG+zdajqQVW9DkA5qP6xAO5X1R2qWkL0S5EvEhEB8HREFhEA+DQiB3cY+yUi8hkR+V8AnxGRtSLy5fhn2q8TkV8J8l0T/xz86xuU9VQR+abp6/+K+/lWEXlp/P5HReR6Eblj/ifgReQPAGwC8GMR+XH83k4RWROn3yEit8f/3m768y4R+Xhc1pUi0jlfnojcGdf7hbj/FMBVAJ4fxk0IIYvFUV3RIISQ5UZ828RJqrozfmsOwEtUdTJeyP1CRK5ApMv7F0RXiAHgFQB+DcCvA9iO6IrxOgB3Afhkg6r+G8BbVfXq+FaF96nq20UkJyInquqDiH7h8TIRySL6sZkXqeqQiLwSwF8jugIPADlVvfAwmrkZwC7zejeAxwFYDWBcVSvm/c1NyjgLwBNVtSAin0N0NfxnIrINkW/7zDjfuQAuBtAN4CYR+VaLuP4cwET8gzoQkcH4/feq6mh8Rf6HInKuqv6riLwDwNNU1VMUisgFAH4nbpMAuFZErkZ0pf5UAK9S1deLyOUAXoroFzbfDeBEVS0Gt6JcD+BJAC5vETchhCQGF9uEkJXGGgDj5rUA+BsReTKAGqLF53pVvUlE1onIJgBrAYyp6i4ReSeAL6pqDcD++auuFhHpBzCgqlfHb30awBfj9OWIFtl/F///SgCnAzgHwPeji89IA9hnirSe7MXiClUtxOlnAjgrjg0A+iT+YRsAX4/zFeK+eCyAm5uU+UxEPx0PAFDVsTj5ChF5A6LPnI2IFvq3tojtiQC+qqozACAiX0G0YL4CwIOqOl//DYi+GCEu77Mi8jUAXzNlHUR0BZ0QQpYELrYJISuNAvyfKH81osX0BapaFpGdZvsXAbwMwAYcuwXvZYh+VOYriO5kuE9EHgXgDlV9fJN9Zg6zjj0AtprXW+L3RgAMiEgmvro9//5CdaYAXKyqczZDvPgOlVWHpbASkRMBvAvARao6JiKfQuOfkG+XoklXAXTG6ecBeDKAFwB4r4g8Ku6DPKI5QQghSwLv2SaErCjiq6lpEZlf0PUDOBgvtJ8G4AST/TJEV2JfBndl+n8BvDS+d3s9gKc2qGMCwJiIPCl+67cAXB1vewDRIvDP4Rbw9wBYK/HPr4tIVkTOPopmXgfgVInMI7m4DfM/4f7juD1A9HPtX2+jvCsBvHX+hYicZ7a9SETyIrIaUV9c16Kc7wOw98EPAuhDtLCfiPvzOSb/FBr/QuVPAbxYRLpEpBvAS+L3GiIiKQBbVfXHAP4E0ZjPX5k/DcDtLWImhJBE4WKbELISuRLRrQgA8FkAF8YPKv42gLvnM6nqHYgWe3tUdf62ji8jutf5TkT3At8IYKJBHa8B8AERuRXAeQCsYm7+Z9Qvj+spIVoA/71EP5F+M4AnLNQIEdkgIrsBvAPAn4nIbhHpi6/YvgXRvdV3Abg8bgsQLTbfISL3I7qH+xML1YPIeHJh/HDhnQDeZLbdimgB/wsA71fVvS3K+X8ABuOHGm9BdD/2LQBuQtTvn0P0ZWaeSwF8N7xVR1VvBPApAL8EcC2A/1TVm1rUmwbwP/EY3wTgX43V5GmIrCSEELIk8EdtCCErDhF5DIA/VNXfOsL9e1R1Or6a+0sAv6Kq+49pkMcBInIJgGlV/ceF8i5H4ivpn1PVZyx1LISQRy68Z5sQsuJQ1RtF5Mcikj5C1/Y3Y6NFDtHV3EfcQnuFsA3AO5c6CELIIxte2SaEEEIIISQheM82IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJwcU2IYQQQgghCcHFNiGEEEIIIQnBxTYhhBBCCCEJ8f8BGKTXYlgupRoAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 864x432 with 2 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"xlab_idx = np.arange(9,len(N),10)\n",
"ylab_idx = np.arange(4,len(S),5)\n",
"plt.figure(figsize=(12,6))\n",
"plt.title('t-test p-values and sample size (stochastic)')\n",
"plt.imshow(sim, origin='lower', aspect='auto')\n",
"plt.xticks(xlab_idx, N[xlab_idx])\n",
"plt.yticks(ylab_idx, S[ylab_idx])\n",
"plt.xlabel('Sample Size')\n",
"plt.ylabel('Standard Deviation')\n",
"plt.colorbar(aspect=50, orientation='horizontal', pad=0.2, label='p-value\\n(avg over 100 replications)')\n",
"plt.savefig('./img/pvalues_big_data_stochastic.svg', dpi=500)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"class welch_ttest:\n",
" \n",
" def __init__(self, a, b):\n",
" self.a = a\n",
" self.b = b\n",
" \n",
" def __driver(self, a_N, b_N, a_mean, b_mean, a_var, b_var):\n",
" a_var_N = a_var / a_N\n",
" b_var_N = b_var / b_N\n",
" tstat = (a_mean - b_mean) / np.sqrt(a_var_N + b_var_N)\n",
" denom = (a_var_N ** 2 / (a_N - 1)) + (b_var_N ** 2 / (b_N - 1)) \n",
" degfree = (a_var_N + b_var_N) ** 2 / denom\n",
" # two-tailed p-value\n",
" pvalue = 2 * stats.t.cdf(-abs(tstat), degfree)\n",
" result = {'tstat': tstat, 'degfree': degfree, 'pvalue': pvalue}\n",
" return(result)\n",
" \n",
" def evaluate(self):\n",
" a_mean = np.mean(self.a)\n",
" b_mean = np.mean(self.b)\n",
" a_N = len(self.a)\n",
" b_N = len(self.b)\n",
" a_var = np.var(self.a, ddof=1)\n",
" b_var = np.var(self.b, ddof=1)\n",
" result = self.__driver(a_N=a_N, b_N=b_N, a_mean=a_mean, b_mean=b_mean, a_var=a_var, b_var=b_var)\n",
" return(result)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"N = np.arange(10,1_001,10)\n",
"S = np.arange(1,31,1)\n",
"wt = welch_ttest(0,0)\n",
"sim = np.zeros((len(S),len(N)))\n",
"for n in range(0,len(N)):\n",
" for s in range(0,len(S)):\n",
" ttest = wt._welch_ttest__driver(N[n], N[n], 2, 2.5, S[s], S[s])\n",
" sim[s,n] = ttest['pvalue']"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAtYAAAFUCAYAAAAAgHuIAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA/40lEQVR4nO3deZhsVXnv8e/b5zDIPBlEQHEer4Ki4o1xvkaNQROJkZgEoxGTmESvxmiMNxIz6L0ZjAlxgGgkiSKIGg0x4ohD4hAmEcFZFBQEBGQ85/Tw3j/2bk6dptbq6n12dVed/n6ep55Ttdfea797VXWf1bt2/SoyE0mSJEnbZ2atC5AkSZJ2BE6sJUmSpB44sZYkSZJ64MRakiRJ6oETa0mSJKkHTqwlSZKkHjixljTxIuK5EfHZta6jq2moPyIujYgn9tznVyLisX322fZ7akQ8o9B2WERkRGzse799iIj/iIjj+l53yXY/FRFfG2G9n42I01bav6QyJ9bSFBll8hMRZ0fEr/ewr8dGxOXb24/Wr8x8QGae3WefEfEg4MHAB3roa9X/4MnMp2TmKX2u2/4hcc+B7T6TmfcZof9/Ax7QjqmkHjixliRNkxcC78wJ+HazST0rvkKnAsevdRHSjsKJtTQlIuKfgbsA/xYRN0XE7w9Z58+AnwJObNc5sV1+34j4aERcGxFfi4hnDWzz1Ii4OCJujIjvR8TvRcTuwH8Ad277uSki7jxkf++IiLe0fd8YEZ+KiLsW6n9zRPzlkmUfiIiXtvdfGRHfavu5OCJ+rtDP7d7qX3qWPiKeFxGXRMR1EXHWYk3ReENEXBURN0TElyPigYX9/Frbx40R8e2IeOFA22Mj4vKIeFnb1xUR8WsD7ftHxAfbfXwRuMewfbTr7hoR/xIRP4qI6yPivyPiwBXU8PsDNTyjfT6/3j7XrxpY/4SIOCMiTmv7Oy8iHlyoaWbg+fhRRJweEfsV1j0gIs5sa782Ij4TETNt223vsLTti6+lm9vn8LC27WkRcUG7zn8tcwb1KcCnBva/ISL+MiKuiYhvAz+zpL69I+Jt7fh8PyL+tN3mfsBbgEe2NV3frr9L29/3IuKH7ev7DkvG/BURcSXwj+24vqd9Dm9sX1P3jog/aJ+XyyLiSQP13PZajfaMebu/6yLiOxHxlMK694zm5+vH7bGe1i7/dLv6l9rj+MVY8m5TRBwaEe+LiKvb5/PEgSE6e+mYSdoOmenNm7cpuQGXAk9cZp2zgV8feLw7cBnwa8BG4AjgGuD+bfsVwE+19/cFHtLefyxw+TL7egdwI/BoYBfgjcBnC+s+uq0jBvZ1K3Dn9vEvAHem+YP/F4GbgYPatucu9gscBiSwcdgxA08Hvgncrz3eVwP/1bb9NHAusA8Q7ToHFer9GZoJcQCPAW5ZMjZzwGuBnYCntu37tu3vBk5vx/6BwPcr4/JC4N+A3YANwEOBvVZQwx+1NbwAuBp4F7An8IB2fO/Wrn8CMAsc067/e8B3gJ2WvraAFwOfBw5pn9e3AqcW6n8dzQR1p/b2UwPP8W19Ltnmz4FPt+sfAVwFPKI9/uPa7XYZst3u7XN/x4FlvwF8FTgU2A/45ODrA3h/W//uwE8AXwReuPR1NdDfG4APtn3t2T43r1sy5v+3HZc7tOO6iea1tRH4p3Zc/3DgeflO4bX63PY5eUF77L8J/GBg/AbXPbXtcwbYFXjUQJ8J3HPg8WNpf3bbfr/UHtfuQ7bdr91+r7X+/ebN245w84y1tON7GnBpZv5jZs5l5vnAe2kmstD8x37/iNgrM6/LzPNW2P+/Z+anM3MzzX/8j4yIQ4es9xma/8B/qn18DPC5zPwBQGa+JzN/kJkLmXka8A3g4SusBZqJ1usy85LMnKOZxB0ezVnrWZrJ0n1pJi+XZOYVwzrJzH/PzG9l41PARwZqp+3rtZk5m5kfAm4C7hMRG4BnAn+UmTdn5kVA7TrZWWB/monRfGaem5k3rKCGP8vMWZrJ/AHAGzPzxsz8CnAxzfXIi87NzDPa9f+aZpJ1VGEM/zAzL2+f1xOAY2L4pQ+zwEHAXdux+ExmFi/TiIhfBH4JeGZbx/HAWzPzC+3xnwJsLtS1T/vvjQPLngX8TWZelpnX0kz0F/d1IM0fPS9pn4uraCaYzy7UFm09/zszr83MG2leP4PrLwCvyczNmXlru+wzmXlW+3p7D3BH4PUDz8thEbEPw303M0/OzHma18lBwIFD1psF7krzh+imzBz12vCH0/zB+vJ2DJZuuziWpfokrYATa2mKtW9TL769/qrCancFHtG+zX59+5b3c4A7te3PpJl8fLd9q/mRKyzjssU7mXkTcC3NJSSvGqjtLe1k693Ase3qvwS8c+BYfnXgcoDrac70HrDCWqA53jcO9HMtzRnfgzPzE8CJwN8DV0XESRGx17BOIuIpEfH59vKG62nGaLCeH7UTqUW3AHvQTKo2Do4L8N1Kvf8MnAW8OyJ+EBH/LyJ2WkEN8+39xUneDwfab21rWjT4XC0Al9NMupa6K/D+gTG8BJhn+ITvL2jeIfhINJervLJ0oBFxBM34/1xmXj2wr5cteX0eWqjr+vbfPQeW3ZnyWN+V5qzxFQN9v5XmzPUwd6R55+DcgfU/3C5fdHVmblqy3dIxv2bI87IHw125eCczb6ms+/s0r+MvRpO28rxCf0sdSjN5nyu0L47l9SP2J6nCibU0XbY5E5iZv5GZe7S3Px+2Ds2k41OZuc/AbY/M/M22j//OzKfTTDb+leYShmH9lNx2djoi9qB5a/kHmfnnA7X9RrvKqTRnPu9K89b/e9vt7gqcDPw2sH9m7gNcRDORWOrm9t/dBpbdaeD+ZTRv9Q8e7x0y87/a4/3bzHwocH/g3sDLl+4gInZpa/tL4MC2ng8V6lnqaprLBQbP2t+ltHJ7lvePM/P+wP+keYfhV7ezhpLB52qG5lKPHwxZ7zLgKUvGcNfM/P6Q+m/MzJdl5t2Bo4GXRsQTlq4XEYuvrxe175oM7uvPluxrt8w8dci+bga+RfO8LbqC8lhfRnP2+4CBvvfKzAcsdrlkF9fQTIQfMLD+3pk5ONFdkw9NZuaVmfmCzLwzzeVDb4qBJJCKy4C7FN5tgOZyqEsX3yWRtH2cWEvT5YfA3Ve4zpnAvSPiVyJip/b2sIi4X0TsHBHPiYi927etb6B5q3uxn/0jYu9l9vfUiHhUROwM/Anw+cy8bNiK7YTqGuAfgLMy8/q2afHa2auh+dAezRnrYX1cTXPN8i9H8yG057HthwPfAvxBRDyg7WvviPiF9v7DIuIR7Rnhm2mujV3g9namuYb2amCu/UDZk4asN6y+eeB9wAkRsVtE3J/muuGhIuJxEfE/2ktIbqB5y39he2qoeGhE/Hw7yXoJzaTz80PWewvwZ7H1Q593jIinF+p/WvvBugB+THNme2HJOhuBM4B/yczTl3RxMvAb7fMSEbF7RPxMROzJcB+iud580enA70bEIRGxL3DbGfNsLvP5CPBXEbFXNB/KvEdELG7/Q+CQ9rW7eBb/ZOAN7R8CRMTBEfHThVpWTUT8QkQc0j68jubnZfBntfR74Ys0f3y8vh3bXSPiJwfaH0PzQWVJPXBiLU2X1wGvbt+m/r3COm+kOSt8XUT8bXud6JNorhP9Ac1bz4sfvgL4FeDSiLiB5tra5wBk5ldpzjB/u93fsLfmofmw3GtoLrl4KPDLyxzDu4Antv/S7uti4K+Az9FMEv4H8J+VPl5Ac6b5RzQf0vuvgb7e3x7fu9tjuogmSQJgL5qJ03U0lwz8iOZShm20Y/a7NJO262guW/ngMsc16Ldp3s6/kuYDnv9YWfdONJPOG2guufgU8M891DDMB2g+GHodzfP+8+0fVEu9sd3XRyLiRprJ9yMKfd4L+BjNNeafA96UmZ9css4hNNeGv2Tg8qCbIuIumXkOzfN5YlvXN2k+1FdyEvCcdiIPzfN5Fs0H9M6j+aNm0K/S/JFycdv/GTTXMQN8AvgKcGVEXNMue0Vbw+fb18/HgGUzoVfBw4AvRMRNNM/NizPz223bCcAp7c/pswY3av/Q+1ngnsD3aC7/+cWBVY6luTxGUg8WP3ksSSsWEe+gSR949VrXorqIOIHmA5LL/eEz8SLiXcDpmfmva13LNIuInwV+JTOftezKkkayI4TbS5LWkcz8pbWuYUeQzTcv/tta1yHtSLwURJIkSeqBl4JIkiRJPfCMtSRJktSDqbjG+oD9NuRhh+40tO3mwhn3TTl8fYAtC+XDns0NQ5fPZflvkLmFctt8ZbuFheFxtJnlmNrqGwyl7arbrF5/1fDd0na1/rq82dLxDZrivrq+4dNhu07H23VfPfdX36bc2PmYO9XR3zadg6b7fgdxVcev3531/tzXTMI7t51LWKXae9/NKo75qj69O+pxlUxEEXUdSrxh/pprMvOOy695e2ObWEfErsCnaSK9NgJnZOZrIuJuNN++tj9wLs0nkrfU+jrs0J344lnDviEZzt08fNOvbjlo6HKA724uf5nbDzbvM3T5NVt2L25z3abdim03bNml2HbL5p2HLt8yW5n4zw6f+AMsbBnellsqb0zMlacAMTt8u5nS93cBUemvy3YxP3Tx8v0NSyYGorLNTGVfpTpq9cVCZcLYZV+FY1qubWZ+eB1d+yuObe14q/11GKfaH1xd6qhMrKIwflCe/NVqqP+xWGksPY/VP3RXXkfMl5+sTvuq1VB7XVRrX8V9rbSGWlv1ddFxgtJlLCqKr90Jqa86hiutYRzbLVReaH3uB+j9Ut5OY9vheKv9jWGi3uG4zrr25Nq35VaN81KQzcDjM/PBwOHAkyPiKJp82Tdk5j1pMkWfP8YaJEmSpFUxtol1Nm5qH+7U3hJ4PE1AP8ApwDPGVYMkSZK0Wsb64cX264YvAK4CPgp8C7g+MxffkL8cOLiw7fERcU5EnHP1jyrvm0uSJEkTYKwT68ycz8zDab7O9uHAfVew7UmZeWRmHnnH/cvXFUuSJEmTYFXi9jLzeuCTwCOBfSJi8dN5hwDfX40aJEmSpHEaZyrIHYHZzLw+Iu4A/C+aDy5+EjiGJhnkOOAD27Of+UKw1Xwlsq4Uqde0Df9boxapV4vim69st1CosRTDB5DzHdoq/UW1rbC8UkMt/aH6qfzivir91dIVOiRrdKqv5+SP2r6qqSUdkjCqSR21JIwu/VVTNyrblY6ruq+V91c73i4JJNWUiY4JKeV0hUp/1TpKaTEdkytKr4su6Rnj2FdNp5SRnlNLamo/P32nZPSdkFKymgkp1TpM8dguk3BMy+n7mJcxzhzrg4BTImIDzZnx0zPzzIi4GHh3RPwpcD7wtjHWIEmSJK2KsU2sM/NC4Ighy79Nc721JEmStMPwK80lSZKkHjixliRJknrgxFqSJEnqwTg/vNibJJkvfKpzoZDIMZvlQ1soJInU+ptbKCeJ1JI/5iupG6X0j4X5cn9Z6a+Y/lH7kG2XhI8O6RkAM7V9dUjxqKdkrGw/y/ZXSI3oku6xXFupji7JH1Ab24799Z2QUk0H6LKvSluH/urpMytPf6illtQ+EV9MvKhuU9tX4YnsmO7SLVmj533VdEn4WMUaOqexdNmm7xSPrnV02WYVUzzK3XXszxSP5Y0jwWMcSSMVnrGWJEmSeuDEWpIkSeqBE2tJkiSpB06sJUmSpB44sZYkSZJ64MRakiRJ6sFUxO0BLBQyr+YL0Xnzlb8ZZivReXOluL3CcoC5StzeQoe2zHIsXT1ub/ji6BKpV9muS7Qb0Cmmr+u+usT3rWZ99ajAvqP9Cv1VI/86ROfVYuk6R/uVxqJjPF6pv1oEVYf4vs6xal2i8zrX3nPEXOlnpGt/fUfC9R0VWFF8jvs+3tp2kxKPV9IlNq/jvtZlPN6kx+CtYgRe5+e/I89YS5IkST1wYi1JkiT1wIm1JEmS1AMn1pIkSVIPnFhLkiRJPXBiLUmSJPVgKuL2ElgoZDnN5i6F5eVIvWpbIYqvFqk3X4nHq7YV4uwWKvF41NpKUXy1uL0OEW5dIvo676trxFwpbq8agbfyiLm+I/Wq++oaj1d8HjtE6tW26xCb1+xr5RFz9fi+DpFwnWPaCv3NlwejGJu37L66RA92aKs9j31H1nUd9y7bdIlA7DsCbxzRacbjbbVa8Xh9x7eNI3puwuPxVjUCr+trtyPPWEuSJEk9cGItSZIk9cCJtSRJktQDJ9aSJElSD5xYS5IkST1wYi1JkiT1YGxxexFxKPBPwIE0oVQnZeYbI+IE4AXA1e2qr8rMDy3X33whmmWeQmRdlv9mqMXtzRW2m6/0N1+L4psvt2Vpu1JsHiwTnTe8re94vL4j9WrbVSPrOrR1idSr91fbZuUReLU+6xFzK6+jU6ReZV/1mLuOtRe2q8a+Vce9FINWqaEW01baV8fYt3pkYaGxa39dIuG6ROd17a+mFD3YOTax58i6TrFvPdc36fF4XePbdsR4vK7ReOOI6SvoPR5vNSPwVjPaj/HmWM8BL8vM8yJiT+DciPho2/aGzPzLMe5bkiRJWlVjm1hn5hXAFe39GyPiEuDgce1PkiRJWkurco11RBwGHAF8oV302xFxYUS8PSL2LWxzfEScExHnXPOjynvWkiRJ0gQY+8Q6IvYA3gu8JDNvAN4M3AM4nOaM9l8N2y4zT8rMIzPzyAP2L18TLUmSJE2CsU6sI2Inmkn1OzPzfQCZ+cPMnM/MBeBk4OHjrEGSJElaDeNMBQngbcAlmfnXA8sPaq+/Bvg54KLl+kqShcLH20tpHdXkj4VaW6G/+co2leSPhUpiyEIh/SMrqRvkyhND+k7q6JKeAd0SPrrvq0MSRofaa8kf1ePtkuLRIT2j6W9l+1m2v0JT5/6qaSJdtqnVXmjrWl+xv/ImxXSPWn+VPqsJKV3a+k7WqKWqdE6GWKUEkq5j0WUbUzwGtum4r0lJ6yjpub7ekzpq+k7xWOWkjqHGlEwyzlSQnwR+BfhyRFzQLnsVcGxEHE4TcHUp8MIx1iBJkiStinGmgnwWhoZML5tZLUmSJE0bv3lRkiRJ6oETa0mSJKkHTqwlSZKkHjixliRJknowzlSQXs0zPJplluExeLW4vdlCRB/AloXhQ1KK4Vuubb4SnVeM1atsE9W2/rapbdd3fF+trb5Nl4i5cn/1OMDh++qyDdTHqVNUYM+RdV1q79pflzi2znF7pbYOMXdQi++rbNSlP+gWMdclLq5jZF2nse0audV3VGCf20D5+e/Y3w4Zj9d3NB50q90IvNGs1nGt5jGNiWesJUmSpB44sZYkSZJ64MRakiRJ6oETa0mSJKkHTqwlSZKkHjixliRJknowFXF7CSwUol5mc/gh1OL25hZqbcP/1qhF6i10bMv5QttCOR6PuZVH8dUj8FYexVePmOvWVuqzFtM2M7fyfXWNx+tUX8d4vCgcV3Wb2nNciqyrRup1aOvYXyFJs9lufviBRS35qUPUXTXmru/+urb1HWdXaKs/V6sYB1jTZV81xuMtbzXrG0cdBb3H4+2IEXiwase1qnGFY7LsxDoi7gi8ADhscP3MfN74ypIkSZKmyyhnrD8AfAb4GFA5BydJkiStX6NMrHfLzFeMvRJJkiRpio3y4cUzI+KpY69EkiRJmmKjTKxfTDO53hQRN7a3G8ZdmCRJkjRNlr0UJDP3XI1CqjUAs4XogC2F9I9aKsjmhfJhzxYSQ2ZLCR7AfKUtK6kblFI8OiR11NpKKROd+6ule1T2VU3xmBv+/HZNICltN1PYz/L9lRIUyttUUzc6pJNUkzqqiRyl5V37K7TV+qt80ru2r+In0Wvj3iXho2t/pbZCmkl1m45tnVM8eqwB6JZa0rWOaU3x6JqQMekpHiZ1bLUDJnXAhKR19Pw6G5eR4vYi4mjg0e3DszPzzPGVJEmSJE2fZS8FiYjX01wOcnF7e3FEvG7chUmSJEnTZJQz1k8FDs9s3h+KiFOA84E/GGdhkiRJ0jQZ9SvN9xm4v/cY6pAkSZKm2ihnrF8HnB8RnwSC5lrrV461KkmSJGnKjJIKcmpEnA08rF30isy8cqxVSZIkSVOmOLGOiPtm5lcj4iHtosvbf+8cEXfOzPPGX14rk/lC1MtsDj+EUmwewFyWr4CZnR++3VxhOcB8JR4vK1F8pbi92hfH16L4ihFztW1q8XjFyLryNl1i5Gr7qkUF1msvRJN1iNSrbVeN76vG2ZXrKG3XKQKv1tZzf9XYvFpUU62O0na1mKku0XnV/jocV8fIuk7ReX3H49V02VfXSLDKvnqPx+sSMbeaEXjVOqY0Hm+aI/BWs/YKI/BG1DXmsqPaGeuXAscDfzWkLYHH1zqOiEOBfwIObNc/KTPfGBH7AacBhwGXAs/KzOtWXLkkSZI0QYoT68w8vr37lMzcNNgWEbuO0Pcc8LLMPC8i9gTOjYiPAs8FPp6Zr4+IV9Jcr/2KTtVLkiRJE2KUVJD/GnHZNjLzisXLRTLzRuAS4GDg6cAp7WqnAM8YqVJJkiRpgtWusb4TzUT4DhFxBE0iCMBewG4r2UlEHAYcAXwBODAzr2ibrqS5VGTYNsfTXIrCwQePmgooSZIkrY3aNdY/TXPZxiHAXw8svxF41ag7iIg9gPcCL8nMGyK2fpAuMzMihl75npknAScBPPhBO03B1fGSJElaz2rXWJ8CnBIRz8zM93bpPCJ2oplUvzMz39cu/mFEHJSZV0TEQcBVXfqWJEmSJskoOdbvjYifAR4A7Dqw/LW17aI5Nf024JLMHDzj/UHgOOD17b8fWLYGYLbQNpvDY/A2L5QPbUslOm92YfhlJ/OF5QALlf5yrhx1V4rOi8o2pVi6pr/C8q7xeIXtqpF1lfi5akxfoa1LpF5tu3rkXyUGrdRfNVKvY1shGWhmvhIZ1CE6rxrt1nd8XzX2beXRed376zmyrtDWKTava1vf++o5Hq9zJNg0x+NNawReTd8Rc+OIilulGLyJiLmDyY+6W+WYu0mx7MQ6It5Cc03144B/AI4BvjhC3z8J/Arw5Yi4oF32KpoJ9ekR8Xzgu8CzVl62JEmSNFlG+Urz/5mZD4qICzPzjyPir4D/WG6jzPwsWz/wuNQTVlKkJEmSNOlGidu4tf33loi4M81VGQeNryRJkiRp+oxyxvrMiNgH+AvgPJpLnv9hnEVJkiRJ02aUDy/+SXv3vRFxJrBrZv54vGVJkiRJ06X2BTGPz8xPRMTPD2ljID5v7JJy4MCWHH4Imxd2Kva3pZIYMltI+Jifr6WClFM8qLYNXzxTSxKpJmsUUkaqKR4rb6vV0CX5o7ZdNWWkmpAyfLsuyR+d++u5rXOKRyFNpGtyRbH2DukeTR2Vxi5JGB366z3Fo5paMobEkOK+Vj4Wvad4dE0G2AFTPMaSJjHpaR2rlNQBE5LWMSlJHTtgIkdOytguo3bG+jHAJ4CfHdKWwKpNrCVJkqRJV/uCmNe0d389MyvnBiVJkiSNkgrynYg4KSKeEIPfRy5JkiTpNqNMrO8LfAx4Ec0k+8SIeNR4y5IkSZKmy7IT68y8JTNPz8yfB44A9gI+NfbKJEmSpCkyyhlrIuIxEfEm4FxgV/wackmSJGkby+ZYR8SlwPnA6cDLM/PmcRe1VBLMFr4dfTaHx+PNZvlvhlKkHsBcIVZvfq7cX1baanF73eLxyv2VI+vK/dXj+4Yv3zDbLbKuHu1XiLPrO26v0l9t3GNueHRRNW6vFg3UYbtaf7X4vlLcVZdtmrZCDV1i8wAKcYDV7WpxgF3i7LrG3BUj5rr2twPG4014BF61hM5jawTeolWNwJuEOLYpjrmblji7STfKNy8+KDNvGHslkiRJ0hQb5VKQO0XExyPiIoCIeFBEvHrMdUmSJElTZZSJ9cnAHwCzAJl5IfDscRYlSZIkTZtRJta7ZeYXlyyrXC0rSZIkrT+jTKyviYh70HyNORFxDHDFWKuSJEmSpswoH158EXAScN+I+D7wHeA5Y61KkiRJmjLLTqwz89vAEyNid2AmM28cf1lLaqAcn7d5Yaehy2+dH74cYPN8+bBLUXzzhRg+ACoReKVIPYCZ2ULcXsd4vGLc3uzKt2naCrFvPUfq1dq69leKkqtFzPXdX1Ri5DpF51UisjpF51VSoTpF53WJzaM+TlMbj1fZxni87d9XpzFcZxF46y7mDqY26m5dxtxN6XM1iurEOiLuAxxP87XmAJdExEmZ+fWxVyZJkiRNkeJp2Ih4JHA2cCPNpSAnAzcDZ0fEUatSnSRJkjQlames/wg4NjPPHlj2rxHxCeA1wFPGWZgkSZI0TWqpIPdYMqkGIDM/Bdx9bBVJkiRJU6g2sa59SPHmvguRJEmSplntUpBDI+JvhywP4OAx1SNJkiRNpdrE+uWVtnOW6zgi3g48DbgqMx/YLjsBeAFwdbvaqzLzQ8v1VYvb21SI29uyUInUWyifqJ+bGx63l3Pd4vZKkXoAMb/ybarxc4VYvfo2lfi5Qn8batv0HLdXj/arRLjNDY/ymanF49Ui8Er7qmwzU6ih2a5DdF4tjqtLdF7HeLxSHbXx6xxnV6pxFePxaoqRZtWx6DjuJX3H4xmBt5UReKOZ4ui0HTbqboqfk2lXnH1m5inb2fc7gBOBf1qy/A2Z+Zfb2bckSZI0UUb5SvNOMvPTwLXj6l+SJEmaJGObWFf8dkRcGBFvj4h9SytFxPERcU5EnHPdtb6lIUmSpMm22hPrNwP3AA4HrgD+qrRiZp6UmUdm5pH77rcW839JkiRpdMVrrCPi72g+NzhUZv7uSneWmT8c6P9k4MyV9iFJkiRNoloqyGLyx08C9wdOax//AnBxl51FxEGZeUX78OeAi0bZLgm2FE6ub87hh7BpbnhaCMCWufJhz80WUkFmy2fNo5IKErWEj0JbKY2ja9vMlpUnf0A5/aNaQ4fkj6bPQipIJcWjmgpS2K6UFrLcvqKQTlFPEuk5FaRLekZtX5X+Ssdb3a6axtFzYkgtraHS38SneHT5JP+kp3iMIwlj0tM6JiFpYkJSIXbI1I0JGVsNWM3EnWUsmwoSEb8JPCoz59rHbwE+s1zHEXEq8FjggIi4nOZr0B8bEYfTnAm/FHjh9pUvSZIkTYbaGetF+wJ7sTXhY492WVVmHjtk8dtGL02SJEmaHqNMrF8PnB8Rn6T51sVHAyeMsyhJkiRp2lQn1hExA3wNeER7A3hFZl457sIkSZKkaVKdWGfmQkT8fWYeAXxglWqSJEmSps4oAdEfj4hnRkQ53kKSJEla50a5xvqFwEuBuYjYRHOddWbmXmOtbMACwaYcHp+3aaGwfL58aJvnhkfqASzMF/5+qETqzWyptM0Vm5jZUlhei7MrbNNst/J4vFKkXm27+ja1+LnadsPji2ZqEXjVuL1CPF4tbq8WC1WK75ufr9RQi5+rxeOtPGKuUzxelxpq/VXj9iY8Hq9rJFhpX5MegdfsbPX2VSxhFSOyVjP2bQLi2HbImDuYiLHdYU1QZN00W3ZinZl7rkYhkiRJ0jQb5Yw1EbEvcC9g18VlmfnpcRUlSZIkTZtlJ9YR8evAi4FDgAuAo4DPAY8fa2WSJEnSFBnlw4svBh4GfDczHwccAVw/zqIkSZKkaTPKxHpTZm4CiIhdMvOrwH3GW5YkSZI0XUa5xvryiNgH+FfgoxFxHfDdcRYlSZIkTZtRUkF+rr17Qvu15nsDHx5rVberIYqxerfOD1++uRK3NztbblvYMjyKL2ZrkXrdovg2FOLsqvF4WyqRdYUovto2pRqaOkrxfZUaOkTqNduV4vFWHqnXbFdoq8X3dYnOq8WP1SLwatF5HeL2qjFJpTp6jserRqetZjxelziuSn9G4G2ndRZzB1MedTchYzi1jKxb94ozzIjYb8jiL7f/7gFcO5aKJEmSpClUO2N9LpA0XwhzF+C69v4+wPeAu427OEmSJGlaFD+8mJl3y8y7Ax8DfjYzD8jM/YGnAR9ZrQIlSZKkaTBKKshRmfmhxQeZ+R/A/xxfSZIkSdL0GSUV5AcR8WrgX9rHzwF+ML6SJEmSpOkzysT6WOA1wPvbx59ul62aBYJNOTz94+a5XYYu3zRXSwUZnvzRNA4/iV9L95jZXGmrpW6UUjw2rzz5A8rpH/VUkNq+CqkgHZM/otZWSN0opYUAUGkrJXwU0z2gW8JHLd2jYypIMUWhus0EpHiMIxWkQ7pC/bg6JB5MQFJHU0aXOnbMRI6pTd0wcWM0Jmtoio0St3ctzbcvSpIkSSpYdmIdEfcGfg84bHD9zHz8+MqSJEmSpssol4K8B3gL8A9A+Rs0JEmSpHVslIn1XGa+eeyVSJIkSVNslLi9f4uI34qIgyJiv8Xb2CuTJEmSpsgoZ6yPa/99+cCyBO7efzmSJEnSdBolFaTTV5dHxNtpvqXxqsx8YLtsP+A0mg9CXgo8KzOvW66vBYJbFobH6t08v/PQ5Zu2DI/nA1jYUo7bi0KsXi1ub8PmYtMybYV4vMo2GytRfMW4vVp8Xy1urxCPV4vAq0Xq1aL4StF5pdi8pq0SyTRX2K7veLzaNrXIqMpx9R6Pt1Aai46RVqXIsEp/naLioDzuffdX0bn2Yg0992fM3VbG2W1lZJ20Jka5FISIeGBEPCsifnXxNsJm7wCevGTZK4GPZ+a9gI+3jyVJkqSpt+zEOiJeA/xde3sc8P+Ao5fbLjM/DVy7ZPHTgVPa+6cAz1hBrZIkSdLEGuWM9THAE4ArM/PXgAcDe3fc34GZeUV7/0rgwI79SJIkSRNllIn1rZm5AMxFxF7AVcCh27vjbC5cLF4EFhHHR8Q5EXHOjdfObe/uJEmSpLEaZWJ9TkTsA5wMnAucB3yu4/5+GBEHAbT/XlVaMTNPyswjM/PIPfcbJbxEkiRJWjujpIL8Vnv3LRHxYWCvzLyw4/4+SBPf9/r23w907EeSJEmaKMtOrCPi45n5BIDMvHTpssp2pwKPBQ6IiMuB19BMqE+PiOcD3wWeNUqRCznDzYW4vVvmVh63l5vLJ+o3bhretmFTOW5vZkuxiQ2bVt5WjdTr0LZhS7cIvJnCdjFbi8CrRfF1iM4rxeZBNbKuWEctbq1LdF6lhmpMW4e4vWpMWy1mrLDduozAM+pue4tY6wpWl5F1klaoOLGOiF2B3WgmxvsCizPLvYCDl+s4M48tNFUn5JIkSdI0qp2xfiHwEuDONNdWL06sbwBOHG9ZkiRJ0nQpTqwz843AGyPidzLz71axJkmSJGnqFC82joiHRcSdFifV7TcufiAi/rb9anJJkiRJrVrc3luBLQAR8WiaDx7+E/Bj4KTxlyZJkiRNj9o11hsyc/EryX8ROCkz3wu8NyIuGHtlA+YJblzYdWjbjVuGp4XMbikfWlRSQWY2D1++sZLusfHWWlv5U+UbNw1vKy0H2LC5/Kn8Ulsp3QNgppbUUWirpnvMVVIDagkfhWSIqKaC1JIwCm19p3jUtlnNFI9askaXZINKf50TOYr76tBfz+kUE5G4ATtm6obJGpLWkdoZ6w0RsTg7fQLwiYE2v7FFkiRJGlCbIJ8KfCoirgFuBT4DEBH3pLkcRJIkSVKrlgryZxHxceAg4CO59f3fGeB3VqM4SZIkaVpUL+nIzM8PWfb18ZUjSZIkTafaNdaSJEmSRuTEWpIkSerBVKR7zOcMP57bbWjbTZuHx+3N37qh2N/GW8t/T2y8JYYu31CL1LulHCe1Uy1u79bh0VobNq08Ug9gZsvwSLjSclgmOq/UVonAq8fjdYjOq8Xj1SLmusTj1erbESPwukbMTXPU3TTH2RlbJ0kTzzPWkiRJUg+cWEuSJEk9cGItSZIk9cCJtSRJktQDJ9aSJElSD5xYS5IkST2Ymri96wpxezdv2nno8rilErd38/BIPYCNtwxfvtNN5airnW+uROrdUo732njr8Hi3DZvKsW/V6LwtcytaDtSj82YL23WM1OsUj1eK4Ws6LLeV4vFqkXpdo/OK3dX66xCdNukxd9MQZWdknSRpjDxjLUmSJPXAibUkSZLUAyfWkiRJUg+cWEuSJEk9cGItSZIk9cCJtSRJktSDNYnbi4hLgRuBeWAuM4+srT+bG7hq055D2zbdPDxub+NN5b8ZdrqpvK+dbxwex7XLjZXYvJvLbTvdUo66m7l1eNvM5vI2tei82Dw7vKEWMVeJzsu5wr5q/dXi8WpRfKUYtNq+KhF41ai7Yn+VbTpEye2wcXZG1kmSNNRa5lg/LjOvWcP9S5IkSb3xUhBJkiSpB2s1sU7gIxFxbkQcP2yFiDg+Is6JiHM2X7dplcuTJEmSVmatLgV5VGZ+PyJ+AvhoRHw1Mz89uEJmngScBLDf/e7oRZ2SJEmaaGtyxjozv9/+exXwfuDha1GHJEmS1JdVP2MdEbsDM5l5Y3v/ScBra9vMLmzgh7cOTwXhhp2GLt75x1Hsb5fryyfAd71+ePLCzj8up3FsqCR/bLhlS7EtNg1P8YjN5W2YLe+rmOJRWg71pI5Swsc4kjpKCRodkzA6JXKsZuqGyRqSJO1w1uJSkAOB90fE4v7flZkfXoM6JEmSpN6s+sQ6M78NPHi19ytJkiSNk3F7kiRJUg+cWEuSJEk9cGItSZIk9cCJtSRJktSDtfqCmBWZm5/hyh8Pj9vb5ZoNQ5ff4epynNluV5fj4na+bnjU3cYbNxe3iVsqbZvK0Xm5ZXhbbhkewwcsE49XaOsSqQed4ueqMXd9x9kZWSdJkiaIZ6wlSZKkHjixliRJknrgxFqSJEnqgRNrSZIkqQdOrCVJkqQeOLGWJEmSejAVcXsLsxu49Yo9hrbtd/nwyLU9v1eOudvl6luKbTM3DG/LW24tF7i5HLe3MDtXbMu5Qlslsq4YqVfTNebOODtJkqSRecZakiRJ6oETa0mSJKkHTqwlSZKkHjixliRJknrgxFqSJEnqwVSkgmy4Bfb70vC/AQ4474bh21xxTbG/hRtvKrbNb5kduryaxlFL3TBZQ5IkaV3wjLUkSZLUAyfWkiRJUg+cWEuSJEk9cGItSZIk9cCJtSRJktQDJ9aSJElSDyKnIA4uIq4GvrvWdayhA4ByfuD64lg0HIetHIutHIutHIutHIutHIutHIutlo7FXTPzjl06moqJ9XoXEedk5pFrXcckcCwajsNWjsVWjsVWjsVWjsVWjsVWjsVWfY6Fl4JIkiRJPXBiLUmSJPXAifV0OGmtC5ggjkXDcdjKsdjKsdjKsdjKsdjKsdjKsdiqt7HwGmtJkiSpB56xliRJknrgxFqSJEnqgRPrNRYRb4+IqyLiooFl+0XERyPiG+2/+7bLIyL+NiK+GREXRsRD1q7y/kXEoRHxyYi4OCK+EhEvbpevu/GIiF0j4osR8aV2LP64XX63iPhCe8ynRcTO7fJd2sffbNsPW9MDGIOI2BAR50fEme3jdTkWEXFpRHw5Ii6IiHPaZevuZwQgIvaJiDMi4qsRcUlEPHI9jkVE3Kd9PSzeboiIl6zHsQCIiP/d/t68KCJObX+frrvfFxHx4nYMvhIRL2mXrZvXRPQ0v4qI49r1vxERxy23XyfWa+8dwJOXLHsl8PHMvBfw8fYxwFOAe7W344E3r1KNq2UOeFlm3h84CnhRRNyf9Tkem4HHZ+aDgcOBJ0fEUcD/Bd6QmfcErgOe367/fOC6dvkb2vV2NC8GLhl4vJ7H4nGZefhA7up6/BkBeCPw4cy8L/BgmtfHuhuLzPxa+3o4HHgocAvwftbhWETEwcDvAkdm5gOBDcCzWWe/LyLigcALgIfT/Gw8LSLuyfp6TbyD7ZxfRcR+wGuAR9CM5WsWJ+NFmeltjW/AYcBFA4+/BhzU3j8I+Fp7/63AscPW2xFvwAeA/7XexwPYDTiv/cG+BtjYLn8kcFZ7/yzgke39je16sda19zgGh7S/BB8PnAnEOh6LS4EDlixbdz8jwN7Ad5Y+t+txLJYc/5OA/1yvYwEcDFwG7Nf+/J8J/PR6+30B/ALwtoHH/wf4/fX2mmA751fAscBbB5Zvs96wm2esJ9OBmXlFe/9K4MD2/uIvjEWXt8t2OO3bcUcAX2Cdjkc0lz5cAFwFfBT4FnB9Zs61qwwe721j0bb/GNh/VQser7+h+U9hoX28P+t3LBL4SEScGxHHt8vW48/I3YCrgX+M5hKhf4iI3VmfYzHo2cCp7f11NxaZ+X3gL4HvAVfQ/Pyfy/r7fXER8FMRsX9E7AY8FTiUdfiaWGKlx7/icXFiPeGy+RNpXWUiRsQewHuBl2TmDYNt62k8MnM+m7d2D6F5C+q+a1vR2oiIpwFXZea5a13LhHhUZj6E5q3LF0XEowcb19HPyEbgIcCbM/MI4Ga2vq0LrKuxAKC9bvho4D1L29bLWLRv0z+d5g+vOwO7c/vLAXZ4mXkJzWUtHwE+DFwAzC9ZZ128JkrGdfxOrCfTDyPiIID236va5d+n+Ytz0SHtsh1GROxEM6l+Z2a+r128bscDIDOvBz5J8/blPhGxsW0aPN7bxqJt3xv40epWOjY/CRwdEZcC76a5HOSNrM+xWDwjR2ZeRXMd7cNZnz8jlwOXZ+YX2sdn0Ey01+NYLHoKcF5m/rB9vB7H4onAdzLz6sycBd5H8ztk3f2+yMy3ZeZDM/PRNNeVf531+ZoYtNLjX/G4OLGeTB8EFj95ehzNtcaLy3+1/fTqUcCPB97SmHoREcDbgEsy868HmtbdeETEHSNin/b+HWiuNb+EZoJ9TLva0rFYHKNjgE+0f41Pvcz8g8w8JDMPo3mb+xOZ+RzW4VhExO4RsefifZrraS9iHf6MZOaVwGURcZ920ROAi1mHYzHgWLZeBgLrcyy+BxwVEbu1/6csvi7W4++Ln2j/vQvw88C7WJ+viUErPf6zgCdFxL7tuyFPapeVrfWF5ev9RvNL8ApgluYMzPNpru/6OPAN4GPAfu26Afw9zbW2X6b51POaH0OPY/EomrdlLqR52+oCmuvC1t14AA8Czm/H4iLgj9rldwe+CHyT5u3eXdrlu7aPv9m2332tj2FM4/JY4Mz1OhbtMX+pvX0F+MN2+br7GWmP73DgnPbn5F+BfdfxWOxOc6Z174Fl63Us/hj4avu785+BXdbp74vP0PxR8SXgCevtNUFP8yvgee3r45vAry23X7/SXJIkSeqBl4JIkiRJPXBiLUmSJPXAibUkSZLUAyfWkiRJUg+cWEuSJEk9cGItSWMWEX8YEV+JiAsj4oKIeMSY93d2RBy5gvWPiogvtLVdEhEntMuPjohXLrO5JKm1cflVJEldRcQjgacBD8nMzRFxALDzGpe11CnAszLzSxGxAbgPQGZ+kOaLEyRJI/CMtSSN10HANZm5GSAzr8nMHwBExB9FxH9HxEURcVL7TXGLZ5zfEBHntGeQHxYR74uIb0TEn7brHBYRX42Id7brnBERuy3deUQ8KSI+FxHnRcR7ImKPITX+BM0XKZCZ85l5cbvtcyPixPb+BQO3WyPiMe03Qb49Ir4YEedHxNPHMH6SNDWcWEvSeH0EODQivh4Rb4qIxwy0nZiZD8vMBwJ3oDmzvWhLZh4JvIXma3dfBDwQeG5E7N+ucx/gTZl5P+AG4LcGd9yeHX818MTMfAjNtxS+dEiNbwC+FhHvj4gXRsSuS1fIzMMz83Dg/7T9/BfwhzRfAf1w4HHAX7RftS5J65ITa0kao8y8CXgocDxwNXBaRDy3bX5ce23zl4HHAw8Y2HTxEowvA1/JzCvas97fBg5t2y7LzP9s7/8L8Kgluz8KuD/wnxFxAXAccNchNb4WOJLmj4BfAj487Fgi4l7AX9BcNjILPAl4Zdv32TRfD32XynBI0g7Na6wlacwyc55m4nl2O4k+LiLeDbwJODIzL2s/MDh4pnhz++/CwP3Fx4u/u3PprpY8DuCjmXnsCDV+C3hzRJwMXD1wVrzpqLmE5HTgBZl5xUD/z8zMry3XvyStB56xlqQxioj7tGd6Fx0OfJetk+hr2knrMR26v0v74UhozjR/dkn754GfjIh7trXsHhH3HlLjzyxe3w3cC5gHrl+y2tuBf8zMzwwsOwv4nYFrw4/ocAyStMPwjLUkjdcewN9FxD7AHPBN4PjMvL49O3wRcCXw3x36/hrwooh4O3Ax8ObBxsy8ur3s5NSI2KVd/Grg60v6+RXgDRFxS1vjczJzfnGuHRF3pZn43zsintdu8+vAnwB/A1wYETPAd9j2OnFJWlcic+k7h5KkSRcRhwFnth98lCRNAC8FkSRJknrgGWtJkiSpB56xliRJknrgxFqSJEnqgRNrSZIkqQdOrCVJkqQeOLGWJEmSeuDEWpIkSeqBE2tJkiSpB06sJUmSpB44sZYkSZJ64MRakiRJ6oETa0mSJKkHTqwlSZKkHjixliRJknrgxFqSJEnqgRNrSZIkqQdOrCVJkqQeOLGWJEmSeuDEWpIkSeqBE2tJkiSpB06sJUmSpB44sZYkSZJ64MRakiRJ6oETa0mSJKkHTqwlSZKkHjixliRJknrgxFqSJEnqgRNrSZIkqQdOrCVJkqQeOLGWJEmSeuDEWpIkSeqBE2tJkiSpB06sJUmSpB5sXOsCRnFA3Cm3sGV4Y8TWu0Pbiw+qfW1fPyOuW2yKkXazzQrLrjtsvZUcx7ayy75XsOtRjy1X1GeXOgbXjdvvc7v7LDcV9zOG4xy6r3GP57B9r9I+h23T6/Pa07a5WuPRadvsvr+uNQ7dbplnbtRfp8vsJ5Z/hTTrLbujrf0s93/W8K5uX8dy+xxW+3L/FSy7zbDlsZJtRj2O+niN2k/p+du66jL7WfbYOtQ5ZPvSc7m4/ah9b7vu6hzbsq+PldR228/t8F7rfQ4uu/3SGPJo6VrnXrj5rMx88tCdL2MqJtZb2MIj4gkQW0+wx8ziqA9ZNrh8m2XtAM7M3G5Zs+7MkGXtNgP7Gdbn0H5K7Yv3h/WzTO1D+xlcd2BZbtM+bD/dtxm6Xqk9hrQP3eb2ywa3H7bP0n5u62ubY+P2684MaR/WT7HO8nrFPrc59pXtZ9k6i7UvLltuvJarbeV1LLsNQ9btcGzbXUeXPlmmvcMYDk5OuvS5tfYh/TBkvW32k0Pbh43HbesWx/D27THqfhjeHkP7rOxnYF+Dk4etv/rq+9n2xzJv377Nr+Bh+7n9PmeGtM9w+/VKfQ7dfrCd+jYzw7YZUseo2wxut237QnG9bftcuG3ZhqF1LNxumw3LtW+zz9vXsYFh2wwsG+h/w23HMVAnQ2ofUtPQ/TB4vLff57bLbj822/Z5+/0MO/YN29S5UDzGbY6tUOeGIc/1Yp8bCq+fxeXbtnP7ZWy14bb2KLRHu2xr++K62yyLwe1nhrTPbNN22/KDvnEAHXkpiCRJktQDJ9aSJElSD5xYS5IkST1wYi1JkiT1wIm1JEmS1AMn1pIkSVIPnFhLkiRJPXBiLUmSJPXAibUkSZLUAyfWkiRJUg+cWEuSJEk9cGItSZIk9SAyc61rWFZEXARsWus6ptgBwDVrXcSUcuy2j+O3fRy/7hy77eP4bR/Hr7tJGLtrMvPJXTbc2HclY7IpM49c6yKmVUSc4/h149htH8dv+zh+3Tl228fx2z6OX3fTPnZeCiJJkiT1wIm1JEmS1INpmViftNYFTDnHrzvHbvs4ftvH8evOsds+jt/2cfy6m+qxm4oPL0qSJEmTblrOWEuSJEkTzYm1JEmS1IOJmlhHxJMj4msR8c2IeOWQ9l0i4rS2/QsRcdgalDmRRhi7R0fEeRExFxHHrEWNk2yE8XtpRFwcERdGxMcj4q5rUeekGmH8fiMivhwRF0TEZyPi/mtR5yRabuwG1ntmRGRETG0M1TiM8Np7bkRc3b72LoiIX1+LOifVKK+/iHhW+/vvKxHxrtWucVKN8Np7w8Dr7usRcf0alDmxRhi/u0TEJyPi/Pb/3qeuRZ0rlpkTcQM2AN8C7g7sDHwJuP+SdX4LeEt7/9nAaWtd9yTcRhy7w4AHAf8EHLPWNU/SbcTxexywW3v/N33trXj89hq4fzTw4bWuexJuo4xdu96ewKeBzwNHrnXdk3Ib8bX3XODEta51Em8jjt+9gPOBfdvHP7HWdU/CbdSf3YH1fwd4+1rXPSm3EV97JwG/2d6/P3DpWtc9ym2Szlg/HPhmZn47M7cA7waevmSdpwOntPfPAJ4QEbGKNU6qZccuMy/NzAuBhbUocMKNMn6fzMxb2oefBw5Z5Ron2Sjjd8PAw90BPzXdGOX3HsCfAP8Xv4F2qVHHT8ONMn4vAP4+M68DyMyrVrnGSbXS196xwKmrUtl0GGX8Etirvb838INVrK+zSZpYHwxcNvD48nbZ0HUycw74MbD/qlQ32UYZO5WtdPyeD/zHWCuaLiONX0S8KCK+Bfw/4HdXqbZJt+zYRcRDgEMz899Xs7ApMerP7jPbt5LPiIhDV6e0qTDK+N0buHdE/GdEfD4iOn3N8w5o5P832ksH7wZ8YhXqmhajjN8JwC9HxOXAh2jO+k+8SZpYSxMvIn4ZOBL4i7WuZdpk5t9n5j2AVwCvXut6pkFEzAB/DbxsrWuZYv8GHJaZDwI+ytZ3PTWajTSXgzyW5qzryRGxz1oWNIWeDZyRmfNrXciUORZ4R2YeAjwV+Of2d+JEm6QCvw8Mnkk4pF02dJ2I2Ejz1sCPVqW6yTbK2KlspPGLiCcCfwgcnZmbV6m2abDS19+7gWeMs6ApstzY7Qk8EDg7Ii4FjgI+6AcYb7Psay8zfzTw8/oPwENXqbZpMMrP7uXABzNzNjO/A3ydZqK93q3k996z8TKQpUYZv+cDpwNk5ueAXYEDVqW67TBJE+v/Bu4VEXeLiJ1pXogfXLLOB4Hj2vvHAJ/I9qr2dW6UsVPZsuMXEUcAb6WZVHuN4bZGGb/B/4h/BvjGKtY3yapjl5k/zswDMvOwzDyM5vr+ozPznLUpd+KM8to7aODh0cAlq1jfpBvl/45/pTlbTUQcQHNpyLdXscZJNdL/uxFxX2Bf4HOrXN+kG2X8vgc8ASAi7kczsb56VavsYGIm1u01078NnEXzi+/0zPxKRLw2Io5uV3sbsH9EfBN4KVCMplpPRhm7iHhYe53SLwBvjYivrF3Fk2XE195fAHsA72mjk/zDpTXi+P12G9V1Ac3P7nHDe1tfRhw7FYw4fr/bvva+RHNt/3PXptrJM+L4nQX8KCIuBj4JvDwz1/07xSv42X028G5PAm5rxPF7GfCC9mf3VOC50zCOfqW5JEmS1IOJOWMtSZIkTTMn1pIkSVIPnFhLkiRJPXBiLUmSJPXAibUkSZLUAyfWkrQORMRzI+LEta5DknZkTqwlSZKkHjixlqQJExGHRcRXI+KdEXFJRJwREbsNtM9ExKURsc/Asm9ExIER8bMR8YWIOD8iPhYRBw7p/x0RcczA45sG7r88Iv47Ii6MiD8e42FK0g7HibUkTab7AG/KzPsBNwC/tdiQmQvAB4CfA4iIRwDfzcwfAp8FjsrMI4B3A78/6g4j4knAvYCHA4cDD42IR/dyNJK0DjixlqTJdFlm/md7/1+ARy1pPw34xfb+s9vHAIcAZ0XEl4GXAw9YwT6f1N7OB84D7ksz0ZYkjcCJtSRNplzyeO+IuKC9HQ18DrhnRNwReAbwvna9vwNOzMz/AbwQ2HVI33O0v/8jYgbYuV0ewOsy8/D2ds/MfFuvRyVJOzAn1pI0me4SEY9s7/8ScObAhPeDmZnA+4G/Bi7JzB+16+4NfL+9f1yh70uBh7b3jwZ2au+fBTwvIvYAiIiDI+InejsiSdrBObGWpMn0NeBFEXEJsC/w5iHrnAb8MlsvAwE4AXhPRJwLXFPo+2TgMRHxJeCRwM0AmfkR4F3A59pLSc4A9tz+Q5Gk9SGakx6SpEkREYfRnKF+4FrXIkkanWesJUmSpB54xlqSJEnqgWesJUmSpB44sZYkSZJ64MRakiRJ6oETa0mSJKkHTqwlSZKkHvx/mK8zIGYHjUEAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 864x432 with 2 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"xlab_idx = np.arange(9,len(N),10)\n",
"ylab_idx = np.arange(4,len(S),5)\n",
"plt.figure(figsize=(12,6))\n",
"plt.title('t-test p-values and sample size (deterministic)')\n",
"plt.imshow(sim, origin='lower')\n",
"plt.xticks(xlab_idx, N[xlab_idx])\n",
"plt.yticks(ylab_idx, S[ylab_idx])\n",
"plt.xlabel('Sample Size')\n",
"plt.ylabel('Standard Deviation')\n",
"plt.colorbar(aspect=50, orientation='horizontal', pad=0.2, label='p-value')\n",
"plt.savefig('./img/pvalues_big_data_deterministic.svg', dpi=500)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Chi-squared test\n",
"\n",
"In this section a similar relationship between p-values and large sample sizes are shown using the Chi-squared test. The test statistic (with Yate's correction for continuity) is computed as follows,\n",
"\n",
"$$\n",
"\\text{stat} = \\sum_{i=1}^{N} \\frac{(|O_i - E_i| - 0.5)^2}{E_i}\n",
"$$\n",
"\n",
"where $O_i$ is the observed frequency, $E_i$ is the expected frequency under the null hypothesis and $N$ is the number of events. $E_i$ is calculated as the product of the marginal sums divided by $N$. The statistic is chi-squared distributed so the p-value for a 2x2 table can be calculated as follows using the chi-squared CDF $F$.\n",
"\n",
"$$\n",
"\\text{pvalue} = 1-F(\\text{stat}, 1)\n",
"$$\n",
"\n",
"`scipy.stats.chi2_contingency` can be used to compute the test statistic and associated p-value. A manual Python implementation of this calculation is provided in the Appendix."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"prob_a = 0.5\n",
"prob_b = 0.55\n",
"sim = []\n",
"for n in range(10,5_001,10):\n",
" rep = []\n",
" for i in range(0,100):\n",
" a = np.random.binomial(1, prob_a, n)\n",
" b = np.random.binomial(1, prob_b, n)\n",
" a_sum = a.sum()\n",
" b_sum = b.sum()\n",
" table = np.array([[a_sum, n-a_sum],[b_sum, n-b_sum]])\n",
" test = stats.chi2_contingency(table)\n",
" rep.append(test[1])\n",
" sim.append(np.mean(rep))"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAEWCAYAAACaBstRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABAMElEQVR4nO3dd3gc1dn38e9P3ZLcey9gQ6gGjDGhhB7TQp7QSYFAwpOHkISSAikOgRTSIISWkITQa14SHDC9d2xjigvYxr3bsiVZVpfu94+ZXY/WK3lt76ren+vaSzszZ2bO7K723lPmHJkZzjnnXCqy2joDzjnnOg4PGs4551LmQcM551zKPGg455xLmQcN55xzKfOg4ZxzLmUeNNJM0jWS7mth+xxJR7VejjJP0l2SftnGeRglySTltGU+tmd7n4/2IHwdd0/zMSskjUnnMcPjviHpgHQft4XztcpnXdKXJT2bQrrvSPptpvMT5UFjJ0g6T9KM8B9htaSnJB2eyr5mtreZvZzhLDrXrphZsZktSucxJZ0KbDazWdtJ165/UCTLn5ndb2YnpLD734AvSxqQuRw25UFjB0m6AvgT8GtgIDACuA04rQ2zlTbt9R/LuSS+Bdzb1ploS2ZWDTwFfK21zulBYwdI6glcC3zbzB4zsy1mVmdm/zWzH0SS5km6R9LmsDpqQuQYSyQd18zxCyTdJ6lEUqmk6ZIGhttGS3olPOZzkm6JVXNIOkrSioRjxc8jaaKkt8Jjrg73zYukNUnflrQAWBCuO0XS++E+b0raL5L+AEnvhXl5GCho4TW7IKxCuEVSmaSPJR3bTNqzJc1IWHe5pKnh85MlzZJULmm5pGtaOG+T1zmxWkjSpPC6SiV90FKVoaSrJH0aXu9cSf+TcH2vS/qDpE2SFks6MbK9yfsG9GvhPP0kPRHmaaOk1yRlpZiHNyTdGO67SNJnw/XLJa2TdH4k/V2S/hJ+jjaH+RvZTJ7yw2tbJmltuF+3ZtLuHh6rTNKG8LMR22bh9iEKSuixR6Uki6S7UNK88LV8poV85QHHAK9E1k1UUANQHub1hnDTq+Hf0vCch0rKkvRTSUvD1+ceBf/fsWMdHvl8LJd0QeT0vSU9Gb5270jaLbLfTWH6ckkzJR2xk/m7QNLrkX33Dt+vjeG+P47k52Xg5GSvU0aYmT9SfACTgXogp4U01wDVwElANvAb4O3I9iXAcc3s+7/Af4HCcN+DgB7htreAG4B84EhgM3BfuO0oYEXCseLnCY8zCcgBRgHzgMsiaQ14DugDdAMOANYBh4T5OD88Xj6QBywFLgdygTOAOuCXzVzTBeFrFkt/NlAG9EmStjC8rrGRddOBcyLXuS/Bj539gLXAF8Nto8LryEn2OofvS+z1GgqUhO9RFnB8uNy/mWs4ExgSpj0b2AIMjlxfHfDN8LX6P2AVoO29b0nO8xvgL+HrlAscETnO9vJQD3w9zMMvgWXAreF5TwjPWxymvytcPjLcfhPwesLnYffw+Y3A1PCz0Z3g8/mbZvL/IPCTMI8FwOHJjpmwz/3Ag+Hz04CFwGcIPqs/Bd5s5lx7A1sS1r0FfDV8XgxMSvbZCNddGJ5rTJj2MeDecNvI8PU5N3wf+gLjI69dCTAxzOP9wEOR434lTJ8DXAmsAQp2In8XxN6T8HVfHR6vIFw+JJL2QGBjq30PttaJOsMD+DKwZjtprgGejyzvBVRFlpfQfNC4EHgT2C9h/QiCL4WiyLoHSDFoJDnPZcC/I8sGHBNZvh24LmGfT4DPEXzRxL8Uw21v0nLQSEz/buyfJ0n6+4Ap4fOx4T9vYTNp/wTcGD5v8o+XeP00DRo/IvyCiGx/Bjg/xc/B+8BpketbGNlWGOZj0PbetyTHvRZ4nCRfrinkYUFk275hHgZG1pXQ9Isv+kVXDDQAwyOfh90BEQSn3SJpDwUWN5One4A7gGFJtm0TNML3YSbQLVx+Crgosj0LqARGJjneYST8LxL8Yv8F0C9hfZPPRrjuBeCSyPIeBME/B7iayP9HwrHuAv4eWT4J+LiF92kTsP9O5O8CtgaNc4FZLZxjLNCQymc3HQ+vntoxJUA/bb/ef03keSVQkGyfhGL6CIL62WeAhyStkvQ7SbkEvzA3mdmWyO5LU820pHFhtccaSeUE7TGJ1STLI89HAleGRfNSSaXA8DAfQ4CVFn5aU8xLsvRDJB0Ruf454bYHCP5JAM4D/mNmleF1HCLpJUnrJZUR1Gk3W93TgpHAmQnXdzgwOFliSV/T1qq6UmCfhPPG3+9YXgm+iHf0ffs9wa/fZ8Mqpqt2IA9rI8+rwrwkriuOLMffbzOrADaG+Y3qTxAEZ0bO+3S4PpkfEgSadxVUy17Y3IUqqML7HkFJsSpcPRK4KXKujeHxhiY5xCaCX9xRFwHjgI8VVO2e0tz5Ca41+l4sJQgYAwk+65+2sG/i/3f8dZX0/bB6rSy8hp5sfZ92JH9R28tPd4LSe6vwoLFj3gJqgC+m42AW9CiJPZZZ0D7yCzPbC/gscApBA9dqgnrUosjuIyLPtxD8cwMgKZum/9i3Ax8TVPv0AH5M8M/YJDuR58uBX5lZr8ij0MweDPMyVFJ0/xG0LFn6VWb2WuT69w63PQf0lzSeIHg8ENnvAYKqkuFm1pOgKifxOmKavCYEv/yj13dvwvUVmdn1iQcJ69T/BlwK9DWzXsDsFs4btb33rQkz22xmV5rZGOALwBWSjt3FPDRneOyJpGKC6qdVCWk2EASbvSOvU08zKyYJM1tjZt80syEEVa23KUnXXUl7AHcDZ5lZ9MfKcuB/E96Xbmb2ZpLTLQwOpXhAMbMFZnYuMAD4LfCv8LW3JPuvIghSMbFS4dowH7sl2adFYfvFD4GzgN7h+1RG+D7tYP6ilhNUozXnM8AHO5rfneVBYweYWRkwBbhV0hclFUrKlXSipN/t6vElHS1p3/BLv5yguNxoZkuBGcAvJOUp6N57amTX+QSlmZPDkslPCeqqY7qHx6uQtCdBvXtL/gZ8K/xlL0lF4bG7EwTOeuC74bV/iaB+tyUDIunPJPiQT0uW0MzqgEcJfnX3IQgi0evYaGbVkiYSlESa8z5wTnjOCQRtLzH3AadK+rykbAUdEI6SNCzJcWL/1OsBJH2d4Ff+dqXwvjWhoPPB7mGALSOoMmrclTy04CQFjb15wHUE7W7RL3DMrJHgs3Cjwi6dkoZK+nwz+T8z8hpuCvPcmJCmB0EV3E/M7PWEQ/wFuFrS3mHanuHnZRtmVgs8T1BlGjv2VyT1D/NdGq5uJHjdGmn6xfsgcLmCjgrFBKXvh82snqCd4jhJZ0nKkdQ3/BGzPd0J/jfWAzmSpgA9djJ/UU8AgyVdpqBjQndJh0S2f46gaq9VeNDYQWb2R+AKgi/m9QS/Ai4F/pOGww8C/kXwBT+PoGdIrEvheQQN0xuBnxPUH8fyVAZcAvwdWEnwKzvam+r74f6bCb4EHqYFZjaDoGH3FoJ//oUEdayxf9YvhcsbCRplH9vOdb1DUO+6AfgVcIaZlbSQ/gHgOODR8J845hLgWkmbCYL3Iy0c42cEvxY3EdQjx0ss4ZfjaQQlrth7+AOS/D+Y2VzgjwTBci1Be8EbLZw3UbPvWxJjCb4IK8Lz3WZmL6UhD8k8EOZnI0FHia80k+5HBO//2wqqNp8nqP9P5mDgHUkVBCXC79m292YcGO5/oyLVswBm9m+CX+APheeaDZxI8/4KfDWyPBmYEx7vJoIOFFVhleGvgDfCqq9JwJ0E/1uvAosJOq98J8zHMoK2iivD1+d9YP8W8hHzDEH13XyC6q5qmlb77kj+4sxsM0FnjVMJqsYWAEdD0OMyzOvdKeQvLWI9M1wHo6C76e5m1tw/e7ugoKviN8wspZsfXeZJuoug48RP2zovu0rSG8Cltp0b/DorSd8hqK79YWud02/kcs51WGZ2WFvnoS2Z2c2tfU6vnnLOOZcyr55yzjmXMi9pOOecS1mnb9Po16+fjRo1qq2z4ZxzHcbMmTM3mFnSmzg7fdAYNWoUM2bM2H5C55xzAEhqduQCr55yzjmXMg8azjnnUuZBwznnXMo8aDjnnEuZBw3nnHMp86DhnHMuZR40nHPOpcyDRjPufH0x0z5a3dbZcM65dsWDRjPuf2cpT37oQcM556I8aDSjb3E+Gypq2jobzjnXrnjQaEa/4jxKttS2dTacc65d8aDRjL5F+ZR4ScM555rwoNGMvsV5bKqso76hsa2z4pxz7YYHjWb0Lc4HYGOlV1E551yMB41m9CvKA6CkwoOGc87FeNBoRu8waGz0xnDnnIvzoNGM4vxgfqqKmvo2zolzzrUfHjSaURQGjS0eNJxzLs6DRjOK8rMBDxrOORflQaMZW6unGto4J84513540GhGt9xssuQlDeeci/Kg0QxJFOXleEO4c85FeNBoQVF+jpc0nHMuot0EDUmTJX0iaaGkq5pJc5akuZLmSHog03kqys9mS60HDeeci8lp6wwASMoGbgWOB1YA0yVNNbO5kTRjgauBw8xsk6QBmc5XcX6ON4Q751xEeylpTAQWmtkiM6sFHgJOS0jzTeBWM9sEYGbrMp2povwcKr16yjnn4tpL0BgKLI8srwjXRY0Dxkl6Q9LbkiY3dzBJF0uaIWnG+vXrdzpTRfneEO6cc1HtJWikIgcYCxwFnAv8TVKvZAnN7A4zm2BmE/r377/TJyzOz/E2Deeci2gvQWMlMDyyPCxcF7UCmGpmdWa2GJhPEEQypig/my3epuGcc3HtJWhMB8ZKGi0pDzgHmJqQ5j8EpQwk9SOorlqUyUx59ZRzzjXVLoKGmdUDlwLPAPOAR8xsjqRrJX0hTPYMUCJpLvAS8AMzK8lkvorzcqitb6TOZ+9zzjmgnXS5BTCzacC0hHVTIs8NuCJ8tIroSLe9CvNa67TOOddutYuSRnvlc2o451xTHjRaUBgOj374b1/ywOGcc3jQaFGsegqgtNKnfXXOOQ8aLSiOBI26BmvDnDjnXPvgQaMFRXlbg0Z1nd+v4ZxzHjRaEC1peNBwzjkPGi2KzRMOUF3n92o455wHjRYUF0RKGvVe0nDOOQ8aLcjPyeb+bxwCQI1XTznnnAeN7RncswDw6innnIM0DyMiqQA4BTgCGAJUAbOBJ81sTjrP1VoKcoN2DW8Id865NAYNSb8gCBgvA+8A64ACgtForw8DypVm9mG6ztkaPGg459xW6SxpvGtmP29m2w3hnN4j0ni+VtEtFjTqvXrKOefSFjTM7MnEdZKygGIzKw/n9M74vN7plp8TNPt4ScM55zLQEC7pAUk9JBURtGfMlfSDdJ+ntWRlCYA/Pb+AdZur2zg3zjnXtjLRe2ovMysHvgg8BYwGvpqB87S6G59b0NZZcM65NpWJoJErKZcgaEw1szqgU4z2t8WHR3fOdXGZCBp/BZYARcCrkkYC5Rk4T6vbUFHT1llwzrk2lfagYWZ/NrOhZnaSBZYCR6f7PG1h0fotAKwqrSKYfdY557qWTDSE50s6T9KPJU2RNAX4cbrP05p+ctJn6F2YS1lVHfPXbuaz17/InW8saetsOedcq8tE9dTjwGlAPbAl8uiwvnnkGC747Giq6hpYsLYCgLc+LWnjXDnnXOtL6zAioWFmNnlHd5I0GbgJyAb+bmbXJ2y/APg9sDJcdYuZ/X0X85qy2DDpm6vrAMj2Ubucc11QJoLGm5L2NbOPUt1BUjZwK3A8sAKYLmmqmc1NSPqwmV2axrymLDYhU3k8aKgtsuGcc20qE0HjcOACSYuBGkCAmdl+LewzEVhoZosAJD1EUMWVGDTaTGxujfKqoNttljxoOOe6nkwEjRN3Yp+hwPLI8grgkCTpTpd0JDAfuNzMlidJg6SLgYsBRoxIz3BXRWFJY1NlLQA5XtJwznVBmehyuxToBZwaPnqF63bVf4FRYYnlOeDuFvJwh5lNMLMJ/fv3T8Opt1ZPrd8c3KuRlSXvduuc63Iy0eX2e8D9wIDwcZ+k72xnt5XA8MjyMLY2eANgZiVmFru77u/AQenJcWqK8oKg8ezctQA89t5Krnz0g9bMgnPOtblM9AG6CDjEzKaY2RRgEvDN7ewzHRgrabSkPOAcYGo0gaTBkcUvAPPSmOft6l6wbU3eY++tTJLSOec6r0y0aQiIjiPeEK5rlpnVS7oUeIagy+2dZjZH0rXADDObCnxX0hcI7v/YCFyQgbw3K9am4ZxzXVkmvgn/Cbwj6d/h8heBf2xvJzObBkxLWDcl8vxq4Or0ZXPHxO7TcM65riztQcPMbpD0MkHXW4Cvm9msdJ+nteXneNBwzrm0tWlI6hH+7UMwyu194WNpuK7Dm//LbXsTNzR6DyrnXNeRzpLGA8ApwEyazp+hcHlMGs/VJvJyto2xmypr6Vec3wa5cc651pfOOcJPCf+OTtcx26PCvGwqa7e285dUeNBwznUdmbhP44VU1nVU7/3s+CbLNfUNzaR0zrnOJ20lDUkFQCHQT1Jvtnaz7UEwTEinUJDbtEG8rqGxjXLinHOtL51tGv8LXAYMIWjXiAWNcuCWNJ6nXamt94Zw51zXkc42jZuAmyR9x8xuTtdx2zsvaTjnupJM3Kdxs6R9gL2Agsj6e9J9rraSn5NFTX0QLOobPWg457qOTDSE/xy4OXwcDfyOYKyoTuPJ7x7BRYcHncS8eso515VkYsDCM4BjgTVm9nVgf6BnBs7TZnYfUMy5E4NBeb16yjnXlWQiaFSZWSNQH94lvo6mw553CrnhJOF1DY3MXLqJ5Rsr2zhHzjmXeZkYsHCGpF7A3wh6UVUAb2XgPG0qGjROv/1NsrPEp78+qY1z5ZxzmZWJhvBLwqd/kfQ00MPMPkz3edpaLGhU1AQ39/kYVM65riATDeH/I6kngJktAZZJ+mK6z9PW8sKg4dVSzrmuJBNtGj83s7LYgpmVAj/PwHnaVE52cO/i0pItbZwT55xrPZkIGsmO2emmvYtVTy31koZzrgvJRNCYIekGSbuFjxsIGsQ7ldywpBGrnirK80manHOdXyaCxneAWuDh8FEDfDsD52lTksjNFnUNQQN4ZV0Do656kr+/tqiNc+acc5mTid5TW4Cr0n3c9ig3O4u6hqD3lIWdp349bR7fOKLDzzflnHNJpXNo9D+Z2WWS/kvTmfsAMLNONZQIxNo1ms6n0Wgwc+kmDhzRC0nJd3TOuQ4qnSWNe8O/f9jZA0iaDNwEZAN/N7Prm0l3OvAv4GAzm7Gz59tVscbwRKff/iZ/PHN/Tj9oWCvnyDnnMiudQ6PPDP++sjP7S8oGbgWOB1YA0yVNNbO5Cem6A98D3tm1HO+6vLAxPCdL1Cfc3Le+oqYtsuSccxmVzuqpj0hSLUUwGZOZ2X7bOcREYKGZLQqP9xBwGjA3Id11wG+BH+xajnddbk5Q0hjWuxtLSpp2ve3VLZctNfUU5Xe63sbOuS4snd9op+zi/kOB5ZHlFcAh0QSSDgSGm9mTkpoNGpIuBi4GGDFixC5mq3k5WUFJY1jvwm2Cxu2vfMpVj33EOz8+loE9CpLt7pxzHU7autya2dLYg6Cb7f7AfkBNuG6XSMoCbgCuTCEvd5jZBDOb0L9//109dbNi3W2H9+m2zbalYRBZWVqVsfM751xry8TYU98A3gW+RDC3xtuSLkxh15U0HUJ9WLgupjuwD/CypCXAJGCqpAnpyPfO2LilFoDR/Yri6+6+cGKTNGY+kKFzrvPIRIX7D4ADzKwEQFJf4E3gzu3sNx0YK2k0QbA4BzgvtjEcz6pfbFnSy8D327L3VEVNPQDjh/eOr9tjYHd6FORQXh1s2xz+dc65ziATd4SXAJsjy5vDdS0ys3rgUuAZYB7wiJnNkXStpHZ9j8feQ3rEn3fLzW7S+B0LLM451xmkXNKQdDgw1sz+Kak/UGxmi5MkXQi8I+lxgt5UpwEfSroCwMxuaO4cZjYNmJawbkozaY9KNe+ZFg0SBXlZFEbGofKShnOuM0kpaEj6OTAB2AP4J5AL3AccliT5p+Ej5vHwb/edz2b7dMt5B1BWVddkXV52VpMgsrm6LnE355zrsFItafwPcADwHoCZrQpvstuGmf0CQFKhmXXqccNP2W/INuskeUnDOddppdqmUWtBNyADkFTUXEJJh0qaC3wcLu8v6bZdzmkHUpQXLWl40HDOdR6pBo1HJP0V6CXpm8DzwN+aSfsn4POEjd9m9gFw5C7ms0MpjFRPlXv1lHOuE0mpesrM/iDpeKCcoF1jipk910L65QkjvDY0l7azOG38ED5aEcxyO7B7fny9lzScc51Jyr2nwiDRbKCIWC7ps4BJyiUYXHDeTuavw7jpnAPiz38weQ9O2HsQ1z4xh/IqL2k45zqPlKqnJG2WVB4+qiU1SCpvJvm3CGbqG0pwk954OuHMfS3Jz8lm4ug+DOnZjdJKDxrOuc4j1eqpeE8pBfVOpxEM49FEOLz5TWb25bTlsAPrXZjH+8tL2zobzjmXNjt8R7gF/kPQ2J24rQEYKSkvDXnr8HoV5VJaWefjTznnOo1Ub+77UmQxi+BGv+pmki8C3pA0FdgSW9nSneCdVe/CPGobGqmsbfB5NZxznUKq32SnRp7XA0sIqqiSid0RnkUnvAt8R/QuzAVgU2UtW2rrGdDd59VwznVsqbZpfD3VA8buCHfQqzCopTv8ty8B8MxlR7LHoC4dR51zHVyLQUPSzSSfwhUAM/tu2nPUifTslttkecG6zR40nHMd2vZKGm02V0VnEJ2cCaCh0RvEnXMdW4tBw8zubq2MdEYDexRw7Wl7M+XxOQDU1DW2cY6cc27XpNp7qj/wI2AvIN6aa2bHRNLkABcRjIgbG/51JcHQ6P8wsy55l1vvwq29jzdV1rZhTpxzbtel2nvqfuBh4GSCO77PB9YnpLkXKAWuAVaE64aFae8Dzt61rHZM0XaNTZV1lFTUxGfzG9m32cGCnXOuXUo1aPQ1s39I+p6ZvQK8Iml6QpqDzGxcwroVwNuS5u9yTjuoHpGgUVpZy0G/fD6+/PRlR7DnoB7JdnPOuXYp1TvCY1VLqyWdLOkAoE9Cmo2SzpQUP6akLElnA5vSkNcOKVrSWFve9H7IyX96jaraTj8AsHOuE0m1pPFLST2BK4GbgR7A5QlpzgF+C9wmaRMgoBfwYritS4oGjdmrth3jsaKmnm6Rmf6cc649SzVovGNmZUAZcHSyBGa2hLDdQlLfcF1JGvLYoXUv2PoSr99cs812L2k45zqSVIPGG5KWEDSGP2ZmSaubJO1JMLzI0HB5JfC4mX2cykkkTQZuArKBv5vZ9QnbY8OuNwAVwMVmNjfFa2gTudkt1wBuqfVJmpxzHUdKbRphA/dPgb2BmZKekPSVaBpJPwIeIqiWejd8CHhI0lXbO0c4rPqtwIkEXXvPlbRXQrIHzGxfMxsP/A7oEIMgfjDlBGb97Pik2yq9pOGc60BSHhrdzN41syuAicBGIPHGv4uAg83sejO7L3xcH6a/KIVTTAQWmtkiM6slCEBNBkU0s2ijQBEtDHHSnvQszKV3UR7/+tah22yr9JKGc64DSXXmvh6Szpf0FPAmsJrgSz6qka039UUNDrdtz1BgeWR5RbguMS/flvQpQUkj6dhXki6WNEPSjPXrE28naTsTRvXhy4eMAOCsCcMA2FBRk7Stwznn2qNU2zQ+AP4DXGtmbzWT5jLgBUkL2PrlPwLYHbh0F/LYhJndCtwq6TyCKrPzk6S5A7gDYMKECe2qNCIFfwvzgpf+8oc/AODTX59EdpbaKlvOOZeSVIPGGAunn5N0ipk9kZjAzJ6WNI6gBBIrIawEpocz+m3PSmB4ZHlYuK45DwG3p5L59qgov2k32zc/3cARY/u3UW6ccy41qTaER3+tX9tCukYze9vM/l/4eNvMGiQVp3Ca6cBYSaPD6WLPAaZGE0gaG1k8GViQSv7bk/ycIFj0KGg6bPr0JV32/kfnXAeyM3OQ7kwdylyCqqpmmVm9pEuBZwi63N5pZnMkXQvMMLOpwKWSjiO4Q30TSaqm2rvvHjsWM/jqoSP5zVNbeyJv2uKDGTrn2r9UR7ktAC4BDgc2SbocuN3MqiNprmhudyCVkgZmNg2YlrBuSuT591I5TnvWs1suU04NehLn52RRUx/0Ebj37aX0LsrjiuMTh+9yzrn2I9Uut/cQ3KNxM3AdwX0U9yak+TXQm2Be8OijeAfO06UUJgwf8ucXOlxtm3Oui0m1emofM4veaPeSpMQ7sd8D/mNmMxN3lvSNnc2gc8659iPVEsB7kibFFiQdwrZTwX4dWNrM/hN2Im+dXo+EOcQB6hp8dj/nXPuVatA4CHhT0pJwDKq3gIMlfSTpQwAz+8TMNiTb2czWpiW3ncxeg7edSyNx+HTnnGtPUq2empzRXHRR3zhiNE/NXkNOlqhvDHo1ry6rZljvwjbOmXPOJZdS0DCz5qqd3C44aGQfFv/mJOat3sylD77HovVbKKnwIUWcc+1XWns1ScoOu+O6FEliryE9uPP8gwGoqvNRb51z7Vdag0Y4XMi56TxmVxGbva+q1hvCnXPtVybun3hD0i2SjpB0YOyRgfN0KgW5QdDYuKWG/8xaSdORW2DR+go2hFVXL328jveW+bAjzrnWtzPDiGzP+PBvdIwqA47JwLk6jW5h0PjDs/MBGNyzgEPG9I1vP+aPr1CUl82cayfz9bumA7Dk+pNbP6POuS4t7UHDzJLOIe5alpeT1aQX1dKNlfQuymPcwO7xNFtqG/w+Dudcm0p79ZSkgZL+EU7YhKS9JKUyc1+XFyttAPzwXx9ywo2vbpNm5aaq1sySc841kYk2jbsIRqqNzeI3n2CCJrcdBQljUQGYGQ2NW9s3lm2sbM0sOedcE5kIGv3M7BHCKV7NrB7wfqQpiJY0YmobGpvMI758UxA0ehVuOwSJc85lWiYawrdI6kvQ+E04ZlVZBs7T6SQLGis3VfHRyq0v3+bqIIAU5Gyb1jnnMi0TQeNKghn3dpP0BtAfOCMD5+l0klVPffuBWcxbXR5frg5v/pNPJ+6cawOZ6D01U9LngD0IJmD6xMzq0n2ezqhb7ra1hdGAAcQnbaqt915UzrnWl4neUx8CPwSqzWy2B4zUFeZtP4bHSho1HjScc20gEw3hpwL1wCOSpkv6vqQW5wd3gWRtGomq64JgUVPvfQucc60v7UHDzJaa2e/M7CDgPGA/YHG6z9MZFaQQNGrCkkZdQ9OuuM451xoy0RCOpJHA2eGjgaC6ym1H4pzhyUSrpWrrG+MDHTrnXGvIRJvGO8C/gWzgTDObaGZ/TGG/yZI+kbRQ0lVJtl8haa6kDyW9EAamTiWVey+qI0OnexWVc661ZaKk8TUz+2RHdpCUDdwKHA+sAKZLmmpmcyPJZgETzKxS0v8BvyMoyXQavQvztpsmOt9GrH3DOedaSyYawtdIukHSjPDxR0k9t7PPRGChmS0ys1rgIeC0aAIze8nMYmNovA0MS3/W21afouRB4/i9Bsafl1dv7YzmJQ3nXGvLRNC4E9gMnBU+yoF/bmefocDyyPKKcF1zLgKeam6jpItjQWv9+vUpZbo96B0GjSE9C/jaoVtr324+9wC+MinogDZv9eb4eu9265xrbZkIGruZ2c/DUsMiM/sFMCZdB5f0FWAC8Pvm0pjZHWY2wcwm9O/fP12nzrg+YfVUZV0D5392VHx9QW42ewzqAUBDo5GXE7xtNV495ZxrZZkIGlWSDo8tSDoM2N543iuB4ZHlYeG6JiQdB/wE+IKZ1aQhr+1K76KgIbyypoG87KZvTX5kuUdBkC5WPbW2vJp3FpW0Ui6dc11ZJoLGt4BbJS2RtAS4Bfjf7ewzHRgrabSkPOAcgvGr4iQdAPyVIGCsS3+2216sTaO2oZH8nKZvTV5kOSscd+rrd02nuq6BE258lbPveLvV8umc67oyMfbUB8D+knqEy+Xb2QUzq5d0KcE8HNnAnWY2R9K1wAwzm0pQHVUMPKpgtL5lZvaFdOe/LcXuCD/n4OFNggQ0HaCwe0EO6zbXsLm6nj1/9nR8/c8fn01NfSPXn74fa8urqaptYFS/olbJu3Oua8jIzX2QWrBISD8NmJawbkrk+XFpylq7JYmPr5tMXnYW1Qk9o0ort/aaKi7I5ZUfHMXnfv9ykzTPzFkbL60c8usXAJ9H3DmXXpmonnK7oCA3m6wskZvQphG987u8qo6RfYvonXAz4Jry6ib3cTjnXLp50GincrKaTphxxoHD+N3p+wEwfngvAHp02/YO8ugsf845l25pr56S9KUkq8uAjzprA3YmKGGWpawscdbBwxk/ohcj+hQCbNPDCqCqtmlJo6HRyM7yGZucc+mRiTaNi4BDgZfC5aOAmcBoSdea2b0ZOGeXMW5g9/jzxMZyCIYZ+WjF1ulhL7xrOudOHMHkfQa1Sv6cc51bJoJGDvAZM1sLIGkgcA9wCPAq4EEjRZcfN44DR/Zqdvuxew5gzqqm/Q3qGoxTb3k9vvzK/PVkKRiKZENFDQN7FGQqu865LiATbRrDYwEjtC5ctxHwWfx2wPeOG8sRY5u/o/2y48bx4pWf482rjmnxOGVVddz0wgIO+fULrC7b3n2WzjnXvEyUNF6W9ATwaLh8eriuCCjNwPm6rKwsMaZ/cZNBDJMprazj5U+C5qQ1ZdUM7tmtNbLnnOuEMlHS+DbBAIXjw8c9wLfNbIuZHZ2B83V5RUnmFr/utL3jzzdW1sZ7Y9X6IIfOuV2QiZLG5cDDZvb/MnBsl0R2luhXnEdhXg7LNgajx+81pEd8e1lVXbwHVXm1d8l1zu28TJQ0ugPPSnpN0qVhQ7jLsIf/91B+OHkPIJgBcFCkCsoMKmqCrrhlVd6s5JzbeWkPGmb2CzPbm6CaajDwiqTn030e19Ru/YsZ0D3oGTWiT2F8mPWDR/UGYHPY7hENGis2VfLUR6tbPO47i0ooqeh0Awo753ZSJu8IXwesAUqAARk8jwsVhkONjOhTSLe8bJZcfzLfOWYssDVYlFXWxtOf/OfX+b/738PMkh7PzDj7jrd9BF3nXFzag4akSyS9DLwA9AW+aWb7pfs8blvdIkEjJjaA4eawLaOsqo5fPTmXV+avjweS6AyAtfWNbKmpb7J+4bqKzGfeOdchZKIhfDhwmZm9n4Fjuxb0LcojLyeLfYZunZK9V8KghiVbanniw9X87bXF8XXVdQ0U5GZjZhzxuxfJz8nm1R8eTbUPfuicS5CJ+TSuBpA0ACiIrF+W7nO5pnoV5vH21cc2Gf02VtKIifWuiqqsbaBXIcxaXsra8pom651zLioTAxaeCtwADCFo1xgJzAP2bmk/lx6JQaJbbjb5OVnxqqYFa7etalpaUklDo1FSsbW9w8x8mHXn3DYy0RD+S2ASMN/MRgPHAt6S2kYk0btwayBJFgjO/dvbHPG7l6io2dqzqqa+cZsRc5/8cLW3bzjXxWUiaNSZWQmQJSnLzF4CJmTgPC5Fie0azYndywHBRE+vLljfZPu3H3iP4298Ja15c851LJkIGqWSiglGtL1f0k3Algycx6VoxaZgkML9h/VsMV2s1xTA1A9W8bunP4kv14TTzzbTO9c510VkImicBlQSDCfyNPApcGoGzuNSVBEGg7MOHt5iuug85L98cl6TbZfc9176M+ac63DS1hAuSRaIlSoagbuTpUnXOV1q7rlwIp+s2czYAd1bTLeqtPlh01/4eOuki4vWV/D6wg187dBR6cqic66DSGdJ4yVJ35E0IrpSUp6kYyTdDZzf3M6SJkv6RNJCSVcl2X6kpPck1Us6I4357vSOHNefbx45hr2H9GBIzwL2HJQ8eLQUNKIm3/QaUx6f4/dxONcFpTNoTAYagAclrZI0V9JiYAFwLvAnM7sr2Y6SsoFbgROBvYBzJe2VkGwZcAHwQBrz3KUU5efwxlXH8NT3jki6fWVpFT27bb/RPDa8+paaeu56Y3HKwcY51/GlrXrKzKqB24DbJOUC/YAqMytNYfeJwEIzWwQg6SGCtpG5keMvCbf5hBC7QFKz21aXVfOZwT1SHgl38YYtXPPfuTw6cwVPfjd5IHLOdS4ZGbDQzOrMbHWKAQNgKLA8srwiXLdTJF0saYakGevXr9/+Di6ub1Eej37r0JTSzlpWCsCmLbUtJ3TOdRqZHOW2zZjZHWY2wcwm9O/f/BzbbqtwjiaK8rMZ3a8opX1+NS3oYZWV1XzpxTnXuWRiwMKdsZJgoMOYYeE6lyHPXn4kjWa8v6yUPQf34O+vLeKJD1ezx6AedC/YsY9FVljl9cK8tdTUNzKsdzf2G9YrA7l2zrW19hI0pgNjJY0mCBbnAOe1bZY6t3EDgx5Uew4KpoW97rR9OGW/wRy/16D41LAA/730cE695XUArj5xT37z1MfbHGvZxkq+++Aspn6wKr5u8W9OarH9xDnXMbWL6ikzqwcuBZ4hGNzwETObI+laSV8AkHSwpBXAmcBfJc1puxx3Pr2L8pi8z+AmAQNg32E9ue6L+/Db0/elML/53xjRgAGw58+ebnZyJ+dcx9VeShqY2TRgWsK6KZHn0wmqrVwr++qkkQA89t6KlPepqW/krL++xb0XHUJBbnamsuaca2XtoqTh2p8rjh/H9V/at8m6sQO6k5stvn/CuJSOMX3JJj5eszkT2XPOtZF2U9Jw7ct3jx27zbp9h/Vkwa9OAuAPz85vdt/hfbqxfGNww9+mSu+O61xn4iUNl3aXHLU7z11+JAAbNtc02baqtIqz/vIWy5PMIOica/88aLid8pevHNTsth4FuQzt3Q0I5iSPeujdZby7ZCNfu/Pd+Oi7MS/MW+vjWTnXznnQcDtl8j6D+H//l/zO8e4FORTm5dAtN5t5q8vjw5KUVdXFe1kt3rCFv7+2KB4kFq2v4KK7Z7Dnz55m1rJNO52vR2YsZ31C6cY5lz4eNNxOO2hkH47dc8A264vCrrlVdQ08/v4qjvr9S8xfu5mTbnqNJSWVXHjYaAD+9PyC+EyAGyMlktte/nSn8rOytIof/utDLrl/5k7t75zbPg8abpfceM547rlwYpN1ifd6bKqsY/KfXmVlaRVfnTSSKafuxYg+hQAs31jFRyvKUh4k8cWP17Juc3XSbbFSy6rS5Nudc7vOg4bbJT0KcjlyXH/OOCi4haZ7QQ67DygG4Pdn7Mfhu/cDoNGC+z1+dOKeQNN5y0+95XV+9p/ZAAzpWRAvdVTXNfCXVz5l3upyHnx3GXNWlXHhXTM4/fY3k+Zlc3V9eC6/qdC5TPEuty4t/nDm/vzhzP2brDtzwnDGDezO6ws3AHD1SXtSmBd85BIHGFlVFpQORvcv4o2FJbz1aQlLS7Zw/VMfc304dMl+4RznyzdWcfE9M7j8+HHMWlbKuROHI4nN1UFpxYOGc5njJQ2XUYN7FsSfxwIGQH1j8i/2kX2DEXbP/dvbRIeuGtKzgHmry+PLz85dy4k3vcaP//0RH6woA7aWNNaW1/DrafOormtg+pKNabsW55wHDZdhfYvzk65vCIPGA988pMn6fpH0n6ypACAvJ4tDxvSlriF5oFlXHpRSYiUNgDteXcQv/juHM//yFktLtiTdzzm34zxouIzKzhK79S/iW5/bLen23oV5vPKDo+LLsalkAe58YzEAH11zAv2K85o9x8X3zuTdxRvjJY2YuauCksnqsp1vGC+rquPmFxbEg5xzXZ23abiMe+HKo7ZZd+PZ4/nbq4sYO6CYnOytv10K87Yd3DA/J7tJCSSZyx6aFW+Mj+8XDpR4zh1v8/DFk5gdBpGLDh/d7HE2VNQ0OdeUx2fz+PurGD+iF0eM9Qm9nPOShmsTnxncgxvOHt8kYHTLzebiI8dwwWdHbZN+e0FjfUUN5QkljdzsrY0iz8xZy3VPzOW6J+aysjQYF6u6roGbnl/Am58GDfXX/ncuE375PO8sKolvf/mTYLpgM/j547N56qPVKV/jnFVlvLvY21Rc5+JBw7ULs3/xeab/9DgKcrO56sQ9GdO/6ZSzfVqongKoazDmRhrKIZgcKmbqB1sngjzs+hepb2jk7jeXcOPz87n95U8xM+56M6gOmxnekf7jxz6K3z/yw399yN1vLU06CVVzTv7z65z117dSTu9cR+BBw7ULxfk5FId3khfkZvNiWKUVm3p278E9GDewmAe/OYlrTt0rvt/Tlx3BdaftDcC7izey79Ce/OUrBwLER9oF2FDRdAysQ69/Md4V+LUFG/jzCwuJNVvMX7OZh6cv47FZWwPNmrCxPT8nix/960NmryyLbyurrOP0299k4bqKpNdW19CYdH1MRU09f33lU283cR2COvvsahMmTLAZM2a0dTbcTpi9soy+xXkM7tltm20NjUZtfSPd8rLZXF3Hvtc8C8DPTtmLI8b244QbXwXggs+O4q43l6R8zr5FeU0GWTz/0JHc/dZSAMYOKGZBGBi+sP8Q/nzuAZRX1/HsnLV8/9EPOHGfQXz9sNG8t2wTdfWN/PG5YPj4w3fvx70XTdxm+tunZ68mOyuLV+ev5963l3LHVw/ihL0Hpf4COZchkmaa2YRk27wh3LVb+wzt2ey27CzRLWw0716Qy6Fj+mIY5xw8nKL8HL5/wjj+8Ox89hjUnbzsLGq382s/5htHjOG3TwdVUC9//yhG9i2MB40Jo3rHg8a0j1Yzqm8hf35xIUVhPuoaGpNWR72+cAOry6p5Zf56JozsTf/u+fQoyOVb973XJF2Vj/DrOgAPGq5TePDiSU2WLz1mLBcdPob8nCz2GdKTWcs3sWJTFXe8ugiAnt1yKauqIz8ni5qwm+/XDxvFRYePZt7qck7ZbzCj+jVtV5k0pi8PvrscCG5O/POLCwHYUht82T8/b12z+fvBvz7gjYUl8eVf/8++26RZG1aBNTYatQ2N5GZnbTOOF4CZsWxjJTOXbuJLB247A/KGihpE8/fItKVHZiwnLzuLLx4wtK2z4naSBw3XacVKIvsO68m+w3rS2Gj0K87jrAnD6dktl+fnrWPG0o389ZVFnD1hOFNO2QtJ/PncA5Ie7+R9B/PrafMozMth8YbkNwxmCW4970AueeA9zOBrh47knreWNgkYAD/+90fb7LtgbQXzVpdz/VMf88r89RTkZvHuT46jR0EuJRU1ZGeJp2ev4TdPfczm6joaDa545AOuPH4c/3hjMW9ffSwFudlceNd0PlxRxtOXHcGeg3ps93VaU1bNoMid+wDl1XX0KMhtZo+d98N/fQjgQaMD86DhuoysLHHxkVtvMjx+r4Hx6WgPGtV7mzaHRDnZWbzyg6OpbzQ+XVfBv2et5IWP17J8YxX9u+fz7aN24/P7DGJwz268/qNj+Ofri7ns+HHcE1ZvJfrxSXvy7Jy1zFga9NZ6dOYKHp25Ir69uq6Rc+94m+G9C3l27hp6dsulsrYhXjKKibWdfLC8lG552XwYDqvyxAer40GjsTHoXbbP0J7c9vJCxvQrZmVpFb265XLlox9w63kHcsS4flwzdQ5m8O9ZK7nnwolsqannxH0HNznfpi21vL+8lNkry8jPzWrymtY1NNLQaBTkbnu/TbSh38xafL3rGhrJydJ23xPX+tpVQ7ikycBNQDbwdzO7PmF7PnAPcBBQApxtZktaOqY3hLuWNDYaL89fx9F7DGj2C2r8tc9SWlnHkutPTrp9ackWCnKzGdijIOn2cT99isZGY2jvbiwtqeRPZ4+nX3E+h4/tF09z2UOz+M/7wQRVp+w3mOM+M5DLHn4/6fH+cf4Eausb2X94L069+fVtZkeMGdyzgOP3Gki/4nw2bqnlrjeXMKRnQXxwyESTxvTh7UVb7yvJzhINjcYNZ+3PaeOH8u7ijfz+mY95b1lpk/3uvGACv3v6E7526ChuemE+3XKzefK7R/DO4hKOGNuf3PBenJWlVRx2/YsAvPez41m0voLNNfWMH9aL6voGKqrrGdm3iJr6Bva95lmuPH4c30kyV/2OWrS+gtcWbOD8z47i6dlrmDi6D32KWu7C3dW11BDeboKGpGxgPnA8sAKYDpxrZnMjaS4B9jOzb0k6B/gfMzu7peN60HC7aktNPY1mdN/J6po1ZdUU5Wdz84sLuePVRbz742MZkBBgqusaMAu69GaF7RhvLNxAyZZadutfxH8/WM1fXgkmp5p77efjgz9e+9+58eFWon40eU8efHcZ5dV1lFamNlcJwOh+RUmr3sb0L2JpSeUOdwse0D2fz+89iGfnrmHvIT158eOg3SfaKy2qW252kw4BE0f1oV/3PIrzc3h/eSmfGdyDS47anSG9Cnh+3lrWb67h8N378+qC9byzqISqugZG9ytit/7FjB3YnX2G9ODgXz1Po8FZE4bxyIwVjOlfxFcnjeTAEb3ZXF3Pv2etpF9xHsP7FDK6XxF7DuoOwLTZaxjQPZ8JI3vzf/e9xxkHDWNYn24M6J7PoJ7d+HB5KbsPKOatRSVM3mcQ+TnZmBkrNlVRsqWW/Yf1pLy6nnmry8nLyeKA4b1oaDSq6hqorG1o8iMj1o4VK6EtXLeZ4vzceLVhtGTW0GjMWVXG6H5FTT6TZkZdg5GXs+t3UnSUoHEocI2ZfT5cvhrAzH4TSfNMmOYtSTnAGqC/tXARHjRce1Hf0MjK0qr4SL47alVpFavLqjhoZJ/4OjNjaUklw3p3Y1VpNWs3V1OQk82+w7b2PHtmzhqenr0GCcqr6jjjoGFc8cgHPHzxoUybvZp15TU0mlFV28C1p+3Np+u38PrC9WRnZXHvW0vYFAadSWOC8769aCMTR/fhR5P34Pw7p1Nd18DwPoUs3rCFn578GV7+ZD2vL9xAXk4WDY22S/efBNdVRUe4hSUnS01Gb87PyaKuoTGe94LcLOoatr4eWQruSSrMy6G+sZHSyjp6FzatgszPyaI4P4fSqjoG9SigrqGRTZW1QXDIzqIoPwgykqirb6S6voHehXnkZImCvGyeuezIeElvR3SUoHEGMNnMvhEufxU4xMwujaSZHaZZES5/GqbZkHCsi4GLAUaMGHHQ0qXJ65Sdcy2rqKknJ0u8t2wT44f3Iicr+CKMTem7pqyaHt1yqKlrZHHJFg4c0ZvSylpmLS/lqHH9Kauq45X56zlq3ADmrCqL/3KuqAmGfNlvWC9WbKrk7UUbmTSmD70L85j20WoOHBmUAvYb2pMVm6oor65j1vJSBnbPZ1VpFeXV9Rw+th8CZizZRGF+Nv2L8+lTlEdtQyNVtQ2sKq1ifUUtR+3Rn5q6RtZtrmaPQd35ZM1msiSyJHKyxRFj+1GQk807izfy1qIS6hoaKczN5rO796Osqpb3l5cxqEcBq8uq2HdoUHpYvKGC7gW5rNxUFe9wkZMlcrKz6F+cR052FovWV1CQm81+w3pRVlXHx6vLKcjNJjc7iy219UhBCaOipoH6hkaG9OrGhooacrODQJGVJTZuqaGu3uhdlMfa8mryc7LoVZhHcX426zfXxAOSYWRJ5OdkUVFTT12DUVPfyM3NdOrYni4XNKK8pOGcczumpaDRnoYRWQkMjywPC9clTRNWT/UkaBB3zjnXCtpT0JgOjJU0WlIecA4wNSHNVOD88PkZwIsttWc455xLr3Zzn4aZ1Uu6FHiGoMvtnWY2R9K1wAwzmwr8A7hX0kJgI0Fgcc4510raTdAAMLNpwLSEdVMiz6uBM1s7X8455wLtqXrKOedcO+dBwznnXMo8aDjnnEuZBw3nnHMpazc392WKpPXAzt4S3g9o9sbBTsqvuWvwa+4advaaR5pZ/2QbOn3Q2BWSZjR3V2Rn5dfcNfg1dw2ZuGavnnLOOZcyDxrOOedS5kGjZXe0dQbagF9z1+DX3DWk/Zq9TcM551zKvKThnHMuZR40nHPOpcyDRhKSJkv6RNJCSVe1dX7SSdKdktaFE1rF1vWR9JykBeHf3uF6Sfpz+Dp8KOnAtsv5zpE0XNJLkuZKmiPpe+H6znzNBZLelfRBeM2/CNePlvROeG0Ph1MQICk/XF4Ybh/VphewCyRlS5ol6YlwuVNfs6Qlkj6S9L6kGeG6jH62PWgkkJQN3AqcCOwFnCtpr7bNVVrdBUxOWHcV8IKZjQVeCJcheA3Gho+LgdtbKY/pVA9caWZ7AZOAb4fvZ2e+5hrgGDPbHxgPTJY0CfgtcKOZ7Q5sAi4K018EbArX3xim66i+B8yLLHeFaz7azMZH7sfI7GfbzPwReQCHAs9Elq8Grm7rfKX5GkcBsyPLnwCDw+eDgU/C538Fzk2WrqM+gMeB47vKNQOFwHvAIQR3BueE6+Ofc4I5bA4Nn+eE6dTWed+Jax0WfkkeAzwBqAtc8xKgX8K6jH62vaSxraHA8sjyinBdZzbQzFaHz9cAA8Pnneq1CKsgDgDeoZNfc1hN8z6wDngO+BQoNbP6MEn0uuLXHG4vA/q2aobT40/AD4HGcLkvnf+aDXhW0kxJF4frMvrZbleTMLm2Z2YmqdP1w5ZUDPw/4DIzK5cU39YZr9nMGoDxknoB/wb2bNscZZakU4B1ZjZT0lFtnJ3WdLiZrZQ0AHhO0sfRjZn4bHtJY1srgeGR5WHhus5sraTBAOHfdeH6TvFaSMolCBj3m9lj4epOfc0xZlYKvERQNdNLUuyHYvS64tccbu8JlLRuTnfZYcAXJC0BHiKoorqJzn3NmNnK8O86gh8HE8nwZ9uDxramA2PDXhd5BPOQT23jPGXaVOD88Pn5BPX+sfVfC3tdTALKIsXeDkFBkeIfwDwzuyGyqTNfc/+whIGkbgRtOPMIgscZYbLEa469FmcAL1pY6d1RmNnVZjbMzEYR/M++aGZfphNfs6QiSd1jz4ETgNlk+rPd1g057fEBnATMJ6gH/klb5yfN1/YgsBqoI6jTvIigLvcFYAHwPNAnTCuCnmSfAh8BE9o6/ztxvYcT1Pt+CLwfPk7q5Ne8HzArvObZwJRw/RjgXWAh8CiQH64vCJcXhtvHtPU17OL1HwU80dmvOby2D8LHnNh3VaY/2z6MiHPOuZR59ZRzzrmUedBwzjmXMg8azjnnUuZBwznnXMo8aDjnnEuZBw3nAEk/CUeE/TAcMfSQDJ/vZUkTtp8ynn5SOBrr+5LmSbomXP8FdbKRmF375sOIuC5P0qHAKcCBZlYjqR+Q18bZSnQ3cJaZfRCOxLwHgJlNpfPffOraES9pOBeMBLrBzGoAzGyDma0CkDRF0nRJsyXdEd5hHisp3ChpRvjL/2BJj4VzGPwyTDNK0seS7g/T/EtSYeLJJZ0g6S1J70l6NBwnK9EAgpsyMbMGM5sb7nuBpFvC5+9HHlWSPhfeNXyngvk1Zkk6LQOvn+tCPGg4B88CwyXNl3SbpM9Ftt1iZgeb2T5AN4ISSUytBXMY/IVgqIZvA/sAF0iKjZi6B3CbmX0GKAcuiZ44LNX8FDjOzA4EZgBXJMnjjcAnkv4t6X8lFSQmsGBOhfHAz8LjvAn8hGCIjInA0cDvwyEnnNspHjRcl2dmFcBBBBPTrAcelnRBuPnosC3hI4JB8PaO7BqrFvoImGNmq8PSyiK2Dgy33MzeCJ/fRzCsSdQkgsm+3giHMj8fGJkkj9cCEwgC3HnA08muRdJY4PcEVVl1BOMRXRUe+2WC4TNGtPByONcib9NwjvhQ4i8DL4cB4nxJDwG3EYzRszxsfI7+wq8J/zZGnseWY/9bieP0JC4LeM7Mzk0hj58Ct0v6G7A+UpoJDhRUaz0CfNO2DkQn4HQz+2R7x3cuFV7ScF2epD3CX+gx44GlbA0QG8Iv5DMS903BiLChHYISwusJ298GDpO0e5iXIknjkuTx5Fh7CsF0nQ1AaUKyO4F/mtlrkXXPAN+JtMUcsBPX4FyclzScg2Lg5nA48XqCkU8vNrPS8Ff9bIIZ0KbvxLE/IZiX/E5gLgnzMpvZ+rAq7EFJ+eHqnxKMshz1VeBGSZVhHr9sZg2xOCJpJEFQGyfpwnCfbwDXEcxo96GkLGAxTdtlnNshPsqtcxmiYHrZJ8JGdOc6Ba+ecs45lzIvaTjnnEuZlzScc86lzIOGc865lHnQcM45lzIPGs4551LmQcM551zK/j//xrwKJCX2twAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plt.plot(sim)\n",
"plt.title('Chi-squared p-value and sample size (stochastic)')\n",
"plt.xlabel('Sample Size')\n",
"plt.ylabel('p-value\\n(avg over 100 replications)')\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"sim = []\n",
"for n in range(10,5_001,10):\n",
" rep = []\n",
" a_sum = prob_a * n\n",
" b_sum = prob_b * n\n",
" table = np.array([[a_sum, n-a_sum],[b_sum, n-b_sum]])\n",
" test = stats.chi2_contingency(table)\n",
" rep.append(test[1])\n",
" sim.append(np.mean(rep))"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEWCAYAAABrDZDcAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAqWUlEQVR4nO3debxdZX3v8c/3zBlO5hOGjEDCEJAxMhR6QYoWUMFaKaBWVK7YQdTiteLVUkR7q7UVteKA1aKiII5NKZahoiLIEGSQJAwhAkkIZJ6nM/zuH+vZYWVzzslOOPvsnL2+71fOK2t41tq/tffa67ef51mDIgIzMyuuhloHYGZmteVEYGZWcE4EZmYF50RgZlZwTgRmZgXnRGBmVnCFTgSSrpR0fT/z50k6bfAiqj5J10n6VI1jmC4pJDXVMo5d2dX+sTdI7+OMAV7nRkkHDuQ603rvlnRMH/NOk7RkoF9zoOzOsWBPjxuS3ibptgrKXSrpM7u7/v7UfSKQ9FZJc9POvUzSzySdUsmyEXF4RPyiyiGa7VUiYmRELBrIdUp6I7AhIh4agHUNeoLenWNBJWV7+zEUEd+NiNdV8BJfB94maWIl8VSirhOBpMuAzwP/D9gHmAp8GTi3hmENmL39F7VZzl8A36l1EDD0vzcRsRX4GfCOgVxpXf4Bo4GNwHn9lLkSuAn4NrABmAfMzs1/Bjijj2XbgOuBVcBa4AFgnzTvAOCXaZ23A18Crk/zTgOWlK1rx+sAxwO/SetclpZtyZUN4K+Bp4Dfp2lvAB5Oy9wDHJkrfwzw2xTL94EbgU/1sU3vBO5Or7kOeBz4oz7Kng/MLZv2N8CcNPx64CFgPbAYuDJXbnrajqbe3uf0uVyfGz8xbdda4BHgtH4+08uBp9P2zgf+pGz7fg38M7AG+D1wVm5+n59bL68zAbg5xbQauAtoqDCGu4Gr07KLgD9I0xcDy4GLcuWvA76a4tmQ4ptWtj/MSMOtadueA15Myw3rI/4ZaV3rgJXA98vXCexP9h0q/W0GIlfu3cCC9F7emo+r7LVagC3A5Ny0YWnb1qT36MPkvhfptX8ErEif0/vT9DOB7UBniumR3Pf9G2TfmaXAp4DGXt7zVWnedWQ/Cn+W1nM3sC/ZD8c1ZPv+MX18R6+kwuMG2fd5Ltn34EXgc2n6c+l9Lr23J6U4f51bz+Hpc1+dlv2/uXlvA+4csOPlQK1ob/tLO0wX6WDTR5krga3A2UAj8I/Avb19oL0s+17gP4HhadnjgFFp3m+Az5F9Mf9X2lkqTQTHkR34msgOmAuAD5Z9SW8HxpF9mY4hO3ickOK4KK2vlewL+CzZAboZeAvZF6i/RNCVK38+2YFiXC9lh6ftmpmb9gBwQW47X0VW6zwy7chvSvOmU2EiACaRfXnPTut6bRrv6GMbziM7iDSk+DcB++W2rxN4T3qv/hJ4HtCuPrdeXucfyQ60zenvD3Pr2VUMXcC7UgyfIjsoXJNe93XpdUem8tel8f+V5n+BnQ8W+URwNTAn7RvtZPvnP/YR/w3Ax1KMbcApva2zbJnvAjek4XOBhcBhZPvqx4F7+nitw4FNZdM+TZY8xwFTgMdI34sU04PAFWT78IFkCfOPy/eP3Pp+AnwNGAFMBO4H3lv2nl+aYi0loZVk37c24OdkCecduc/lzj6+o1dS4XGDbJ/68zQ8Ejixt+9ALs5fp+F2sqT2oRRfO3BCruyxwOoBO14O1Ir2tj+yjPnCLspcCdyRG58FbOntA+1l2XdT9us7TZ+adroRuWnfo8JE0MvrfBD4SW48gNNz418BPlm2zBPAqWQHjx0HujTvHvpPBOXl7y/tyL2Uvx64Ig3PJDtgDe+j7OeBq9PwTl+C8u1n50TwEeA7Zeu6ldyv5l18xg8D5+a2b2Fu3vAUx767+tx6We9VwH/QywGzghieys17VYphn9y0VcDRafg64MbcvJFANzAltz/MAESWcA7KlT2JVGvsJaZvA9eS+5Veto/NKJv2EbKD87A0/jPg4tz8BrIaw7Re1ncyZd9FsgP7mbnxS3gpEZwAPFdW/qPAv5fvH2l8H2AbudoPcCHpQJ7e8/L1XQd8PTd+KbCg7HNZmxvfsY+yG8cN4FfAJ4AJZa8/nf4TwYXAQ/3sUzOB7kq+A5X81XMfwSpgQgXtgS/khjcDbb0tkzqbS39Tydo7bwVulPS8pH+S1Ez2S3BNRGzKLf5spUFLOljSzZJekLSerH9jQlmxxbnhacCHJK0t/ZH9wto//S2NtOdUGEtv5feX9Ie57Z+X5n2PbIcFeCvw04jYnLbjBEl3SlohaR1ZG3H5dlRiGnBe2fadAuzXW2FJ75D0cK7sEWWvu+PzLsVKdnDd3c/ts2S/iG+TtEjS5bsRw4u54S0plvJpI3PjOz7viNhI1lSwf1k8HWSJ7cHc6/53mt6bvyVLHvens1ze3deGSjoL+ABZjW5LmjwN+ELutVan9U3qZRVryH7R5u3Pzvtx/r2eRrbP5T/z/0t2wO/NNLJa2bJc+a+R1QxKFveyXPl73t9nUK6i4wZwMXAw8LikByS9oZ915k0ha17sSztZbX1A1HMi+A3Zr4Q3DcTKIjuTovT3XER0RsQnImIWWRvvG8iqlcuAsZJG5BafmhveRPaFBUBSIzt/Wb9C1j45MyJGkX0BVB5Obngx8A8RMSb3NzwibkixTJKUX34q/eut/PMRcVdu+w9P824HOiQdTZYQvpdb7ntkzRRTImI0WTNK+XaU7PSekP1Cz2/fd8q2b0REfLp8JZKmkZ1R8T5gfESMIWty6Ot183b1ue0kIjZExIci4kDgHOAySX/0CmPoy5TSgKSRZM0pz5eVWUl28Do89z6NjoheD2YR8UJEvCci9idr5vxyb6ehSjoE+BbwZxGRP5guJmt6yX8uwyLinl5ebmG2KuWTxLL8drHze72YrCaTX3d7RJxdCr9s/YvJvusTcuVH5fbT3pYZFBHxVERcSJaUPgP8MO1ju4pnMVmTWF8OI+svGxB1mwgiYh1ZG+M1kt4kabikZklnSfqnV7p+Sa+R9Kp0IF9P1vbcExHPknUOfUJSSzpV9Y25RZ8k+/Xw+lSD+DhZ229Je1rfRkmHkrVj9+frwF+kX+CSNCKtu50sGXYB70/b/mayzqv+TMyVP49sh7ult4IR0Qn8gOzX8TiyxJDfjtURsVXS8WQ1hr48DFyQXnM2WV9GyfXAGyX9saRGSW3Kzjmf3Mt6Sl+wFQCS3kX2a3yXKvjcdiLpDZJmpKS5jqy5pueVxNCPsyWdIqkF+CRZe/ROv3AjoodsX7i6dFqhpEmS/riP+M/LvYdrUsw9ZWVGkTV/fSwifl22iq8CH5V0eCo7Ou0vLxMR24E7yJorS25Ky49NcVyam3c/sEHSRyQNS5/7EZJenea/CEyX1JDWvwy4DfgXSaMkNUg6SFL+9WpC0tsldaTPZ22a3EO2f/TQ98H+ZmA/SR+U1CqpXdIJufmnkjXPDYi6TQQAEfEvwGVkB9sVZFn2fcBPB2D1+wI/JDtoLyA7A6N0etxbydo5VwN/T9YeW4ppHfBXwL+Rnd2wCchfSPN/0vIbyL7Y3+8viIiYS9b5+SWyL/RCsrbG0hfwzWl8NVnH5Y93sV33kbU/rgT+AXhLRKzqp/z3gDOAH0REV276XwFXSdpAlpBv6mcdfwcclOL/BLmaRTrgnUtWMyp9hh+ml303IuYD/0KWAF8ka+e9u5/XLdfn59aLmWQHt43p9b4cEXcOQAy9+V6KZzVZ5+bb+yj3EbLP/15lzYp3AIf0UfbVwH2SNpLV3D4QL7924Ni0/NXKNY0CRMRPyH7h3phe6zHgrH624WvAn+fGP0HWHPR7soP4jlNLI6KbrIZ9dJq/kuz7MjoV+UH6f5Wk36bhd5B1LM8n249+SB/Nh4PsTGBeet++QHYyxZbULPkPwN2pOevE/EIRsYHsxIg3kjVDPQW8BkBSG1lH9bcGKsjSWQ5WRZKuJOt86+sLvFeQ9E7gf0dERRfcWfVJuo6sE/XjtY7llZJ0N/C+GICLyopM0qVkTa5/O1DrHNIXVpjZ0BERJ9c6hnoQEf860Ous66YhMzPbNTcNmZkVnGsEZmYFN+T6CCZMmBDTp0+vdRhmZkPKgw8+uDIier3AcMglgunTpzN37txah2FmNqRI6vNKeTcNmZkVnBOBmVnBORGYmRWcE4GZWcE5EZiZFVzVEoGkb0paLumxPuZL0hclLZT0qKRjqxWLmZn1rZo1guvI7rzXl7PI7uA4k+zpRF+pYixmZtaHqiWCiPgV2W1z+3Iu8O3I3AuMkVS128Y+8MxqPnfbE8x/fn21XsLMbEiqZR/BJHZ+fNwSen/MHZIukTRX0twVK1bs0Yv99tk1fPHnC/nKL/t7+puZWfEMic7iiLg2ImZHxOyOjr4ewdq/9556EIftN4qtnd0DHJ2Z2dBWy0SwlJ2fWTo5TaualqYGtnX17LqgmVmB1DIRzAHekc4eOhFYl549WjWtTQ1sc43AzGwnVbvpnKQbgNOACZKWkD1ztRkgIr5K9kD0s8mesboZeFe1YilpbWpg47auXRc0MyuQqiWCiLhwF/MD+OtqvX5vWpsaWLXRTUNmZnlDorN4oLQ2NbKty01DZmZ5BUsE7iw2MytXrETQ3MB2JwIzs50UKhG0NLpGYGZWrlCJoLXZfQRmZuWKlQiasqah7IQlMzODgiWClsYGegK6epwIzMxKCpUIWpuzzXU/gZnZS4qVCJoaAXzmkJlZTqESQUtTqUbgDmMzs5JCJYLWUiLodI3AzKykYIkgNQ11OxGYmZUUKhG0uEZgZvYyhUoEre4jMDN7mYImAtcIzMxKCpUIhrVkfQRbtrtGYGZWUqxE0Jwlgs1+XKWZ2Q7FSgSpRrDVNQIzsx0KlQiGt2RP5ty83c8tNjMrKVQiKDUNbfHpo2ZmOxQqEbQ1NyDBFtcIzMx2KFQikMSw5kY2u4/AzGyHQiUCgOEtjWzxWUNmZjsULhG0NTf6OgIzs5zCJYLhLW4aMjPLK1wiGNbS5KYhM7Oc4iWC5gY3DZmZ5RQuEQx3jcDMbCeFSwTDWhp9ZbGZWU7hEsFwX0dgZraTwiWC9rZmNm51jcDMrKSqiUDSmZKekLRQ0uW9zJ8q6U5JD0l6VNLZ1YwHYGRbExu3d9HTE9V+KTOzIaFqiUBSI3ANcBYwC7hQ0qyyYh8HboqIY4ALgC9XK56SUW1NRMBG9xOYmQHVrREcDyyMiEURsR24ETi3rEwAo9LwaOD5KsYDwMjW7FbUbh4yM8tUMxFMAhbnxpekaXlXAm+XtAS4Bbi0txVJukTSXElzV6xY8YqCam9rBmCDE4GZGVD7zuILgesiYjJwNvAdSS+LKSKujYjZETG7o6PjFb1ge1uqEWzrfEXrMTOrF9VMBEuBKbnxyWla3sXATQAR8RugDZhQxZgYmRLBetcIzMyA6iaCB4CZkg6Q1ELWGTynrMxzwB8BSDqMLBG8srafXRiVEoGbhszMMlVLBBHRBbwPuBVYQHZ20DxJV0k6JxX7EPAeSY8ANwDvjIiqntc5sjXrI3BnsZlZpqmaK4+IW8g6gfPTrsgNzwdOrmYM5dp31AjcR2BmBrXvLB50w1saaRBs3OYagZkZFDARSGJka5P7CMzMksIlAsiuJVjvpiEzM6CwiaDJncVmZklhE4GbhszMMoVMBCNbm9xZbGaWFDIRtLc1+/RRM7OkoInATUNmZiWFTAQj25rY4KYhMzOgoIlgVFsz27t62NblZxebmRUyEZQeTuPmITOzgiaCMcOzG8+t3ewOYzOzQiaCcSNaAFizeXuNIzEzq71CJoKxw7NEsHqTE4GZWSETwY4agROBmVkxE8GOGoGbhszMipkIhrU00tbc4M5iMzMKmggAxg1vcR+BmRkFTgRjR7S4j8DMjAIngnEjWtxHYGZGgRPB2OGuEZiZQaETQTNr3FlsZlbgRDCihXVbOunq7ql1KGZmNVXYRFC6qGztFtcKzKzYCpsISheVuZ/AzIqusImgVCPwtQRmVnSFTQS+8ZyZWaawiaCjvRWAFRu31TgSM7PaKmwiGD+ihcYGsXy9E4GZFVthE0FDg5gwsoXlG7bWOhQzs5qqaiKQdKakJyQtlHR5H2X+TNJ8SfMkfa+a8ZTraG9l+QbXCMys2JqqtWJJjcA1wGuBJcADkuZExPxcmZnAR4GTI2KNpInViqc3E9vbeHG9awRmVmzVrBEcDyyMiEURsR24ETi3rMx7gGsiYg1ARCyvYjwvM9E1AjOzqiaCScDi3PiSNC3vYOBgSXdLulfSmb2tSNIlkuZKmrtixYoBC3BieyurNm6juycGbJ1mZkNNrTuLm4CZwGnAhcDXJY0pLxQR10bE7IiY3dHRMWAv3tHeSk/AKp9CamYFVs1EsBSYkhufnKblLQHmRERnRPweeJIsMQyKjvY2ADcPmVmhVTMRPADMlHSApBbgAmBOWZmfktUGkDSBrKloURVj2snEUemiMicCMyuwihOBpFMkvSsNd0g6oL/yEdEFvA+4FVgA3BQR8yRdJemcVOxWYJWk+cCdwIcjYtWebMiemJiuLva1BGZWZBWdPirp74HZwCHAvwPNwPXAyf0tFxG3ALeUTbsiNxzAZelv0E0YmRKBry42swKrtEbwJ8A5wCaAiHgeaK9WUIOlrbmRscObecHXEphZgVWaCLanX+8BIGlE9UIaXPuPGcbza7fUOgwzs5qpNBHcJOlrwBhJ7wHuAL5evbAGT5YIXCMws+KqqI8gIv5Z0muB9WT9BFdExO1VjWyQTBozjHufHrT+aTOzvU7F9xpKB/66OPjn7T+mjQ3buli/tZNRbc21DsfMbNBV1DQkaYOk9elvq6RuSeurHdxg2G/0MACWuXnIzAqq0qahHWcISRLZzeNOrFZQg2n/MVkieH7tFg7Zd8ifCGVmttt2+8riyPwU+OOBD2fwTUqJYKnPHDKzgqr0grI350YbyC4uq4u2lI72Vpoa5FNIzaywKu0sfmNuuAt4hpc/W2BIamwQ+45uY9m6ushrZma7rdI+gndVO5Ba2n/MMJas2VzrMMzMaqLfRCDpX0lXE/cmIt4/4BHVwJSxw/n1woF74I2Z2VCyqxrB3EGJosamjR/Oj367ja2d3bQ1N9Y6HDOzQdVvIoiIbw1WILU0bfxwAJ5bvZmD9/EppGZWLJWeNdQBfASYBbSVpkfE6VWKa1BNHZclgmdXORGYWfFUeh3Bd8keLnMA8Amys4YeqFJMg27a+Oxmqs+tdoexmRVPpYlgfER8A+iMiF9GxLuBuqgNAIwd3kx7axPPrdpU61DMzAZdpdcRdKb/l0l6PfA8MK46IQ0+SUwdP5xnXSMwswKqNBF8StJo4EPAvwKjgL+pWlQ1MHXccJ54YUOtwzAzG3SVJoL7ImIdsA54TRXjqZlp40dwx4IX6eruoalxt2/BZGY2ZFV6xLtb0m2SLpY0tqoR1ciBHSPo7A6WrPE9h8ysWCpKBBFxMPBx4HDgQUk3S3p7VSMbZAd1jATg6RUbaxyJmdngqrgNJCLuj4jLgOOB1UBdXWw2w4nAzAqq0ieUjZJ0kaSfAfcAy8gSQt0YPbyZCSNbeXq5TyE1s2KptLP4EeCnwFUR8ZvqhVNbB3WMYKFrBGZWMJUmggMjIgAkvSEibq5iTDVz0MSR/Nejy4gIsidympnVv0o7i/O3or6qSrHU3EEdI1m3pZPVm7bXOhQzs0GzJyfM1+1P5RkTsw7jJ19085CZFUelncVtki6T9GNgjaS/kdS2ywWHmMP2ze48+vgL62sciZnZ4Km0j+DbwAay20sAvBX4DnBeNYKqlY72VsaNaGHBMicCMyuOShPBERExKzd+p6T51QioliRx2H7tLFjmew6ZWXFU2kfwW0knlkYknUAFj7GUdKakJyQtlHR5P+X+VFJIml1hPFVz2L6jePLFDXR199Q6FDOzQVFpIjgOuEfSM5KeAX4DvFrS7yQ92tsCkhqBa4CzyJ5sdqGkWb2Uawc+ANy3B/EPuEP3G8W2rh6e8bMJzKwgKm0aOnMP1n08sDAiFgFIuhE4FyhvUvok8Bngw3vwGgPusP2yDuP5yzYwY6IfW2lm9a/S6wie7e+vj8UmAYtz40vStB0kHQtMiYj/6u/1JV0iaa6kuStWrKgk5D02Y+JImhrE4+4wNrOCqNmN9yU1AJ8je9hNvyLi2oiYHRGzOzo6qhpXa1MjMyaO9JlDZlYY1UwES4EpufHJaVpJO3AE8IvU73AiMGdv6DA+dF+fOWRmxVHNRPAAMFPSAZJagAuAOaWZEbEuIiZExPSImA7cC5wTEbs8G6naZu0/ihfWb2Xlxm21DsXMrOqqlggiogt4H3ArsAC4KSLmSbpK0jnVet2BcPSU7CFsDz+3traBmJkNgkrPGtojEXELcEvZtCv6KHtaNWPZHa+aNJrGBvHw4rWcMWufWodjZlZVfkp7L4a1NHLovu08vHhtrUMxM6s6J4I+HD1lDI8sXktPT+y6sJnZEOZE0Iejp4xhw7YuP8PYzOqeE0EfjpmadRg/5OYhM6tzTgR9OHDCCNrbmnjIZw6ZWZ1zIuhDQ4M4esoYHnpuTa1DMTOrKieCfrx6+jieeHEDazf7GcZmVr+cCPpx0kHjiYB7F62udShmZlXjRNCPIyePpq25gXsXrap1KGZmVeNE0I/WpkZmTxvnRGBmdc2JYBdOPHAcj7+wgdWb3E9gZvXJiWAXTjpoPAD3uVZgZnXKiWAXjpw8hmHNjfzGicDM6pQTwS40NzZwwoHj+NWT1X1EpplZrTgRVOD0QyfyzKrNLPJ9h8ysDjkRVOA1h0wE4M4nXCsws/rjRFCBKeOGM2PiSO58fHmtQzEzG3BOBBU6/dCJ3Pf7VWza1lXrUMzMBpQTQYVec8hEOruDXy9cWetQzMwGlBNBhWZPH0t7WxO3zXux1qGYmQ0oJ4IKNTc28LpZ+3L7/BfY3tVT63DMzAaME8FueP2R+7J+axd3u3nIzOqIE8FuOGVGB+1tTfzX75bVOhQzswHjRLAbWpqy5qHb5rl5yMzqhxPBbio1D931lC8uM7P64ESwm06Z0cG4ES386LdLah2KmdmAcCLYTS1NDfzJMZO4ff6LfkaBmdUFJ4I9cN7syXR2B//x8NJah2Jm9oo5EeyBQ/cdxasmjeYHc908ZGZDnxPBHjpv9mTmL1vPo0vW1joUM7NXpKqJQNKZkp6QtFDS5b3Mv0zSfEmPSvofSdOqGc9AetMxkxjR0sh19zxT61DMzF6RqiUCSY3ANcBZwCzgQkmzyoo9BMyOiCOBHwL/VK14Btqotmbectxkbn5kGSs2bKt1OGZme6yaNYLjgYURsSgitgM3AufmC0TEnRGxOY3eC0yuYjwD7h1/MJ3t3T3ccP9ztQ7FzGyPVTMRTAIW58aXpGl9uRj4WW8zJF0iaa6kuStW7D0Xch3UMZJTD+7gO/c+y7au7lqHY2a2R/aKzmJJbwdmA5/tbX5EXBsRsyNidkdHx+AGtwvv+cMDWbFhGz980GcQmdnQVM1EsBSYkhufnKbtRNIZwMeAcyJiyDW2nzxjPEdPGcNXfvE0nd2+/5CZDT3VTAQPADMlHSCpBbgAmJMvIOkY4GtkSWBIPhBYEpeePoMla7bwHw8/X+twzMx2W9USQUR0Ae8DbgUWADdFxDxJV0k6JxX7LDAS+IGkhyXN6WN1e7XTD53IrP1G8aWfP+VagZkNOU3VXHlE3ALcUjbtitzwGdV8/cEiiQ+97mAu/tZcbrz/Of78pOm1DsnMrGJ7RWdxPTj90ImccMA4Pn/HU2zc1lXrcMzMKuZEMEAk8dGzD2PVpu187ZdP1zocM7OKOREMoKOnjOGco/bna79axDMrN9U6HDOzijgRDLCPvf4wWhobuGLOPCKi1uGYme2SE8EA22dUG5e99mB+9eQKP+TezIYEJ4IqeMdJ03jVpNH83U8fY/mGrbUOx8ysX04EVdDU2MDV5x/F5u3dXP6j37mJyMz2ak4EVTJjYjuXn3UoP398OTfcv3jXC5iZ1YgTQRVddNJ0Tpkxgatunsf859fXOhwzs145EVRRQ4P43PlHMXpYM++9fi5rNm2vdUhmZi/jRFBlE9vb+Orbj+PFddt4/40P0eV7EZnZXsaJYBAcM3Usn3zT4dz11Er+3tcXmNlepqo3nbOXnP/qqfx+5Wa++sunGT+ihcted0itQzIzA5wIBtVHzjyENZu288WfL2T08BYuPuWAWodkZuZEMJgk8Q9/cgTrt3byyZvn09ndw1+celCtwzKzgnMfwSBramzgixcewxuP2p9P/+xx/vnWJ9xnYGY15RpBDTQ3NvD5849mREsjX7pzIcvWbeX/vfkIWpsaax2amRWQE0GNNDaIf3zzq9h3dBufv+Mpnlm1ia++/Tg62ltrHZqZFYybhmpIEh8842CueeuxzHt+Ha//4l3cvXBlrcMys4JxItgLvP7I/fjxX55Me1sTb//GfXz6Z4+zrau71mGZWUE4EewlZu0/iv+89BTOnz2Fr/7yac76wl3cu2hVrcMyswJwItiLDG9p4tN/eiTfevfxdHb3cMG19/LhHzzC8vV+poGZVY8TwV7o1IM7uO2Dp/KXpx3ETx5ayqmf/QWfvfVx1m/trHVoZlaHNNTOYZ89e3bMnTu31mEMmmdWbuJfbn+S/3zkeUYPa+aiP5jORSdNY/xIn11kZpWT9GBEzO51nhPB0PDY0nV8/o6nuGPBi7Q2NfCW4ybzthOmMWv/UbUOzcyGACeCOrJw+Ub+7a5F/PihpWzv6uFVk0bzZ7Mnc85Rkxg9vLnW4ZnZXsqJoA6t3bydnz60lO/PXcKCZetpbhR/cNAEzjxiX147ax8muOnIzHKcCOrcY0vXMeeR5/nvx17gudWbaRAcOXkMJ88Yz8kzJnDs1LG0Nfv2FWZF5kRQEBHBgmUbuHXeC9z11AoeWbKO7p6gtamBY6eO5agpYzh6ymiOmjKGfUe1IanWIZvZIHEiKKgNWzu5b9Fq7n56JQ8+u4YFy9bT2Z193h3trRyx/yhm7tPOjI6RzNhnJDMmjmRUm/sZzOpRf4nAN52rY+1tzZwxax/OmLUPAFs7u1mwbD2PLlnHI4vX8vgLG7jn6VVs63rpOcr7jGplytjhTBo7jMljhzFpzPDs/7HD6Ghvpb21yTUJszpT1UQg6UzgC0Aj8G8R8emy+a3At4HjgFXA+RHxTDVjKrK25kaOmTqWY6aO3TGtuydYsmYzT724kaeWb+TpFRtZsmYzDz67hpsfXUZ3z841xtamBiaMbGXCyJb0fysT2lsYO7yFUcOaGdXWzKhhTYxqa2Z0Gh/Z1kRjg5OH2d6qaolAUiNwDfBaYAnwgKQ5ETE/V+xiYE1EzJB0AfAZ4PxqxWQv19ggpo0fwbTxI3bUHEq6unt4ccM2lq7ZwtK1m1mxYRsrN25n5YZtrNi4jefXbeXRpetYvWn7yxJGuZGtTQxraWRYc/pr2fn/tuZGhrdk4y2NDTQ3NtDcpJeGGxtobhQtTWXjjQ00peGmhgYaGqBBorFBO/5vlJCybc1PbxA0pPmNDanMjmEnLiuOatYIjgcWRsQiAEk3AucC+URwLnBlGv4h8CVJiqHWcVGnmhobmDRmGJPGDAPG9VmupyfYsK2L9Vs6Wb+1k/VbutiwtZP1W3eetqWzm62d3Wze3sWWzh62bu9mxYZtbN7exdbOHrakedu6eqj1HiBlCQVAaVyI9I9SnhBK87Lbiqu0ALlpZfNfyjH5eS+tK79sKZbdyUti95LY7q17d9ZbeendTrtDMeYB8P4/mskbj9p/wNdbzUQwCVicG18CnNBXmYjokrQOGA/sdFN+SZcAlwBMnTq1WvHaHmpoEKOHZU1BA6W7J+js7kl/2fD2rrLx7h46u14a7+oJunuCiKA7suGeCHp6oDuCnp7Y8X9PsGN+dxovDZemR0BQ+p8d46QklU3rZX4apjQ/Ny92WjYbinj5sr29XiV2N3/uzm+u3Vn37iTyvSXm3Xufa/NLZSC/Y3lDorM4Iq4FroXsrKEah2ODIGvGafT1D2aDoJp3H10KTMmNT07Tei0jqQkYTdZpbGZmg6SaieABYKakAyS1ABcAc8rKzAEuSsNvAX7u/gEzs8FVtaah1Ob/PuBWstNHvxkR8yRdBcyNiDnAN4DvSFoIrCZLFmZmNoiq2kcQEbcAt5RNuyI3vBU4r5oxmJlZ//yEMjOzgnMiMDMrOCcCM7OCcyIwMyu4IXcbakkrgGf3cPEJlF21XADe5mLwNhfDK9nmaRHR0duMIZcIXglJc/u6H3e98jYXg7e5GKq1zW4aMjMrOCcCM7OCK1oiuLbWAdSAt7kYvM3FUJVtLlQfgZmZvVzRagRmZlbGicDMrOAKkwgknSnpCUkLJV1e63gGiqRvSlou6bHctHGSbpf0VPp/bJouSV9M78Gjko6tXeR7TtIUSXdKmi9pnqQPpOl1u92S2iTdL+mRtM2fSNMPkHRf2rbvp1u+I6k1jS9M86fXdAP2kKRGSQ9JujmN1/X2Akh6RtLvJD0saW6aVtV9uxCJQFIjcA1wFjALuFDSrNpGNWCuA84sm3Y58D8RMRP4nzQO2fbPTH+XAF8ZpBgHWhfwoYiYBZwI/HX6POt5u7cBp0fEUcDRwJmSTgQ+A1wdETOANcDFqfzFwJo0/epUbij6ALAgN17v21vymog4OnfNQHX37eyZq/X9B5wE3Job/yjw0VrHNYDbNx14LDf+BLBfGt4PeCINfw24sLdyQ/kP+A/gtUXZbmA48FuyZ4CvBJrS9B37OdlzQE5Kw02pnGod+25u5+R00DsduJnsefF1u7257X4GmFA2rar7diFqBMAkYHFufEmaVq/2iYhlafgFYJ80XHfvQ2oCOAa4jzrf7tRM8jCwHLgdeBpYGxFdqUh+u3Zsc5q/Dhg/qAG/cp8H/hboSePjqe/tLQngNkkPSrokTavqvj0kHl5vey4iQlJdniMsaSTwI+CDEbFe0o559bjdEdENHC1pDPAT4NDaRlQ9kt4ALI+IByWdVuNwBtspEbFU0kTgdkmP52dWY98uSo1gKTAlNz45TatXL0raDyD9vzxNr5v3QVIzWRL4bkT8OE2u++0GiIi1wJ1kTSNjJJV+0OW3a8c2p/mjgVWDG+krcjJwjqRngBvJmoe+QP1u7w4RsTT9v5ws4R9PlfftoiSCB4CZ6YyDFrJnI8+pcUzVNAe4KA1fRNaGXpr+jnSmwYnAulx1c8hQ9tP/G8CCiPhcblbdbrekjlQTQNIwsj6RBWQJ4S2pWPk2l96LtwA/j9SIPBRExEcjYnJETCf7vv48It5GnW5viaQRktpLw8DrgMeo9r5d646RQeyAORt4kqxd9WO1jmcAt+sGYBnQSdY+eDFZ2+j/AE8BdwDjUlmRnT31NPA7YHat49/DbT6FrB31UeDh9Hd2PW83cCTwUNrmx4Ar0vQDgfuBhcAPgNY0vS2NL0zzD6z1NryCbT8NuLkI25u275H0N690rKr2vu1bTJiZFVxRmobMzKwPTgRmZgXnRGBmVnBOBGZmBedEYGZWcE4EVpckfSzdpfPRdBfHE6r8er+QVPFDxSWdmO6S+bCkBZKuTNPPUR3dHdeGBt9iwuqOpJOANwDHRsQ2SROAlhqHVe5bwJ9FxCPp7riHAETEHOr7YkfbC7lGYPVoP2BlRGwDiIiVEfE8gKQrJD0g6TFJ16arlEu/6K+WNDf9Qn+1pB+n+79/KpWZLulxSd9NZX4oaXj5i0t6naTfSPqtpB+keyKVm0h2ISAR0R0R89Oy75T0pTT8cO5vi6RT05Wn31T2bIKHJJ1bhffPCsaJwOrRbcAUSU9K+rKkU3PzvhQRr46II4BhZDWHku2R3f/9q2SX8P81cATwTkmlO1keAnw5Ig4D1gN/lX/hVPv4OHBGRBwLzAUu6yXGq4EnJP1E0nsltZUXiOx+9EcDf5fWcw/wMbLbJxwPvAb4bLoVgdkecyKwuhMRG4HjyB7UsQL4vqR3ptmvSW3zvyO7kdnhuUVLTTK/A+ZFxLJUq1jESzf2WhwRd6fh68lud5F3ItnDj+5Ot4y+CJjWS4xXAbPJktZbgf/ubVskzQQ+S9aM1El275nL07p/QXZrhan9vB1mu+Q+AqtLkd2y+RfAL9JB/yJJNwJfJrsfy+LUQZv/Jb4t/d+TGy6Nl74r5fdkKR8XcHtEXFhBjE8DX5H0dWBFrtaRrShrUroJeE+8dCMxAX8aEU/sav1mlXKNwOqOpEPSL+mSo4FneemgvzIdZN9SvmwFpqbOaMh+yf+6bP69wMmSZqRYRkg6uJcYX1/qnyB7zGA3sLas2DeBf4+Iu3LTbgUuzfVtHLMH22C2E9cIrB6NBP413ba5i+yOlJdExNr06/sxsqc8PbAH636C7BnJ3wTmU/aM2IhYkZqhbpDUmiZ/nOzOt3l/DlwtaXOK8W0R0V3KDZKmkSWqgyW9Oy3zv4FPkj2561FJDcDv2bmfw2y3+e6jZhVS9ljMm1NHs1ndcNOQmVnBuUZgZlZwrhGYmRWcE4GZWcE5EZiZFZwTgZlZwTkRmJkV3P8HZpd1eyKgxxQAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plt.plot(sim)\n",
"plt.title('Chi-squared p-value and sample size (deterministic)')\n",
"plt.xlabel('Sample Size')\n",
"plt.ylabel('p-value')\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Appendix\n",
"Calculation of chi-squared test with Yate's correction for continuity."
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [],
"source": [
"def chi_squared_test(table):\n",
" E = np.outer(table.sum(axis=1), table.sum(axis=0))/table.sum()\n",
" c = np.sum(np.square(np.abs(table - E)-0.5)/E)\n",
" dof = table.size - sum(table.shape) + table.ndim - 1\n",
" pvalue = 1-stats.chi2.cdf(c, dof)\n",
" return(c, pvalue)"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"(0.0026399847052020935, 0.9590220987022877, 1, array([[ 9.40298507, 20.59701493],\n",
" [11.59701493, 25.40298507]]))\n"
]
}
],
"source": [
"table = np.array([[10, 20],[11,26]])\n",
"print(stats.chi2_contingency(table))"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(0.0026399847052020935, 0.9590220987022877)"
]
},
"execution_count": 13,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"chi_squared_test(table)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The square of the two proportion Z-test statistc is the same as the test statistic associated with the Chi-squared test."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [],
"source": [
"def z_prop_test(table):\n",
" g1_sum = table[0,0]\n",
" g2_sum = table[1,0]\n",
" g1_n = table[0].sum()\n",
" g2_n = table[1].sum() \n",
" g1_mean = g1_sum/g1_n\n",
" g2_mean = g2_sum/g2_n\n",
" pool_mean = (g1_sum + g2_sum) / (g1_n + g2_n)\n",
" z = (g1_mean - g2_mean) / np.sqrt(pool_mean * (1-pool_mean) * ((1/g1_n)+(1/g2_n)))\n",
" pvalue = stats.norm.cdf(loc=0, scale=1, x=-abs(z))\n",
" return((z, pvalue*2))"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"(0.09997575214966518, 0.7518587343512207, 1, array([[ 9.40298507, 20.59701493],\n",
" [11.59701493, 25.40298507]]))\n"
]
}
],
"source": [
"print(stats.chi2_contingency(table, correction=False))"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.099975752149665"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"z_prop_test(table)[0]**2"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# References\n",
"* https://en.wikipedia.org/wiki/Welch%27s_t-test\n",
"* https://en.wikipedia.org/wiki/Pearson%27s_chi-squared_test\n",
"* https://en.wikipedia.org/wiki/Yates%27s_correction_for_continuity"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.7"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment