Last active
November 12, 2019 12:09
-
-
Save LuxXx/819f9b3441b8f5cc9e549a7eb4fbb740 to your computer and use it in GitHub Desktop.
Kernel Destiny Estimation using gauss kernel and silvermans rule and drawing new samples with the acceptance-rejection method
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| { | |
| "cells": [ | |
| { | |
| "cell_type": "code", | |
| "execution_count": 1, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "import numpy as np\n", | |
| "import matplotlib.pyplot as plt" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 2, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "X = np.append(np.random.normal(100, 100, 1000), np.random.normal(300, 40, 1000))" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 3, | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "text/plain": [ | |
| "array([166.20962345, 100.62571625, 51.98102988, ..., 316.92602048,\n", | |
| " 380.78154065, 292.82678533])" | |
| ] | |
| }, | |
| "execution_count": 3, | |
| "metadata": {}, | |
| "output_type": "execute_result" | |
| } | |
| ], | |
| "source": [ | |
| "X" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 4, | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAD8CAYAAAB3u9PLAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAFAFJREFUeJzt3XGMXeV55/HvrzaQKmlJAm6V2GbtCKeqWbVJNPJGSrXKhiQ4IcKtBMpk211ri2TtFrSpulLXbCR2y8pS2JWSalWiCBW0lE1ivLRRp4lbQpagqFJjPGwIwXa8mQArZkHFEYQ2WgE7zrN/3Jfm5vbOzJmZOzN3xt+PNOKc97zvOc+xL/Pc933PeZ2qQpKkn1rvACRJ48GEIEkCTAiSpMaEIEkCTAiSpMaEIEkCTAiSpMaEIEkCTAiSpGbregewFJdffnnt2rVrvcOQpA3j0Ucf/X5VbetSd0MlhF27djE9Pb3eYUjShpHkf3et65CRJAkwIUiSGhOCJAkwIUiSGhOCJAkwIUiSGhOCJAkwIUiSGhOCJAnYYG8qS7ow7Tr85b/bfvqT165jJJubPQRJEmAPQdKY6u8VaG3YQ5AkASYESVJjQpAkASYESVJjQpAkASYESVJjQpAkASYESVJjQpAkASYESVJjQpAkASYESVJjQpAkASYESVLj8teSxoZLXq8vewiSJKBjQkiyP8nZJDNJDg85fkmS+9rxE0l29R27pZWfTXJNX/nTSb6d5LEk06O4GUnS8i06ZJRkC3AH8AFgFjiZZKqqTvdVuxF4saquTDIJ3A58NMleYBK4Cngr8NUkb6+q863dP6mq74/wfiRJy9Slh7APmKmqJ6vqVeAocGCgzgHgnrZ9P3B1krTyo1X1SlU9Bcy080mSxkyXhLAdeKZvf7aVDa1TVXPAS8Bli7Qt4CtJHk1yaOmhS5JGqctTRhlSVh3rLNT2PVX1bJKfAx5M8p2q+vrfu3gvWRwCuOKKKzqEK0laji49hFlgZ9/+DuDZ+eok2QpcCrywUNuqeu2/zwNfZJ6hpKq6s6omqmpi27ZtHcKVJC1Hl4RwEtiTZHeSi+lNEk8N1JkCDrbt64GHqqpa+WR7Cmk3sAd4JMnrk/wMQJLXAx8Enlj57UiSlmvRIaOqmktyM/AAsAW4u6pOJbkNmK6qKeAu4N4kM/R6BpOt7akkx4DTwBxwU1WdT/LzwBd7885sBT5fVX+xCvcnSeqo05vKVXUcOD5Qdmvf9svADfO0PQIcGSh7EvjlpQYrSVo9vqksSQJMCJKkxoQgSQJMCJKkxuWvJW0o/UtkP/3Ja9cxks3HHoIkCTAhSJIah4wkrYn5hnr8V9LGhz0ESRJgQpAkNSYESRJgQpAkNSYESRJgQpAkNSYESRJgQpAkNSYESRJgQpAkNSYESRLgWkaSVpHrFG0s9hAkSYAJQZLUmBAkSYAJQZLUOKksaaScSN647CFIkgATgiSpMSFIkoCOCSHJ/iRnk8wkOTzk+CVJ7mvHTyTZ1XfsllZ+Nsk1A+22JPlmki+t9EYkSSuzaEJIsgW4A/gQsBf4WJK9A9VuBF6sqiuBTwO3t7Z7gUngKmA/8Jl2vtd8HDiz0puQJK1clx7CPmCmqp6sqleBo8CBgToHgHva9v3A1UnSyo9W1StV9RQw085Hkh3AtcAfrvw2JEkr1SUhbAee6dufbWVD61TVHPAScNkibX8f+F3gR0uOWpI0cl3eQ8iQsupYZ2h5ko8Az1fVo0neu+DFk0PAIYArrrhi8WglXTD633l4+pPXrmMkm0OXHsIssLNvfwfw7Hx1kmwFLgVeWKDte4DrkjxNbwjqfUn+27CLV9WdVTVRVRPbtm3rEK4kaTm6JISTwJ4ku5NcTG+SeGqgzhRwsG1fDzxUVdXKJ9tTSLuBPcAjVXVLVe2oql3tfA9V1W+M4H4kScu06JBRVc0luRl4ANgC3F1Vp5LcBkxX1RRwF3Bvkhl6PYPJ1vZUkmPAaWAOuKmqzq/SvUiSVqDTWkZVdRw4PlB2a9/2y8AN87Q9AhxZ4NwPAw93iUOStHp8U1mSBJgQJEmNy19LWnMukT2e7CFIkgATgiSpMSFIkgATgiSpMSFIkgATgiSpMSFIkgATgiSpMSFIkgDfVJa0CfkP5yyPPQRJEmBCkCQ1JgRJEmBCkCQ1TipL2hRcUnvl7CFIkgATgiSpcchI0rL4rP/mYw9BkgSYECRJjQlBkgSYECRJjZPKkja1wfcTnACfnz0ESRJgQpAkNZ0SQpL9Sc4mmUlyeMjxS5Lc146fSLKr79gtrfxskmta2euSPJLkW0lOJfm9Ud2QJGl5Fk0ISbYAdwAfAvYCH0uyd6DajcCLVXUl8Gng9tZ2LzAJXAXsBz7TzvcK8L6q+mXgHcD+JO8ezS1JkpajSw9hHzBTVU9W1avAUeDAQJ0DwD1t+37g6iRp5Uer6pWqegqYAfZVzw9b/YvaT63wXiRJK9AlIWwHnunbn21lQ+tU1RzwEnDZQm2TbEnyGPA88GBVnRh28SSHkkwnmT537lyHcCVJy9ElIWRI2eC3+fnqzNu2qs5X1TuAHcC+JP9w2MWr6s6qmqiqiW3btnUIV5K0HF0Swiyws29/B/DsfHWSbAUuBV7o0raqfgA8TG+OQZK0TrokhJPAniS7k1xMb5J4aqDOFHCwbV8PPFRV1con21NIu4E9wCNJtiV5I0CSnwbeD3xn5bcjSVquRd9Urqq5JDcDDwBbgLur6lSS24DpqpoC7gLuTTJDr2cw2dqeSnIMOA3MATdV1fkkbwHuaU8c/RRwrKq+tBo3KGn1+a+VbQ6dlq6oquPA8YGyW/u2XwZumKftEeDIQNnjwDuXGqwkafX4prIkCTAhSJIaE4IkCTAhSJIaE4IkCTAhSJIaE4IkCTAhSJIaE4IkCTAhSJIaE4IkCTAhSJIaE4IkCTAhSJIaE4IkCTAhSJIaE4IkCTAhSJIaE4IkCTAhSJIaE4IkCTAhSJKaresdgKSl2XX4y3+3/fQnr13HSLTZ2EOQJAEmBElSY0KQJAEmBElS46SyNAJLneh1YljjqFMPIcn+JGeTzCQ5POT4JUnua8dPJNnVd+yWVn42yTWtbGeSryU5k+RUko+P6oYkScuzaEJIsgW4A/gQsBf4WJK9A9VuBF6sqiuBTwO3t7Z7gUngKmA/8Jl2vjng31TVLwLvBm4ack5J0hrqMmS0D5ipqicBkhwFDgCn++ocAP5D274f+IMkaeVHq+oV4KkkM8C+qvor4DmAqvrbJGeA7QPnlFaFwzXScF2GjLYDz/Ttz7ayoXWqag54CbisS9s2vPRO4MSwiyc5lGQ6yfS5c+c6hCtJWo4uPYQMKauOdRZsm+QNwB8Dv11VfzPs4lV1J3AnwMTExOB1pZHp0nPor7Ma5xzVdaXl6NJDmAV29u3vAJ6dr06SrcClwAsLtU1yEb1k8Lmq+pPlBC9JGp0uCeEksCfJ7iQX05sknhqoMwUcbNvXAw9VVbXyyfYU0m5gD/BIm1+4CzhTVZ8axY1IklZm0SGjqppLcjPwALAFuLuqTiW5DZiuqil6v9zvbZPGL9BLGrR6x+hNFs8BN1XV+SS/Avwz4NtJHmuX+ndVdXzUNyitlMM1ulB0ejGt/aI+PlB2a9/2y8AN87Q9AhwZKPtLhs8vSJLWiW8qS2NkPR+Jne/aPqZ74XAtI0kSYEKQJDUOGUmbhEM7S+ef2U+yhyBJAuwh6AKxkjeMpQuFPQRJEmBCkCQ1DhlJY8phK601ewiSJMAegjaZcXiMcC2/2Xe9Vpe3kCV7CJIkwIQgSWocMpLUmUNMm5s9BEkSYA9Bm5jfZqWlsYcgSQJMCJKkxoQgSQJMCJKkxkllrYuVvFE8yreRN+vE82a9r9U0Dm+5rzd7CJIkwIQgSWocMtLY6tKFd2hk6fwz03zsIUiSABOCJKkxIUiSgI4JIcn+JGeTzCQ5POT4JUnua8dPJNnVd+yWVn42yTV95XcneT7JE6O4EUnSyiw6qZxkC3AH8AFgFjiZZKqqTvdVuxF4saquTDIJ3A58NMleYBK4Cngr8NUkb6+q88B/Bf4A+KNR3pA0rjbSZO5GilWj06WHsA+Yqaonq+pV4ChwYKDOAeCetn0/cHWStPKjVfVKVT0FzLTzUVVfB14YwT1Ikkagy2On24Fn+vZngX80X52qmkvyEnBZK//GQNvty45Wm5JviErjoUsPIUPKqmOdLm0XvnhyKMl0kulz584tpakkaQm6JIRZYGff/g7g2fnqJNkKXEpvOKhL2wVV1Z1VNVFVE9u2bVtKU0nSEnQZMjoJ7EmyG/g/9CaJ/+lAnSngIPBXwPXAQ1VVSaaAzyf5FL1J5T3AI6MKXpJWw4U6jLloD6Gq5oCbgQeAM8CxqjqV5LYk17VqdwGXJZkBfgc43NqeAo4Bp4G/AG5qTxiR5Av0EsgvJJlNcuNob02StBSpWtKQ/rqamJio6enp9Q5DI+BjjdooNnoPIcmjVTXRpa5vKkuSABOCJKlx+WtJ6mizTzbbQ5AkAfYQtAo2+7coabOyhyBJAkwIkqTGISMtyXzvD3QZGvLdA2m82UOQJAEmBElSY0KQJAEmBElS46SyRsIJY21WF9Jn2x6CJAkwIUiSGhOCJAkwIUiSGieVN5mVLCy3GovSXUgTctJGZw9BkgTYQ7jgzfcN3iWspQuPPQRJEmBCkCQ1DhlpUU4MSxcGewiSJMAewqaw1G/wfuOXVs9GfiDDHoIkCTAhSJIah4zWyVL/beLB+hutKyptNht5aGg+nXoISfYnOZtkJsnhIccvSXJfO34iya6+Y7e08rNJrul6TknS2kpVLVwh2QL8L+ADwCxwEvhYVZ3uq/NbwC9V1b9MMgn8WlV9NMle4AvAPuCtwFeBt7dmC55zmImJiZqenl76XbJ+2Xy+63aZ2F1qfUkbQ5ffQaP6nZXk0aqa6FK3Sw9hHzBTVU9W1avAUeDAQJ0DwD1t+37g6iRp5Uer6pWqegqYaefrck5J0hrqkhC2A8/07c+2sqF1qmoOeAm4bIG2Xc4pSVpDXSaVM6RscJxpvjrzlQ9LREPHrpIcAg613R8mOTtQ5XLg+8Pazie3L6X2SFwOfH+p113jOJf857hONkKcxjg6GyHOVf8dtMLfBf+ga8UuCWEW2Nm3vwN4dp46s0m2ApcCLyzSdrFzAlBVdwJ3zhdckumu42PrxRhHZyPEaYyjsxHi3AgxdtVlyOgksCfJ7iQXA5PA1ECdKeBg274eeKh6s9VTwGR7Cmk3sAd4pOM5JUlraNEeQlXNJbkZeADYAtxdVaeS3AZMV9UUcBdwb5IZej2Dydb2VJJjwGlgDripqs4DDDvn6G9PktRVpxfTquo4cHyg7Na+7ZeBG+ZpewQ40uWcyzTvcNIYMcbR2QhxGuPobIQ4N0KMnSz6HoIk6cLgWkaSJGADJYQk/zHJ40keS/KVJG9t5UnyX9oSGI8neVdfm4NJvtt+Ds5/9pHG+Z+TfKfF8sUkb+w7NhbLeCS5IcmpJD9KMjFwbCxiHBLz2Cx1kuTuJM8neaKv7M1JHmyftQeTvKmVz/v5XOUYdyb5WpIz7e/64+MWZ5LXJXkkybdajL/Xyne3JXC+m96SOBe38nmXyFmDWLck+WaSL41rjCNRVRviB/jZvu1/DXy2bX8Y+HN67zy8GzjRyt8MPNn++6a2/aY1iPODwNa2fTtwe9veC3wLuATYDXyP3oT6lrb9NuDiVmfvKsf4i8AvAA8DE33lYxPjQLzrev0h8fxj4F3AE31l/wk43LYP9/29D/18rkGMbwHe1bZ/ht5SMXvHKc52rTe07YuAE+3ax4DJVv5Z4F+17d/q+/9+ErhvDf/Ofwf4PPCltj92MY7iZ8P0EKrqb/p2X8+PX2Q7APxR9XwDeGOStwDXAA9W1QtV9SLwILB/DeL8SvXe1gb4Br13LF6LcyyW8aiqM1U1+ILfWMU4YL2v/xOq6uv0nqbr1798yz3Ar/aVD/t8rnaMz1XV/2zbfwucobcawNjE2a71w7Z7Ufsp4H30lsAZFuOwJXJWVZIdwLXAH7b9jFuMo7JhEgJAkiNJngF+HXjtKadxXh7jN+l962KBeMYhzteMa4zrff0ufr6qnoPeL2Pg51r5usfehi3eSe8b+FjF2YZiHgOep/el7XvAD/q+VPXHMd8SOavt94HfBX7U9i8bwxhHYqwSQpKvJnliyM8BgKr6RFXtBD4H3PxasyGnWmjZjFWPs9X5BL13Lz63HnF2iXFYs7WMcQnW+/orsa6xJ3kD8MfAbw/0sv9e1SFlqx5nVZ2vqnfQ60nvozecOV8cax5jko8Az1fVo/3FC8SxkT+r4/UP5FTV+ztW/TzwZeDfM//yGLPAewfKH15xkCweZ3oT2B8Brq42mLhAnCxQvmoxzmNNY1yCLsunrLe/TvKWqnquDbU838rXLfYkF9FLBp+rqj8Z1zgBquoHSR6mN4fwxiRb2zfs/jjmWyJnNb0HuC7Jh4HXAT9Lr8cwTjGOzFj1EBaSZE/f7nXAd9r2FPDP21MS7wZeal3hB4APJnlTe5Lig61stePcD/xb4Lqq+r99hzbCMh7jGuN6X7+L/uVbDgJ/2lc+7PO5qtq49V3Amar61DjGmWRb2lN4SX4aeD+9uY6v0VsCZ1iMw5bIWTVVdUtV7aiqXfQ+dw9V1a+PU4wjtd6z2l1/6H3TeQJ4HPgzYHv9+EmFO+iNPX6bn3xq5jfpTYzOAP9ijeKcoTeG+Fj7+WzfsU+0OM8CH+or/zC9p0C+B3xiDWL8NXrfZF4B/hp4YNxiHBLzul5/IJYvAM8B/6/9Od5Ib5z4fwDfbf9982Kfz1WO8VfoDVU83vdZ/PA4xQn8EvDNFuMTwK2t/G30vojMAP8duKSVv67tz7Tjb1vjv/f38uOnjMYyxpX++KayJAnYQENGkqTVZUKQJAEmBElSY0KQJAEmBElSY0KQJAEmBElSY0KQJAHw/wF59rRaLA80kgAAAABJRU5ErkJggg==\n", | |
| "text/plain": [ | |
| "<Figure size 432x288 with 1 Axes>" | |
| ] | |
| }, | |
| "metadata": { | |
| "needs_background": "light" | |
| }, | |
| "output_type": "display_data" | |
| } | |
| ], | |
| "source": [ | |
| "p = plt.hist(X, bins=100, density=True)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 5, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "def K(x, mu=0, std=1):\n", | |
| " return 1 / (np.sqrt(2*np.pi*std)*std) * np.exp(-((x-mu)*(x-mu)) / (2*std*std))" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 6, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "# Silverman rule\n", | |
| "h = 1.06 * np.sqrt(X.var()) * np.power(len(X), -1/5)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 7, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "# estimated density\n", | |
| "def f_n(x):\n", | |
| " s = 0\n", | |
| " for xi in X:\n", | |
| " s += K((x-xi)/h)\n", | |
| " return s / (len(X)*h)" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 8, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "xvals = np.arange(np.floor(X.min()), np.ceil(X.max()))" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 9, | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAD8CAYAAAB3u9PLAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XuclHXd//HXZ2YPgAgIrsoZhAUFRUFC0BISD5gmappgvzS1zNLKuq0bzSwtS7M7u7u1zFLzlICaySnxFJKm4HKUo6wcBFFBUQ4qe5j5/P64LnBcdtlZdmavmd338/GYx17znevwHpjdz3yvw/cyd0dERCQWdQAREckNKggiIgKoIIiISEgFQUREABUEEREJqSCIiAiggiAiIiEVBBERAVQQREQkVBB1gIY48MADvVevXlHHEBHJG/PmzXvX3UvSmTevCkKvXr0oKyuLOoaISN4ws3XpzqtdRiIiAqggiIhISAVBREQAFQQREQmpIIiICKCCICIiIRUEEREBVBBERCSkgiAiIkCeXaksIi1TrwnTd0+vvfn0CJM0b+ohiIgIoB6CiOSo1F6BNA31EEREBFBBEBGRkAqCiIgAKggiIhJSQRAREUAFQUREQioIIiICqCCIiEhIBUFERAAVBBERCakgiIgIoIIgIiIhFQQREQFUEEREJKThr0UkZ2jI62iphyAiIkCaBcHMxpjZSjMrN7MJtbxebGaTwtfnmFmvlNeuCdtXmtmpKe1rzexVM1toZmWZeDMiIrLv6t1lZGZx4A7gZGAD8IqZTXH3ZSmzXQq87+59zWwccAtwvpkNAMYBA4EuwDNm1s/dE+Fyn3f3dzP4fkREZB+l00MYBpS7+2p3rwQmAmNrzDMWuC+cfhQYbWYWtk909wp3XwOUh+sTEZEck05B6AqsT3m+IWyrdR53rwa2Ap3qWdaBp8xsnpld1vDoIiKSSemcZWS1tHma8+xt2ePdfaOZHQQ8bWYr3H32HhsPisVlAD169EgjroiI7It0eggbgO4pz7sBG+uax8wKgPbAlr0t6+67fm4CHqeOXUnufpe7D3X3oSUlJWnEFRGRfZFOQXgFKDWz3mZWRHCQeEqNeaYAF4XT5wLPubuH7ePCs5B6A6XAXDPbz8z2BzCz/YBTgCWNfzsiIrKv6t1l5O7VZnYlMBOIA/e4+1IzuxEoc/cpwN3AA2ZWTtAzGBcuu9TMJgPLgGrgCndPmNnBwOPBcWcKgL+5+5NZeH8iIpKmtK5UdvcZwIwabdenTO8Ezqtj2ZuAm2q0rQaOamhYERHJHl2pLCIigAqCiIiEVBBERARQQRARkZCGvxaRvJI6RPbam0+PMEnzox6CiIgAKggiIhLSLiMRaRJ17erRXdJyh3oIIiICqCCIiEhIBUFERAAVBBERCakgiIgIoIIgIiIhFQQREQFUEEREJKSCICIigAqCiIiEVBBERATQWEYikkUNGafoYLbwnYLHGR1fQGsqWJDsy+3VZzHP+2cxoaRSD0FEIjfMlvNk8QTOiz9PWbIf0xPDGRBbxyNFN/L1uAa/ayrqIYhIpAbaWu4pupW3vSPnVN3AGu8MwC+qv8KthX/iusKH+JBWPJwYHXHS5k8FQUSi8/H7/LnoN3xAW75SeS3v0PGTl2jF96quZH8+5mcF97MgWcoK7xFh2OZPu4xEJDozfkQJW7m88qpPFYNdEsT5QdW3+ID9+J/CPxIjGUHIlkM9BBHJqLQPJK+eBa9O5o7EOSzxQ+uc7T3a8/Oqr3J70f9xbvx5Jic+n5mgsgf1EESk6SUTMPM6aN+DP1afWe/s05LDKUv24+qCRyimsgkCtkwqCCLS9JY8Bu+8Cif9lAqK0ljAuLXqfA6yDzg3Pjvr8VoqFQQRaWIOL9wGJYfDwHPSXmqOH8aCZF8ui08jTiKL+VqutAqCmY0xs5VmVm5mE2p5vdjMJoWvzzGzXimvXRO2rzSzU2ssFzezBWY2rbFvRETyw+djC2HTMvjsVRBryHdS487qL9IztonRsflZy9eS1fu/YWZx4A7gNGAAMN7MBtSY7VLgfXfvC9wG3BIuOwAYBwwExgB/CNe3y/eA5Y19EyKSP75VMAXad4cjvtTgZZ9JDuFtP4Bx8X9lIZmkU56HAeXuvtrdK4GJwNga84wF7gunHwVGm5mF7RPdvcLd1wDl4fows27A6cBfGv82RCQfHG7rGBZbCcdeDvHCBi+fIM6kxChGxRbRhXezkLBlS6cgdAXWpzzfELbVOo+7VwNbgU71LPs74EegE4tFWorx8eeo8EI4+oJ9Xsfk6lEAOricBelch2C1tHma89TabmZnAJvcfZ6Zjdrrxs0uAy4D6NFDVymK5KvW7OSs+AtMTx7LOW32vAgtXW9Swiveny/GX6LXhGns+jOz9ubTM5S05Uqnh7AB6J7yvBuwsa55zKwAaA9s2cuyxwNnmtlagl1QJ5rZg7Vt3N3vcveh7j60pKQkjbgikovOiL9MO/uYv1Wf2Oh1TU2MoDT2Jv1tff0zS9rSKQivAKVm1tvMiggOEk+pMc8U4KJw+lzgOXf3sH1ceBZSb6AUmOvu17h7N3fvFa7vOXf/fxl4PyKSo86Nz+b1ZGfKMjCc9T8Tw0i48cX4SxlIJrvUWxDCYwJXAjMJzgia7O5LzexGM9t1ieHdQCczKwd+AEwIl10KTAaWAU8CV7i7TiAWaWE68x7HxlbweOKz1L4nuWHeoz0vJo/gjNjL7LkHW/ZVWmMZufsMYEaNtutTpncC59Wx7E3ATXtZ9yxgVjo5RCQ/nRF+k5+aHJGxdU5PDueWwj8zwNaxzHtlbL0tma5UFpGsGxv/DwuTfVjnh2Rsnc8lBgNwYmxBxtbZ0qkgiEhW9bE3OSK2limJ4zK63s10YGGyD6PjKgiZouGvRSSrzoy/RNKNqYnhu9sacq/lvXk2MZjvFzzGgWzNyPpaOvUQRCSrTovNYU7ycDZzQMbX/WxyCDFzPq9eQkaoIIhI1vS2t+gXe5OZyaFZWf8y78lG78hoHUfICBUEEcmak2NlADyVyE5BAOP5xFEcF1sKieosbaPlUEEQkaw5JT6PV5O92MiBWdvGi8kjaGcfwVsLs7aNlkIFQUSyooQPGGKrstg7CPwnOTCYWK0hsRtLBUFEsmJ0fD4xc57K0vGDXbbQjiXJXvD6rKxupyVQQRCRrDglVsa65EGs9O71z9xILySPgPVzoPLDrG+rOVNBEJGM24+POT62JOwdNH7sovq8mDwCklWwToPdNYYKgohk3GdjSyi2ap5NDmmS7c1NHgbxYh1HaCRdqSwiGTcytpBt3pqyZL8m2V4FRdDjWFj9PLDnldC6eU561EMQkQxzRsUX8WLyCKqb8jtnz8/CO0vg4w+abpvNjAqCiGRUP9tAF9vCrOTRTbvhnscBHhxcln2igiAiGTUqFlwg9nxiUNNuuNtQiBXCuhebdrvNiAqCiGTUqNgiViS78zadmnbDha2h6xCdadQIKggikjkV2xkaW8ms5FFNvuleE6bzhzUHU7V+Hq2oaPLtNwcqCCKSOaufp8gSPB9BQQCYm+xPoSUYHCuPZPv5TgVBRDKn/Bl2eCvKkv0j2fy8ZH+SbgyzFZFsP9/pOgQR2Sep5/qvvfl0cIfyZ3gxeQRVEf1p2U4blnlPhsVWQCKSCHlNPQQRyYzNK2Hr+kiOH6R6JdmfIbFVFKL7IzSUCoKIZEb50wA8n4i2IMxJHk5rq+QIWxNpjnykgiAimbHqaSg5LKs3w0nHvHC4jCGx1yLNkY9UEESk8Sp2wBsvQd+Tok7CZjqwPlnCkNiqqKPkHR1UFpHGWzMbEpVQejL8K/p7Esz3Uo6NLQccsD0PgEut1EMQkcYrfxoK94MeI6JOAsC8ZCmH2Pt04b2oo+QVFQQRaSSHVc/AoSOhoDjqMADMT5YCaLdRA6VVEMxsjJmtNLNyM5tQy+vFZjYpfH2OmfVKee2asH2lmZ0atrUys7lmtsjMlprZDZl6QyLStPrYRtj6Rk4cP9hlhffgYy9SQWigeguCmcWBO4DTgAHAeDMbUGO2S4H33b0vcBtwS7jsAGAcMBAYA/whXF8FcKK7HwUcDYwxs+GZeUsi0pRGxRYFEzlUEKopYLEfqjONGiidHsIwoNzdV7t7JTARGFtjnrHAfeH0o8BoM7OwfaK7V7j7GqAcGOaBHeH8heHDG/leRCQCI2OL4MB+cEDPqKN8yrxkPwbaOoqpjDpK3kinIHQF1qc83xC21TqPu1cDW4FOe1vWzOJmthDYBDzt7rXe1cLMLjOzMjMr27x5cxpxRaSptKKCY2MroO/JUUfZw/xkKYWW4EhbHXWUvJFOQbBa2mp+m69rnjqXdfeEux8NdAOGmdkRtW3c3e9y96HuPrSkpCSNuCLSVEbEllFsVVCaO7uLdlmQ7AvAMTqOkLZ0CsIGoHvK827AxrrmMbMCoD2wJZ1l3f0DYBbBMQYRySMjY4v4yIuhx3FRR9nDe7RnbfJgHVhugHQKwitAqZn1NrMigoPEU2rMMwW4KJw+F3jO3T1sHxeehdQbKAXmmlmJmXUAMLPWwEmAxqsVyTOjYot4KTkACltFHaVW87w0LAg6RJmOeq9UdvdqM7sSmAnEgXvcfamZ3QiUufsU4G7gATMrJ+gZjAuXXWpmk4FlQDVwhbsnzKwzcF94xlEMmOzu07LxBkUkO3ra2/SKvcPdVadxacqVwLlkQbKUL8VfoJttZoMfFHWcnJfW0BXuPgOYUaPt+pTpncB5dSx7E3BTjbbFwOCGhhWR3LHrdNOo7o6Wjl0XqB1jr6kgpEFXKovIPhkZW8Tq5CG84QdHHaVOK707O7yVjiOkSQVBRBquaicjYstyuncAkCDO4uShusdymlQQRKTh1r1Aa6vM+YIAwcinh9sbtKIi6ig5TwVBRBqu/Fl2emFwhlGOW5DsG16gpjuo1UcFQUQabtVTzEkeTgVFUSep1wKNfJo2FQQRaZh3V8F75TyTHBJ1krRsoR1rdIFaWlQQRKRhVv4TgGcT+VEQABZ4aXhgWReo7Y0Kgog0zMp/wsFHspEDo06StvnJUg6yD+hm70YdJaepIIhI+j7aAutfhv6nRZ2kQXYNdDfYtNtob1QQRCR9q54CT+ZdQVjhPfjIi3UcoR4qCCKSvpUzoO0h0PnoqJM0SII4i/1QBqsg7JUKgoikp7oCyp+F/mMgln9/OhYk+zLQ1kHVx1FHyVn5978qItFY82+o3AH98mt30S677qDGW4uijpKzVBBEJD3L/gFF+8Oho6JOsk92XaDGhleiDZLDVBBEpH6JKlgxLTiYnKM3w6nPu7TnjWQJrJ8bdZScpYIgIvVbMxs+fh8GnhV1kkaZ76WwoSzqGDlLBUFE6rdrd1Gf0VEnaZQFyVLYvhG2bog6Sk5SQRCRvUtUwfJpwdlFebq7aJddd1DTbqPaqSCIyN6tmQ0fb4EB+b27CIIL1Chopd1GdVBBEJG9WzwJWrWHvidFnaTRqiiALoNhg3oItVFBEJG6VWyH5VNh4Dl5v7tot26fCa5FqNYd1GpSQRCRui17Aqo+gqMviDpJ5nT7DCQq4a3FUSfJOQVRBxCRhus1Yfru6bU3n569DS18GDr2Cf6INhe73suGudC9Gb2vDFAPQURq9/5aWPcCHDUezKJOkzntOkP77rpiuRYqCCJSuwUPAQZHnR91kszr9hlYr4JQkwqCiOypugLm3Qv9ToUOPaJOk3ndh8G2DbBtY9RJcooKgojsadkU+HAzDPtG1EmyY/dxBPUSUumgskgGNPQgb5MdFN5Xc+8KDiYfemLUSbLjkEEQLw6uWB4wNuo0OSOtHoKZjTGzlWZWbmYTanm92Mwmha/PMbNeKa9dE7avNLNTw7buZvYvM1tuZkvN7HuZekMi0khvzAnOwBn2jby8EU5aCoqg81G6YrmGev+3zSwO3AGcBgwAxpvZgBqzXQq87+59gduAW8JlBwDjgIHAGOAP4fqqgf9y98OB4cAVtaxTRKLw799Am04w5MKok2RX92GwcQFUV0adJGekU/6HAeXuvtrdK4GJQM0+1ljgvnD6UWC0mVnYPtHdK9x9DVAODHP3t9x9PoC7bweWA10b/3ZE6tdrwvTdD6lh40JY9RQM/zYU7Rd1muzq9hlIVMA7r0adJGekUxC6AutTnm9gzz/eu+dx92pgK9ApnWXD3UuDgTm1bdzMLjOzMjMr27x5cxpxRWSfPf9rKG7ffA8mp9p1YFmnn+6WzkHl2q5I8TTn2euyZtYWeAy4yt231bZxd78LuAtg6NChNbcrkjHpHOhtaK+ioevM1Hb3yZrZsHI6nPiTYDC75q59V2jXDda/DMMvjzpNTkinh7AB6J7yvBtQ8+Td3fOYWQHQHtiyt2XNrJCgGDzk7n/fl/AikiHJBMy8Ftr3gBFXRJ2m6fQcAeteAtd3TUivILwClJpZbzMrIjhIPKXGPFOAi8Lpc4Hn3N3D9nHhWUi9gVJgbnh84W5gubv/NhNvREQaoeweePtVOPlnUNg66jRNp+dxsONt2LI66iQ5od5dRu5ebWZXAjOBOHCPuy81sxuBMnefQvDH/QEzKyfoGYwLl11qZpOBZQRnFl3h7gkz+yzwVeBVM1sYbupad5+R6Tco0ljN/uDzu+Xw1E+C22MOPCfqNE2rx3HBz3X/gU59os2SA9K6MC38Qz2jRtv1KdM7gfPqWPYm4KYabS9Q+/EFEWlK1ZXw+DehoBjG3t68BrFLR0n/4BTbN16CIV+NOk3kdKWySAPFSdDZ3uMgPuBge5/29iHFVFFMJYUk4IXXIF4EsUIuiK+gkgI+9FZQXgxFbaGoLd3tHT7yVuygNRUUsuv7UZNewewO074Pb5bBeX+l1y8XAAv22HbOX1XdGGbQYwSsezHqJDlBBUGkHgfxPsfFljIitoyBsbWU2psUW1XdCzzzyeQvC1PaH/zf3ZP/Lv6kudpjfEQrtvj+bKIDm7wDm70DvLQu+AZbchh1n7S3j9zhmZ/Cwgdh5H/DwLOBZr5rrC49j4MV02Drm8GZRy2YCoJIbba/DUv+Dq8+wtxW8wH4wPdjcfJQ7kuewuvehbe9I+/4AXzg+7GTIioopJoCVt14MiSrIFHFsJ8/SbFV0oYKZn5rCFTugMoP+a8HX6SN7aQtO8OfH9PRtnMQH3C4vcGo2CKYOXN3nEXFbZifLKUs2Z/5XsqCZF92UlxX+r2r2gkzroYFD8DQS2HkHqPRtCw9w+MIb7wER54bbZaIqSCIpBhkr3NpwT/htrmQrIZDBnFL1ThmJ49kmffE0zkxr6jN7slNHPDJlTc9jt3d/lgyXs9KnLXXHQubV8C7K5n+xDSGxl7j84WLAKjyOIu8D3OShwW7oroPh+K29e/eWfsCzPgRbFoKn7saTryu5R03SNFrwnTiJFhY3Jr91/1HBSHqACI54a1F3Fd4MyPji9nurWHYN+GYr0FJP/4YyVlGBm1Lgkfvz3HtY50BaMcOhsRWcWxsBcfGlnNZfDo8OAUsDl2OZkJBZ5Yme/GmHwjv9gNPwtYN8PbiYEjrjfODi7HGT4L+YyJ4X7knQZz5yVJGrvtP1FEip4IgLUKdp45ufxueug5efYRBsbb8qmo8DyVGs2NWG5i1CljVpDnrs422zEoOZlZyMABt2MmybxwAa1+EdS9ycfxJiguqg5lv/9mnFz74SDjtVhj8leY/TlEDzUkexsjNk+HD92C/TlHHiYwKgrRIRpIL4s/B7d+C6p3wuas54enD2E6b+hfOIR/RCvqcGDyAQRMep6e9Q1d7l3vH9Q92B+3fGQ4shbYHRZw2d81NHhZMvPESHH5GtGEipIIgLU5XNnNb0R8YFlsJXU6AM34Hnfqw/en8P8umgiJe8+685t1hUDM7RTSLFnsfKGgFa/+tgiDSUpwRe4lfFt6N4Vxd9U1+c+EtOXtQtdlfIZ1DKimEHsNh9fNRR4lUM70dksinFVHFrwr+zO1F/0e5d+G0yl/xaGJkzhYDicCho2Dzctj+TtRJIqMegjQrtZ52ue0tJhX9nMGxcm6vHstt1eeSoL7TPjOToSmks726TkdVLyRF75HBzzWzYVCtI/E0e+ohSPP2xhy4ayT9bD3frLyK31Sfn9ViIHms81HQqgOsnhV1ksiohyDN1/Kp8Ogl0K4rZ1f+IDjQKlKXWBx6fy4oCO4tcneiCoI0S+fFZ8Hkv0DXY+CCybx240tRR2oWmv0upkNHBV8ktqxukcNha5eRNDvfiE/j1sK7gl/uC5+ANh2jjiT5oveo4GcL3W2kHoI0H+78d8HDfKtgKtMSw/n+0q9Rdf2sqFNJPunUJxjaY/Us+MylUadpcuohSPOQTMDU7/Ktgqk8WD2a71ZdSZW+70hDmQU9yzWzg89UC6OCIPmvugIeuQjm38/vq8/iuupLSOqjLfuq72jY+QFsKIs6SZPTb43kt4rt8NB5wYHAU3/Fb6u/jO7OKo3S58Rg9NhVM+uft5lRQZD89eF7cN+ZwRj/Z90JI74ddSJpDlp3CIaxeO2pqJM0Oe1klcg05l69wyfczwNFN9PDNlE8/kE47AuNztDcNOf3lg2f+jyecUpwi9FtG6FdlwhTNS31ECT/vLuKR4tv4GDbwoWVE/a5GIjUqd+pwc9VLauXoIIg+WXjQrhnDK2oZHzlT5jjh0edSJqjksOgfY8Wt9tIu4wkp32qG//NdvDweGjdgfMqr2GNd95jHqmf/r3SYAalJ8OiicFZbAXFUSdqEuohSF44JfYKPPglaN8VLpm5uxiIZE2/U6HqQ1j3YtRJmowKguS8c+PP88fC38EhR8DF/wyKgki29T4BCtvAipbTo1JBkNzlzjfjU/lN4Z/4T3IgXDhF4xJJ0ylsHew2Wj4Vksmo0zSJtAqCmY0xs5VmVm5mE2p5vdjMJoWvzzGzXimvXRO2rzSzU1Pa7zGzTWa2JBNvRJqZZBKeuo5rCh9mamI4l1b9EIrbRp1KWprDz4Qd78D6OVEnaRL1HlQ2szhwB3AysAF4xcymuPuylNkuBd53975mNg64BTjfzAYA44CBQBfgGTPr5+4J4K/A7cD9mXxD0gwkquCJK2DxJO6tPpUbq7+KN4PObD4dzM2nrFnV71SIF8OyJ6DniKjTZF06v2XDgHJ3X+3ulcBEYGyNecYC94XTjwKjzczC9onuXuHua4DycH24+2xgSwbegzQnlR8GZxItngQnXscN1Rc2i2Igeap4/2Bso+VTg5vmNHPp/KZ1BdanPN8QttU6j7tXA1uBTmkuK0KvCdM5esJEFvziBHj9Wfji/8IJP0TjEknkDj8Ttm2AN+dHnSTr0ikItf1G1iyVdc2TzrJ737jZZWZWZmZlmzdvbsiikke62SYeLbqBAbYOvvwAHPO1qCOJBPqfBrFCWPJY1EmyLp2CsAFIvRltN2BjXfOYWQHQnmB3UDrL7pW73+XuQ919aElJSUMWlTxxtJXzj6LrOdC28tXKCXD4GVFHEvlE6w7BsYQlj0KiOuo0WZVOQXgFKDWz3mZWRHCQeEqNeaYAF4XT5wLPubuH7ePCs5B6A6XA3MxEl+ZgTGwuE4t+zofeinMqb2CuhqKQHNBrwvTdDwAGnR+cbdTMb61Zb0EIjwlcCcwElgOT3X2pmd1oZmeGs90NdDKzcuAHwIRw2aXAZGAZ8CRwRXiGEWb2MPAS0N/MNphZy7tfXUvmzmXxqdxZ9DuWei/OrryR1d5yRpWUPNPvVGjVARZPjDpJVpnn0ZHzoUOHellZy7uLUbOTqIIZP4R59zItMZz/qrqcCoqiTiVSq91Ds0/7Pix8GH64Kjj7KE+Y2Tx3H5rOvDqfT5rWh+/BA2fDvHu5o/pMvlN1pYqB5IdB46D6Y1hWc49586GCIE3nrcVw1yhYPxfOupNbq8fpGgPJH92HQae+MP+++ufNU/ptlKbx6qNw9yngCbjkSTh6fNSJRBrGjJ+/PRzWz+G0a/4QdZqsUEGQ7Eom4Onr4bFLofNRcNks6Dok6lQi++TRxAns9EL+X/yZqKNkhW6QIxm361S9jmxj/sDJ8PpzMPQSGHMLFOh4geSvrbRlamIEZ8VfgJ3boFW7qCNllHoIkhXDbDkziq+BtS8Gw1CccZuKgTQLDyROZj+rCO6m1syoIEhmJZN8O/4EDxf9go+8GL7xrIahkGZlsfdhYbIPvHxHs7tyWbuMpEHqGhZ57c2nB6eUPn4ZPyp8hqmJ4VxT9XWWHHJkvcuK5Js/Vp/Jn96/DZb9A448N+o4GaMegmTGqmfgjyNgzWyuq7qY71R9hx20iTqVSFY8lTwGDuwHL/yuWQ2LrYIgjdKKCm4ouBce+hK07ghff5YHEyejYaulOXNicPxV8M6rUN58zjhSQZB9doStZnrRtVxU8DQM/3ZwSmnnQVHHEmkag74M7bvDrF81m16CCoI0WAHVXBH/B48X/ZQ2VsEFldfCmF9BYauoo4k0nXghjJoAb84LbrHZDOigsjTIQFvDrwvvYmBsHVMSI7iu6mK20VYHjKXZ2utn+6jx8NId8OwNcNjpQZHIY+ohSHqqPoanf8oTRT+hxLbyzcqr+G7Vd9hG26iTiUQnFoeTfgZbVsO8v0YcpvHUQ5D6rX0BpnwXtrzOY4lR3FR9gQqByC6lp0Cvz8Fzv4ABY6HtQVEn2mfqIUjdtr8Nf78M/no6JKvhq//gv6svUzEQSWUGp/8Wqj6CJydEnaZRVBBkT9WV8OLv4f+OgaWPw+euhm+/BH0+H3UykdxU0g9O+CEseQxeeyrqNPtMu4yamdQDYLvv9JT2stM4JVbGjwom0Te2EfqNgVN/CZ36ZCSPSLN2/FVBQZj6Xbj8BdjvwKgTNZh6CBJY+yKPFf2Mu4puw3AYPwkumNSoYiDSohQUwTl/ho+2BLtak8moEzWYeggtmTvnX/s/XFHwD06Iv0oX68iPqr7BY4kTSNybAIJv9w3taYi0WJ0HBdfkTP8BvPA/wW6kPKKC0BIlk/Dak/DCb5lU/AqbvR2/rBrPfYlTdX9jkcYaegms+09w1lHHQ+GIL0WdKG0qCC3Jjs2w8EGYdx+8vwY69OS6qot5JDFShUAkU8xg7O0sJmwMAAAJyUlEQVSw7U14/HLYrwR6nxB1qrSoIDR31RVQ/iy8OhmWT4NkFfQ8Hk68DgacxYM/nlnvKnRgWKSBClvD+IfhntPg4fEw7iE4dFTUqeqlgtAMFVPJ8bEl8PcpsHIGVGyD1gfAZ74OQy+Gkv5RRxRp/lofABf+Ax44Gx46D750Nww4M+pUe6WC0Az0mjCNvvYmI2OLub9wEcNiK2hlVfBa++ADOPBs6D1y9zgr+sYvkj17nPr9tenwty/D5AvhhKth1DXBkBc5SAUhHyWq4O3F8MYcWP8yc4qf52D7AIDyZBceSpzE88lB3H/d1bqPsUjU2nSEC6fAjB/C7Fth3Utw5u9z8pRuFYRcl0zAe68HBeCtRfDm/GC43eqPg9fb9+Dl5ABeTh7O84mj2EjKxTAqBiK5oagNnHUH9DwuGN7iDyPgs1fBiCuhVbuo0+2mghChmrtuOrCdvvYmj36pI7z9alAE3lkajJECVHgBK7wH85IjuWTc+dD9WGjfle9pF5BIk9unUQEGfwX6joaZ18Lzt8CcP8Gxl8MxF0G7LllKmr60CoKZjQH+F4gDf3H3m2u8XgzcDxwDvAec7+5rw9euAS4FEsB33X1mOutsthJVsHU9bFnNJfF/0tfepE9sI31sIwfatmCeaUBxOzjkSBhyIRwyCDoP4ojfraYq/C+75AhdLCaSl/Y/BM69J+gdzL4Vnr8ZZv8a+pwYXLPQZzTsf3Ak0czrufWbmcWB14CTgQ3AK8B4d1+WMs+3gUHufrmZjQPOdvfzzWwA8DAwDOgCPAP0Cxfb6zprM3ToUC8rK2v4u6RxY/ykzR12boUdm2DHO8F5yO+v45FnXqB7bDPdbDPdYlvAP7mk/X1vS7l34fVkF8q9azDtXdjgJTixT2XVwWCR5uFTf4Peex0W/g0WPRz8zQA4aCAPbzyIV/1QliZ78sRPL97nXUtmNs/dh6Yzbzo9hGFAubuvDlc+ERgLpP7xHgv8LJx+FLjdzCxsn+juFcAaMysP10ca68ycN+cz2FYRI0kMhzXtgj/Kngh+JpM1noc/dz2SieC1yo+CUzgrtkHF9k8eH777SRFIVOyx+c/FD2C9lzA3eRjdRh4HB/SEA3oz5M432ML+6Ib0Ii1Ypz4w+ifw+R/DO0ug/BlYM5sx8VcYb/8K5vnLg3DlK1mPkk5B6AqsT3m+ATi2rnncvdrMtgKdwvaXayzbNZyub52Zc+8XeLz440+e39fI9cWLoXj/8NEW2hwInfoGN8Zoe3D4OCjYJ9i+O8N/8uzuRc858ZNvBlv4oJFBRKTZiMWCsZA6D4LP/YDBE6bRzTYzwNZx14mDmyRCOgWhtq+vNfcz1TVPXe21jbJa674rM7sMuCx8usPMVtaY5UDg3dqWza7NDV3gQOBdu6VhCzV0/kaK6N+yQZQxM/IhI+RHzgZnTPf3eh3wIvDnXzU4U6qe6c6YTkHYAHRPed4N2FjHPBvMrABoD2ypZ9n61gmAu98F3FVXODMrS3f/WJTyIacyZoYyZk4+5MyHjOlK534IrwClZtbbzIqAccCUGvNMAS4Kp88FnvPgaPUUYJyZFZtZb6AUmJvmOkVEpAnV20MIjwlcCcwkOEX0HndfamY3AmXuPgW4G3ggPGi8heAPPOF8kwkOFlcDV7h7AqC2dWb+7YmISLrSug7B3WcAM2q0XZ8yvRM4r45lbwJuSmed+6jO3Uk5Jh9yKmNmKGPm5EPOfMiYlnqvQxARkZZB91QWEREgzwqCmf3czBab2UIze8rMuoTtZma/N7Py8PUhKctcZGarwsdFda89YxlvNbMVYY7HzaxDymvXhBlXmtmpKe1jwrZyM5vQBBnPM7OlZpY0s6E1XsuJjLVkjnT7NbLcY2abzGxJSltHM3s6/Jw9bWYHhO11fjaznLG7mf3LzJaH/9ffy7WcZtbKzOaa2aIw4w1he28zmxNmnBSeeEJ4csqkMOMcM+uV7YwpWeNmtsDMpuVqxoxw97x5AO1Spr8L3BlOfwH4J8F1D8OBOWF7R2B1+POAcPqALGc8BSgIp28BbgmnBwCLgGKgN/A6wQH1eDh9KFAUzjMgyxkPB/oDs4ChKe05k7FG3ki3X0ueE4AhwJKUtl8DE8LpCSn/77V+NpsgY2dgSDi9P8FQMQNyKWe4rbbhdCEwJ9z2ZGBc2H4n8K1w+tspv/PjgElN+H/+A+BvwLTwec5lzMQjr3oI7r4t5el+fHIx21jgfg+8DHQws87AqcDT7r7F3d8HngbGZDnjU+5eHT59meAai10ZJ7p7hbuvAXYN47F7aBB3rwR2DeORzYzL3b3mBX45lbGGqLf/Ke4+m+BsulRj+eQa+PuAs1Laa/tsZjvjW+4+P5zeDiwnGCUgZ3KG29oRPi0MHw6cSDAETm0Zd2V/FBhtZlkf98XMugGnA38Jn1uuZcyUvCoIAGZ2k5mtB74C7DrTqbbhNbrupb2pXELwrYu9ZIk6Y6pczRj19tNxsLu/BcEfY+CgsD3y7OFui8EE38BzKme4K2YhsIngC9vrwAcpX6pSc3xqiBxg1xA52fY74EfArlEpO+VgxozIuYJgZs+Y2ZJaHmMB3P3H7t4deAi4ctditaxqb0NnZDVjOM+PCa69eChXM9a2WFNmbICot98YkWY3s7bAY8BVNXrYe8xaS1vWc7p7wt2PJuhJDyPYnVlXjibPaGZnAJvcfV5q815y5PNnNfdukOPuJ6U569+A6cBPqXuIjA3AqBrts7Kd0YKD12cAoz3cmbiXjOylPWsZ69CkGRsgneFTovaOmXV297fCXS2bwvbIsptZIUExeMjd/56rOQHc/QMzm0VwDKGDmRWE37BTc9Q1RE42HQ+caWZfAFoB7Qh6DLmUMWNyroewN2ZWmvL0TGBFOD0FuDA8U2I4sDXsDs8ETjGzA8KzKU4J27KZcQzw38CZ7v5Rykv5MIxHrmaMevvpSB2+5SLgiZT22j6bWRXut74bWO7uv83FnGZWYuFZeGbWGjiJ4FjHvwiGwKktY21D5GSNu1/j7t3cvRfB5+45d/9KLmXMqKiPajfkQfBtZwmwGJgKdPVPzla4g2D/46t8+syZSwgOjpYDFzdBxnKCfYgLw8edKa/9OMy4Ejgtpf0LBGeBvA78uAkynk3wTaYCeAeYmWsZa8kc6fZrZHkYeAuoCv8dLyXYT/wssCr82bG+z2aWM36WYFfF4pTP4hdyKScwCFgQZlwCXB+2H0rwRaQceAQoDttbhc/Lw9cPbeL/91F8cpZRTmZs7ENXKouICJBnu4xERCR7VBBERARQQRARkZAKgoiIACoIIiISUkEQERFABUFEREIqCCIiAsD/B29q2/hTl9TZAAAAAElFTkSuQmCC\n", | |
| "text/plain": [ | |
| "<Figure size 432x288 with 1 Axes>" | |
| ] | |
| }, | |
| "metadata": { | |
| "needs_background": "light" | |
| }, | |
| "output_type": "display_data" | |
| } | |
| ], | |
| "source": [ | |
| "b = plt.hist(X, bins=100, density=True)\n", | |
| "a = plt.plot(xvals, f_n(xvals))" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 10, | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYQAAAD8CAYAAAB3u9PLAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xd8FNX6x/HPs5tG70gn9I6gMYCAIF1FsXEF9ap0Bf1hFxQbV669oCICIldFBcQWQelNEIHQm0DooQYIgRBItpzfH7tojAlsIJvZ3Tzv12tf2Z2d2fkuCfvsnHPmjBhjUEoppWxWB1BKKRUYtCAopZQCtCAopZTy0oKglFIK0IKglFLKSwuCUkopQAuCUkopLy0ISimlAC0ISimlvMKsDpAbZcuWNdHR0VbHUEqpoLF69epjxphyvqwbVAUhOjqa+Ph4q2MopVTQEJG9vq6rTUZKKaUALQhKKaW8tCAopZQCtCAopZTy0oKglFIK0IKglFLKSwuCUkopQAuCUkopLy0ISimlgCA7U1kpVTBFD5v55/09r91kYZLQpkcISimlAD1CUEoFqMxHBSp/6BGCUkopQAuCUkopLy0ISimlAC0ISimlvLQgKKWUArQgKKWU8tKCoJRSCtCCoJRSyksLglJKKUALglJKKS8tCEoppQAtCEoppby0ICillAK0ICillPLS6a+VUgFDp7y2lh4hKKWUAnwsCCLSTUS2iUiCiAzL5vlIEZnqfX6FiERnem64d/k2EemaafkeEdkoIutEJD4v3oxSSqlLd9EmIxGxA2OAzkAisEpE4owxWzKt1g9INsbUFpFewOvAXSLSEOgFNAIqAfNEpK4xxuXd7npjzLE8fD9KKaUukS9HCLFAgjFmlzEmA5gC9MiyTg/gM+/96UBHERHv8inGmHRjzG4gwft6SimlAowvBaEysD/T40TvsmzXMcY4gRSgzEW2NcAcEVktIgNzH10ppVRe8mWUkWSzzPi4zoW2bW2MOSgi5YG5IvKHMWbJP3buKRYDAapVq+ZDXKWUUpfClyOERKBqpsdVgIM5rSMiYUAJ4MSFtjXGnP95FPieHJqSjDHjjTExxpiYcuXK+RBXKaXUpfClIKwC6ohIDRGJwNNJHJdlnTjgfu/9O4EFxhjjXd7LOwqpBlAHWCkiRUSkGICIFAG6AJsu/+0opZS6VBdtMjLGOEXkYWA2YAc+NcZsFpGRQLwxJg6YCHwhIgl4jgx6ebfdLCLTgC2AExhijHGJyBXA955+Z8KAr4wxs/zw/pRSSvnIpzOVjTE/Az9nWfZCpvvngJ45bDsKGJVl2S7gytyGVUop5T86dYVSyiepGakcOnOIw2cOczTtKKmOVNIcaZxxnCHDnYFd7ITbwgmzhVEkvAilo0r/eatSrAqlokpZ/RbURWhBUEr9jTGGnSd3svHYRrYnb2db8ja2J28nJT0l2/ULhRUi3BaOy7hwuV043U6cxvmP9UpGliS6eDQ1S9akcdnGNC3blNola2O32f39lpSPtCAopdh/ej+/Jv5K/JF4Vh9ZzYlzJwDPh32dknXoVK0T1YtXp2KRilQoUoEKRSpQLKIYUfaobD/Q0xxpJKcnc+LsCY6dPca+0/vYc2oPe1L2sGDfAr7b8R0AhcMK06x8M1pXak2bKm3IebS6yg/iGQwUHGJiYkx8vE57pNTlMsawLXkb8/fNZ8G+BWxP3g5AxSIVuabCNcRcEUOz8s2oVqxann+DN8aw//R+1ietZ33SelYdXsWulF0AuDNK4TzdEMepZrjPVSG74rDntZvyNE+oE5HVxpgYn9bVgqBUwXH87HFm7prJDzt/YEfyDmxio1m5ZnSo1oHrq15PteLWnPx5IPUAyw4s48W507EX2Y7YXLgzyuBIuRJHSgzGUTrb7bQ4XFxuCoI2GSkV4owxrDm6hi+3fsnCfQtxGidNyjbhuRbP0bl6Z8oUKmN1RCoXrcy/6v2LpycVAdtZwoptIrz4eiLKLiSi7EJcZ+qSkdwSV2o9dNZ+/9GCoFSIcrgczNoziy+2fMHWE1spEVmCfzf8Nz1q96BWyVpWx8uZuxDOlGtwplyDhJ0kvORKwkutonDVz3BnlCbjRFscJ2PAhFudNORoQVAqxDhcDn7c+SMTNkzg4JmD1CxRkxdavUD3mt0pFFbIslyZr4aWuannQldJM86SZBzrQsaxjoQV20JE6V+JqvAjEWXn4zjRmtMZ11EsophfcxckWhCUChEOt4O4hDgmbJzAgdQDnmahls/RpnIbbBLszSx2nKeb4DzdGHvhXUSUWURk+dnc8N3vDGgygF71exFpj7Q6ZNDTgqBUkDPGsHD/Qt5Z/Q57T+2lcZnGPNfCUwi808OEEMGVVouzabWwRSVybes1vBX/FpO3TmbwlYO5pdYtel7DZQj2rw1KFWjbTmxjwJwBDF04FJvY+KDDB3x101e0rdI2BIvB37nPVWFc53F80uUTykaV5YXfXqDnjJ7EH9aRiJdKjxCUCkIp6SmMXjOa6dunUzyyOMNjh9OzXk/CbQWvo7VFxRZ8ddNXzN07l7fj36bP7D7cWONGnoh5gvKFy1sdL6hoQVAqiBhjmLVnFq+vfJ3k9GTuaXAPD175ICUiS1gdzVIiQpfoLrSt0paJGycyadMkFu1fxOBmg7m3wb3ajOQjbTJSKkgknk7koXkP8fSSp6lQpAJTbprCM7HPFPhikFmhsEI83PxhfujxA1dfcTVvxb/Fv3/5NwnJCVZHCwpaEJQKcMYYpm2bxu1xt7P26FqGxQ7jyxu/pEGZBlZHC1hVi1dlTMcxvHndmySeTqTnjJ6MWz8Oh9thdbSApgVBqQB2NO0oD81/iP/8/h+uLHclP976I/c0uEebQHwgInSr0Y0fbv2BTtU68eG6D7ln5j3sPLnT6mgBSwuCUgFq1p5Z3B53O6sPr+bZFs8yrvM4KhSpYHWsoFM6qjRvtnuT99q/x+Ezh7lrxl1M2zaNYJrHLb9oQVAqwKQ50hixdARPLX6KasWqMe3mafSu3zsETi6zVsfqHfmux3dcfcXV/Of3/zB04VCSzyVbHSug6F+YUgEkITmBu2feTdzOOAY1HcTnN3xOjRI1rI4VMsoWKsvYTmN5KuYplh5Yyh1xd7Di0AqrYwUMHXaqVID4IeEHRv0+isLhhRnXeRytKrWyOtJlu9A8RVldwQkeCfuejva1FCKdte7afOi8ldWmXp5msomN+xrdR2zFWJ5e8jQD5w7kkeaP0Ldx3wJ/FFaw371SAeCc8xwjlo7g+WXP07RcU6bfPD0kikFuxMpWZkUOo6d9MfHuusx0taShbS/fRIykv933opIb9UvXZ8pNU+hSvQuj14xm6MKhnMo45Zd9BQs9QlDKQofPHGbowqFsOb6FB698kAebPljgRhA1kj18GvEmh01pbne8zG5TEYBXnPfwZvg4RoR/yRmi+NrVMc/3XTi8MG9c9wbNyjfjrVVv0WtGL95p/w71S9fP830FAz1CUMoiq4+s5q4Zd7H31F7ev/59hjQbUuCKAWeTmRDxFicpyj0Zz/5ZDADOEsVQx8MsdjXlpbDPqS/7/BJBRLinwT1M6jaJdFc69/58L7/s/sUv+wp0WhCUymfGGKb+MZX+s/tTPKI4X930FddXu97qWNb4+WnKkcKDGY9yhH9eJtOFnccdD3GSIrwdPhYbbr9FaVa+GdO6T6NRmUY8veRpPlz7IW7jv/0FIm0yUiofOdwO/rviv0zfPp22ldvy2nWvUTyiuNWx8pTPHcm7FsHGaYxx3c4mUzPH1Y5Tgv84/s2HER9wp30x01z+K55lCpVhQpcJ/Of3/zBuwzh2pezildavUDi8sN/2GUj0CEGpfHI64zSD5w1m+vbp9G/Snw86fBByxcBnbhfMHgElqjHWectFV5/hbkm8uy5Phn1DJBl+jRZhj2DktSN5MuZJ5u+bzwOzHuDwmcN+3Weg0IKgVD44mHqQ+365j/jD8Yy8diRDrxpa8PoLMtv0LRzZCJ1eJJ0IHzYQ3nTcRXk5yZ32JX6PJyLc3+h+PujwAftO76P3zN5sOb7F7/u1mhYEpfxs8/HN3PPzPRw+c5ixncdyW53brI5kMQNL34VyDaDR7T5vtcLUZ627NgPtM7Dj8mO+v1xX5Tom3zCZcFs4D8x6gKUHlubLfq3iU0EQkW4isk1EEkRkWDbPR4rIVO/zK0QkOtNzw73Lt4lI1yzb2UVkrYjMuNw3olQgWrR/EX1m9SHcFs4XN3xBy4otrY5kuett6+DoFmjzKNhy851U+Nh5M9VtR+loW+O3fFnVLlWbyTdOpnrx6jw8/2G+3/F9vu07v130tyEidmAMcAPQEOgtIg2zrNYPSDbG1AbeBV73btsQ6AU0AroBH3lf77yhwNbLfRNKBaJvtn/D0IVDqVmiJl/d9BW1S9W2OlJAeCgsDkpUhcZ35Hrbee6rOGxK0cu+0A/Jcla+cHkmdZ1EbIVYXvjtBcauHxuSk+P5Up5jgQRjzC5jTAYwBeiRZZ0ewGfe+9OBjuK5oGsPYIoxJt0YsxtI8L4eIlIFuAn45PLfhlKBwxjD+A3jGbl8JK0rtebTrp9StlBZq2MFhAayl1jbNmjxINhzf7lPF3amutrT3raeShzzQ8KcFY0oyphOY7il1i18tO4jXl7+Mk63M18z+JsvBaEysD/T40TvsmzXMcY4gRSgzEW2fQ94Gvw4sFipfOY2bt5Y9QYfrP2A7jW7M7rD6AIzZNEXve0LSDfh0OzuS36Nac72APnSuZxVuC2cV1q/woAmA/h2x7c8tugx0l3p+Z7DX3wpCJLNsqzHSjmtk+1yEekOHDXGrL7ozkUGiki8iMQnJSVdPK1SFnG4HTy39Dkmb53MvQ3uZVSbUQXyovc5KcQ5brUvZaa7BRT+50lovjpAOVaZetxsX070sBlED5uZq0n0LpeI8H9X/R/PtniWRfsXMWTeEM44zuTb/v3Jl4KQCFTN9LgKcDCndUQkDCgBnLjAtq2BW0RkD54mqA4iMjm7nRtjxhtjYowxMeXKlfMhrlL576zzLI8ufJQZu2bwSPNHePqapwv8zJlZdbf/TnE5y1fODpf9Wj+5WlHHdoB6sv/iK/tJ7/q9+W+b/xJ/JJ6BcwaSkp5iWZa84stf7CqgjojUEJEIPJ3EcVnWiQPu996/E1hgPD0ucUAv7yikGkAdYKUxZrgxpooxJtr7eguMMffmwftRKt+dzjjNg3Mf5NfEX3m+5fMMbDoQTxeayuxO+xJ2uisSnwfTWf/iisVlhJvty/Mg2aW7udbNvN3+bbae2MoDsx4gKS24WzEuWhC8fQIPA7PxjAiaZozZLCIjReT8KYYTgTIikgA8DgzzbrsZmAZsAWYBQ4wx+TOAWKl8kJKewsA5A9lwbANvtnuTf9X7l9WRAlJFjtPC9gffu9qQfUty7hynBMvcjelu+51/tmDnr47VOvJRp484kHqA+2fdz4HUA5bmuRw+HdMaY342xtQ1xtQyxozyLnvBGBPnvX/OGNPTGFPbGBNrjNmVadtR3u3qGWP+MYWgMWaRMaZ7Xr0hpfJL8rlk+s/pz7bkbbzX/j26Rne9+EYFVHfvN/mf3Hl3nYeZ7pZE247QUPbm2WteqpYVW/JJl09ISU/hvl/uY3fKbqsjXRJt5FTqEhw7e4y+s/uyO2U3H3T4gHZV21kdKaD1sP/GOnct9poKefaaC1zNAehgW5tnr3k5mpZryqRuk3C6nfSd3ZddJ3ddfKMAowVBqVw6cuYIfWb14UDqAT7q+BGtK7e2OlJAqyUHaGzbQ5zr2jx93SRKss5di472wCgIAHVL1WVS10kA9J3dl4TkBIsT5Y4WBKVy4VDqIfrM7kPS2SQ+7vQxsRVjrY4U8G6xL8dthJ9cf03bcX6o6OUOF53vas6VspOyBM4In5ola/Jp10+xiY1+c/qxPXm71ZF8pgVBKR/tP72fB2Y9wMlzJxnfeTxXXXGV1ZGCwg22FaxwNyCJUnn+2vPdV2ETw/UBdJQAUKNEDSZ1m0SYLYx+s/ux7cQ2qyP5RAuCUj44kHqAfrP7ccZ5hk+6fkLTck2tjhQUasgh6toOMNsd45fX32Kqc9CUpmOA9CNkVr14dSZ1nUSkPZJ+c/qx9XjgT9umBUGpiziUeshTDBxnmNB5Ag3LZJ3bUeWksy0egDku/xQEEBa7ruRa22ZwBd68QtWKV2NSt0kUDitMvzn92Hxss9WRLkgLglIXcPjMYfrO7sup9FOM7zKeBmUaWB0pqHSxr2ajO5qD+G9yv2XuxhSXNDi0zm/7uBxVi1VlUrdJFI8ozoC5AwL6SEELglI5OJp2lP5z+nMy/STjOo+jUZlGVkcKKuU4yVWyw49HBx6/ub2/l135OyV2blQuWpmJXSdSNLwoA+cOZEfyDqsjZUsLglLZOHb2GP1m9yMpLYmxncbSpFwTqyMFnY72NdjEMMdP/QfnnaA4m9zRsHORX/dzuSoXrczELhOJsEXQf05/dqUE3nkKWhCUyuL42eP0n92fI2lHGNtpLM3KN7M6UlDqYotnr7s820zVi698mZa6G8P+FZAR2LOOVi1elU+6foIg9J/dn32n9lkd6W+0ICiVSfK5ZAbMHcCB1AOM6ThGh5ZeoiKcpbVtk/fowP8T/S1zNwa3A/ZaO9mdL2qUqMGELhNwuB30m9MvoOY+0oKglNepjFMMnDuQfaf28UHHD7imwjVWRwpabWybiBQn8935U1BXuuuDPTKg+xEyq1OqDhO6TCDNkUa/2f04fOaw1ZEALQhKAZDmSGPIvCEknExg9PWjaVmx5cU3UjlqZ1vHKVOIeHfdfNlfOhFQrQXsWgz8/Uzo/Lx4Tm7UL12fcZ3HkZKeQv85/QNi6mwtCKrAc7gcPL74cTYc28DrbV/XuYkum6G9fT3L3I1xEpZ/u63eBo5sgrMn82+fl6lx2caM7TT2zxFtJ86dsDSPFgRVoLncLp5d+izLDizjhZYv0CW6i9WRgl5dSaSSnGCRO58746tfCxhP53IQaVa+GWM6juFg6kEenPsgpzNOW5ZFC4IqsIwxjFoxill7ZvH41Y9zR907rI4UEtrbPCeILXbl8/QeVWLAFg57l+XvfvPANRWu4d3r32XHyR0MmT+ENEeaJTm0IKgC64O1H/DN9m/o17gffRr3sTpOyGhvW88f7qocpkz+7ji8EFS+KihGGmWnTeU2vN72ddYnreexRY+R4crI9wxaEFSB9L9N/2PCxgncWfdOhl411Oo4oSP9NDG2bSxyX5nvu44eNpOPdl+BY/9qokjP9/3nhS7RXXip1Uv8dvA3nlnyDE53/s7PpAVBFTjf7fiOt1e/TdforoxoMQIR/4+TLzB2LSZCXCy2oCAArHTXI1xcNLcF14VpMrutzm08c80zzNs3jxd/exG3cefbvvNxCIBS1pu7dy4vL3+Z1pVa82qbV7Hb7FZHCi0J80g1UcS761my+9XueriNECt/sJzgnXvq3ob3kupIZcy6MRQNL8qw2GH58sVFC4IqMM4fhjcp24R32r9DuD3c6khBLfP4/j2v3QTGQMI8lrkb47Doo+U0hdliqhNr+wNclkTIM4OaDiI1I5XPtnxG0YiiPNL8Eb/vUwuCKhA2JG3g0YWPEl0imjEdx1A4vLDVkUJP0jZI2c8it7VDd1e569HLvpBwnJYVprwgIjwR8wSpjlS2n9iO0+0kzObf9xO8/1pK+WhH8g4emvcQZQuVZVyncZSILGF1pNCUMBeAxS5r+g/OW+FuQJ+w2TSW3aw1dSzNcrlEhOdbPo8bt9+LAWinsgpxiacTGTR3EJH2SMZ3Hk+5wuWsjhS6dsyFcvX9ejEcX6z2TpdxlS14Lm5/IXabnXBb/jRvakFQISspLYkBcwaQ4c5gfOfxVClWxepIoSs9FfYth9qdrE5CEiXZ7y7HVbbAvAhNINMmIxWSUtJTGDRvEMfPHeeTLp9Qu1RtqyOFtt1LwJUBdTrDQuuvSbDG1KGFbStgAPlnB7jKlh4hqJCT5khjyPwh7EnZw+jrR9O0XD5PoVAQJcyF8CJQrZXVSQBY7a5DBUmmEsetjhJUtCCokOJwOXh80eNsPLaRN657g1aVAuMDKrQZ2DEParaDsEirwwCwxu3pTNZmo9zxqSCISDcR2SYiCSIyLJvnI0Vkqvf5FSISnem54d7l20Skq3dZlIisFJH1IrJZRF7OqzekCi6X28XwpcNZdnAZL7V6iU7VrW/PLghqyUFI2RcQ/Qfn/WGqcdZEaEHIpYsWBBGxA2OAG4CGQG8RaZhltX5AsjGmNvAu8Lp324ZAL6AR0A34yPt66UAHY8yVQDOgm4joFUnUJTPG8MqKV5i9ZzZPxjzJbXVuszpSgdHett5zJ4AKgpMwNpiaITPSKL/4coQQCyQYY3YZYzKAKUCPLOv0AD7z3p8OdBTPedY9gCnGmHRjzG4gAYg1Hqne9cO9N3OZ70UVYKPXjGb69un0b9Kf+xvdb3WcAqWdbT2UrQulqlsd5W9Wu+vSSPYSSf7PGhqsfCkIlYH9mR4nepdlu44xxgmkAGUutK2I2EVkHXAUmGuMyfaqFiIyUETiRSQ+Kcn6S8ypwPPppk+ZuGkiPev25P+a/5/VcQqUKNJpYfsDane2Oso/rHHXIVxcNJFdVkcJGr4UhOxmVMr6bT6ndXLc1hjjMsY0A6oAsSLSOLudG2PGG2NijDEx5crpSUXq777d/i3vrn6XbtHdeK7FczpzaT5rZdtCpDigTuA0F5231u0Zany19iP4zJeCkAhUzfS4CnAwp3VEJAwoAZzwZVtjzElgEZ4+BqV8NmfPHEb+PpLWlVvz3zb/1ZlLLdDOtp40EwnVrrU6yj8cpwR73Fdox3Iu+FIQVgF1RKSGiETg6SSOy7JOHHC+4fZOYIExxniX9/KOQqoB1AFWikg5ESkJICKFgE7AH5f/dlRB8dvB33jm12e4styVvNv+XZ251CLtbetZ7m4I4VFWR8nWalPHWxC0i9IXFz1T2RjjFJGHgdmAHfjUGLNZREYC8caYOGAi8IWIJOA5Mujl3XaziEwDtgBOYIgxxiUiFYHPvCOObMA0Y8wMf7xBFXrWHV3HowsfpWaJmnzY8UMKhRWyOlKBVF0OE207wkTHDfTLdCZwIFnrrsMd9qVUkSQSTXmr4wQ8n6auMMb8DPycZdkLme6fA3rmsO0oYFSWZRuA5rkNq9T25O0MmT/EM3Np53EUjyhudaQC6/xwU6uujuaL8yeoXS3btSD4QM9UVkFj/+n9DJo7iCh7FOM7j6dsIWtn1Szo2tnWs8tdgX3mCquj5GibqUqqidJ+BB9pQVBBISktiYFzBuJwOxjXeZzOXGo1xzla2bYE9NEBgAs7G9w1g/oay/lJC4IKeCnpKQycO5Dj544ztuNYnbk0EOxdSiHJCPiCAJ6ZTxvIPqJItzpKwNOCoAJamiONwfMHs/fUXt7v8D5NyjWxOpICSJjPORPuGWEU4Na6a3tPUNttdZSApwVBBawMVwaPLnyUTcc28cZ1b9Cyok53FTB2zGGFuwHpRFid5KLW6synPtOCoAKSy+1i+K/DWX5ouc5cGmiO7YDjCcxzX2V1Ep+coDi79QQ1n2hBUAHHGMPI30cyZ+8cnbk0EG37BYD5ruAoCABrTR1vx7KeoHYhWhBUQDHG8Hb823y34zsGNh2oM5cGom2/wBVNOEjwDPtd465DeTlJFTlmdZSApgVBBZQJGyfw2ZbP6F2/Nw83e9jqOCqrtBOw/3eod4PVSXLl/ER3zUWbjS5EC4IKGF9t/YoP1n7AzTVvZljsMJ25NBDtmAPGHXQF4Q9TjTQTqf0IF6EFQQWEn3b+xKsrX+X6qtczsvVIbKJ/mgFp289QtAJUbGZ1klxxYWeDqUlzLQgXpP/rlOUW7FvA88uep0WFFrzZ7k3CbD5NsaXymzMdEuZDvW5gC76PjrXu2jSSveA4a3WUgBV8v1UVUlYcWsFTi5+iYZmGjO4wmkh7pNWRVE52/woZqVA3uJqLzjt/BTUOrbc6SsDSgqAsszFpI48seIRqxasxttNYioQXsTqSupAtP0BEMajZ3uokl+T8CWokrrI2SADTgqAssSN5Bw/Nf4gyUWUY13kcJSJLWB1JXYjLAX/M8HQmB+jFcC7mGCXY5y4H+1daHSVgaUFQ+e78NNYRtggmdJlA+cI6T33A270EziZDo1utTnJZ1pg6kBhvdYyApQVB5asjZ44wYM4AMtwZjO88XqexDhbnm4tqdbQ6yWVZ664Dpw9CSqLVUQKSFgSVb46dPUb/Of1JPpes01gHE5cDts7wjC4K0uai885fQU2bjbKnBUHli5PnTjJw7kCOpB3ho04f6TTWwWT3Ejh7AhoGd3MReE5QIyxKm41yoAO+ld+dyjjFwLkD2ZuylzGdxnD1FVdbHUnlxoapEFUCagf/jLMOwqBSc0jUI4Ts6BGC8qszjjM8NO8hdpzcwbvXv6vXNAg26adh60/Q6Pagby76U5VrPOciOPUKallpQVB+c9Z5liHzh7D52Gbeuu4trqtyndWRVG5t+REcadDsbquT5J0q14ArAw5tsDpJwNEmI+UX6a50hi4Yypoja3it7Wt0rB7co1MCTfSwmX/e3/PaTf7b0bqvoXQtz4doqDj/XhJXQtUQel95QI8QVJ5zuBw8segJlh9azsjWI7mx5o1WR1KXInkP7F0KV/aGUJp5tnhFKFFVz1jOhhYElaecbifP/PoMixMXM6LFCG6tHfwjUwqstV8CAlfeZXWSvFflGtivBSErLQgqzzjdTp799Vnm7p3LUzFPcVf9EPwgKSic6bB6EtTtCiWrWZ0m71WNhVOJcOqg1UkCihYElSfOF4Nf9vzCY1c/xn2N7rM6krocW+LgTBLEDrA6iX/82Y+gRwmZaaeyumxOt5Phvw5n1p5ZPH714/Rp3MfqSPkut528+dYpfKlWjvd0JtfsYHUS/6jQFOyRnjOWG/awOk3A8OkIQUS6icg2EUkQkWHZPB8pIlO9z68QkehMzw33Lt8mIl29y6qKyEIR2Soim0VkaF69IZW/tBiEoH0rPCNwYgcE5YU9zmVHAAAeuklEQVRwfBIWARWv1DOWs7job1tE7MAY4AagIdBbRBpmWa0fkGyMqQ28C7zu3bYh0AtoBHQDPvK+nhN4whjTAGgJDMnmNVWAc7qdDPt1GLP2zOKJq5/QYhAqfn0LCpeBq0K82a9qLBxcC84Mq5MEDF/KfyyQYIzZZYzJAKYAWY+xegCfee9PBzqK5wrpPYApxph0Y8xuIAGINcYcMsasATDGnAa2ApUv/+2o/HK+GMzeM5snY57kgcYPWB3JZ9HDZv55U1kcXAc75kDLwRAR4hcsqnINuNLhyEarkwQMXwpCZWB/pseJ/PPD+891jDFOIAUo48u23ual5sCK7HYuIgNFJF5E4pOSknyIq/zN6XbyzJJn/iwG9ze63+pIKq8sfgMiS4RuZ3Jm5zuWdfjpn3zpVM7ujBTj4zoX3FZEigLfAo8aY05lt3NjzHhgPEBMTEzW/ap85nA5GPbrMObsnRNyxcCXjt7cHlXk9jXzar+XZPcS2DYTOjzvmcwu1JWoDMWrwP7foeWDVqcJCL4cISQCVTM9rgJkHbz75zoiEgaUAE5caFsRCcdTDL40xnx3KeFV/kp3pfPYoseYs3cOT8U8FVLFoMBzu2D2s1CiGrQaYnWa/FO9FexdDka/a4JvBWEVUEdEaohIBJ5O4rgs68QB5z8d7gQWGGOMd3kv7yikGkAdYKW3f2EisNUY805evBHlX2mONIbMH8KSxCU83/J5Pc8g1MR/Coc3QueXILyQ1WnyT/VrIfUwnNhldZKAcNEmI2OMU0QeBmYDduBTY8xmERkJxBtj4vB8uH8hIgl4jgx6ebfdLCLTgC14RhYNMca4RKQN8G9go4is8+7qWWPMz3n9BtXlO51xmiHzh7A+aT2j2ozi5lo3Wx0pX4V85/OxBJjzvOfymI1utzpN/qp2refn3t+gTC1rswQAn05M835Q/5xl2QuZ7p8Deuaw7ShgVJZlS8m+f0EFmORzyQyaO4gdJ3fwVru36Fy9s9WRVF5yZsD3gyAsEnp8GFqT2PmiXD3PENt9y+Gqf1udxnJ6prLK0bGzxxgwZwD7Tu1j9PWj9XoGXnZcVJTjlOckV0gyJeQMkTiIJINwXLB0O9gjwBbO3fY/yCCMMyYKEiIhoihEFKWqHCHNRJFKIdIJ5/z3o3w9g9kYmPEYHIiHnv8j+r9rgbX/2HfAn1V9OUSgWivYu8zqJAFBC4LK1qHUQ/Sf05+ks0l81OkjWlRsYXUky5QnmWttm2ll20Ij2x7qyAEixZHzBvP+uvvf8EzLJ4/+8+6vkX8tdhobaURxwhTjKCU5akqSZErC8r2eb7Dl6pPzoL1LZAzMexHWTYZ2z0Cj24AQbxrLSfVr4Y8ZkHLAM/KoANOCoP5hV8ouHpz7IKkZqYzvPJ5m5ZtZHSn/nT4Mm76Djd+wMmoNACdNETa4a/KZuws7TSUOm9IcMaU4aYpwjgjSCcdJGDtGdga3A1wOYv8zi0jJoDDpzH7oKshIhYwzPDF5GYXlHEU55/15ltJymvKcpIHso71tPcye/Wec9ZGFWeOuQ7y7HmtMHda6a3OOyJzSX5jjHPz8JKz9AmL6Qbt/zEZTsFT39iPsWw5N7rQ2i8W0IKi/2ZC0gSHzh2ATGxO7TqRBmQZWR8pXTWUn/cJ+gXdXgtsJFZryuqMXS9xN2GKqY3wZmBdR+M+7Ryn115k31f46yvrWbb/Iixj2jGgBSX/AsW3M/HEGMbbtXB++HgCHsbPe1GKFu76nKapqS4gsevHmnT1L4een4ehmaPskdBhR8PoNMokeNhM7LtZFFqLY3t+0IFgdQAWOZQeW8diixygTVYZxncdRrXgIzoOfk0Pr+Sz8NdrZN3DaFILYQXD1A1CuLmMtGWUkULSc51ajLc9+WxGA4qRylW0HLWx/0MK2lYH2mTA5DsQOlZoxLKwim93RHDBl4VhdMG5ISYTDGzxTWh9c4zkZq/dUqNfNgvcVeFzYWeOuQ7u9v1kdxXJaEBQAM3bN4Pmlz1O7VG3GdhpL2UJlrY6Up3IcOnr6MMwZARu/oamtKK86evOlqyOpiwrDoh3AjnzNeTGnKMoid3MWuZsDUJhzbBlQCvYsg73L6GOfRWSY07Pyhy/9feMrmsANb0Lze0J/nqJcWuGuT7ukaXDmOBQpY3Ucy2hBUHyx5QveWPUG11S4htHXj6ZYRDGrI/md4OZu+wL48CFwnoO2T3Ld3PqcpvDFNw4gaURBrQ6eG9B02PdUlyNUlmNM6lXP0xxUrCKUrQNFy1ucNnCtdNf33Nm3HBp0tzaMhbQgFGDGGEavGc3ETRPpXL0zr7Z9lUj7JXZUBpHKJPFuxEfE2rZBpeug+3tQphan5wb/KJt0IthuqrLdVIWmITZE1I82mFoQFgV7ftWCoAoeh8vBS8tfIm5nHP+q+y+ebfEsdtvFOjqDX3fbcv4bPhHB8KRjEG/d93rAdqqG/BnSASSDcKjWEnYttjqKpUL0ckjqQlLSUxg0bxBxO+MY0mwII1qOCPliEIGDV8Mm8GHEBySYStyQ8SrTXe0CthgoC9RsD0lb4fQRq5NYRo8QCpjE04kMnj+YxNOJvNr2VbrXDK3D42yHXZ46xNSI/9DclsCHzh6867wTF/4rgPn9zd6X/eU0HFWPQjKp0c7zc/cSaJrtTDwhTwtCAbIhaQOPLHgEp9vJ+M7jiakQY3Uk/9u3Aqb9m7pykkEZjzLbHWt1IhWoKl4JUSVh16ICWxC0yaiAmLt3Ln1n96VwWGEm3zi5YBSDrT/BZ90hvDC3ZYzUYqAuzGaHGm09BaGAXh9BjxBCnDGG/23+H++ufpem5Zryfof3KR1V2upYftfTvgimfQKVr4a7p7F95HKrI4WEkG9iqtne80XixK4COR22FoQQds55jpeXv8yMXTPoGt2VV1q/QlRYlNWx/G6AfQbPhX8FNTvAXZP1JCzluxrtPT93LdKCoELHkTNHGLpwKJuPb+aR5o8woMkAJNRH1BjDM2Ff81DYT8xwteSxzQ/geGGR1alUMClTyzO1x65FcE0/q9PkOy0IIWh90noeXfgoaY40Rl8/mg7VOlgdyf/cLpjxKA+F/cRkZ0decPbBrV1kKrdEPM1Gf8zw/E2F+HDsrPR/TIj5IeEH+szqQ5Q9ii9v/LJgFANnOnxzP6z5nPedtzLC2VeLgbp0tTvCuZOQGG91knynRwghwuF28E78O0zeOpkWFVvw1nVvUTKqpNWx/C/9NEy5B3Yvhq6v8s6P1a1OpIJdrQ6e2WN3zP7blOUFgX6NCgFHzhyh76y+TN46mXsb3MvHnT4uGMXgzHH47BbPHP+3fgytBludSIWCQiU901hsn2N1knynRwhB7vdDv/PMkmc46zzLG9e9wQ01brA6ks8u51q9LYd9zhcRr1FNjhLZezLUv/GyM4SaUH5v/vC3v8fuXTyXGD11EIpXsjBV/tIjhCDlNm7GbxjPoLmDKBlZkik3TQmqYnBZju1geuTLXCEnuC9j2CUXA6VyVLer5+eOgnWUoAUhCKWkp/DIgkf4YO0HdI3uytc3fU3NkjWtjpU/Dq6DT7sRRQa9M55nhSlYl/hU+aRcfShRrcA1G2mTUZCJPxzPsF+HcfzccZ5t8Sy96vUK6fML/nYYP6g4fN0bCpWkZ8ZwdpuK/1hHXZz+e/lABOp0hvVTPKPYwkL/OiGgRwhBw+l28uHaD+k3px+R9ki+uOELetfvHdLFILMutlUw+Q4oURn6zv6zGCjlN3W7guMM7F1mdZJ8owUhCBxIPUCfWX0Yt2Ec3Wt2Z9rN02hctrHVsfLNnfbFjA1/Dyo0hj6/eIqCUv5W4zoILwx/FJwjKi0IAW7W7ln0jOtJwskEXm/7OqPajKJIeAGZm8cYBtl/4q3wcfzmbgT3xUHh0J+YTwWI8EKeZqOtP4HbbXWafOFTQRCRbiKyTUQSRGRYNs9HishU7/MrRCQ603PDvcu3iUjXTMs/FZGjIrIpL95IqElJT2H4r8N5aslT1ChZg29u/oYbaxag0TRuN8wZwfDwr/nJ1ZJ+jqcgsqjVqVRB0+AWSD0C+1dYnSRfXLRTWUTswBigM5AIrBKROGPMlkyr9QOSjTG1RaQX8Dpwl4g0BHoBjYBKwDwRqWuMcQH/Az4EPs/LNxQKfk38lZd+e4kT507w0JUPMaDpAMJt4VbHyj8uB/w4BDZMZZKzKyOd/8aEwMFsMHXmBlNWv6rbFeyRsOVHqN7K6jR+58v/slggwRizyxiTAUwBemRZpwfwmff+dKCjeHo7ewBTjDHpxpjdQIL39TDGLAFO5MF7CBmpGam8+NuLDJ4/mOKRxfnypi8Z3GxwwSoGGWc8I4k2TIUOI3jZeV9IFAMVpCKLeeY22vpTgbhoji/DTisD+zM9TgSyTvDx5zrGGKeIpABlvMt/z7Kt9ghmY8WhFTy/7HmOpB2hX+N+DG42mAh7hNWx8k30sJmU5DSTIt6kuX0X3Dwarn4AftZvqspiDW6BbT/DgTVQ5Wqr0/iVLwUhu3GNWUtlTuv4su2Fdy4yEBgIUK1atdxsGhRS0lN4O/5tvk/4nuji0Xx+w+dcWe5Kq2PluypylP+Fv0FVSYJ/fQENulsdSSmPejeALRw2fRvyBcGXY/FEoGqmx1WAgzmtIyJhQAk8zUG+bHtBxpjxxpgYY0xMuXLlcrNpQDPGMHPXTG754RbidsbRp3Efpt08rUAWg2aSwA8RL1BWUvh3xjAtBiqwFCrp6UvYNB1cTqvT+JUvRwirgDoiUgM4gKeT+O4s68QB9wPLgTuBBcYYIyJxwFci8g6eTuU6wMq8Ch+s9p/ezyu/v8JvB3+jSdkmjO88nnql61kdyxLdbCt5L3wMR0wp+jieZpcpOBOJqcD1j4kXm97luWjOrkVQp5N1wfzsogXB2yfwMDAbsAOfGmM2i8hIIN4YEwdMBL4QkQQ8Rwa9vNtuFpFpwBbACQzxjjBCRL4G2gNlRSQReNEYMzHP32EAcbgcfL7lcz5e/zF2m53hscO5q95d2AvYVZkAMIaB9p94NvxrVrvrMCDjCU5Q3OpUSmWvbleIKgkbpoR0QRATRD3nMTExJj4+OK9itCRxCW+seoO9p/bSoWoHhrcYToUiFayOZQ2XA35+ClZPYoarJU84HiSdgtOBroLLn1Ozz3gM1n0NT+3wjD4KEiKy2hgT48u6Ormdn+1O2c0bq95g6YGlRBePZmynsbSp3MbqWNY5c9xzucs9vzLGeQtvOf+lw0pVcGjaC+I/hS1x0Pweq9P4hRYEP0nNSGXchnFM3jqZKHsUT8Y8yd317ybcXoDOKcjq0AbP5S5Tj8CtH/PmFG0iUkGkaiyUqQ1rPgvZgqBfzfKYw+Xgy61fctP3N/G/zf/j5po389NtP3F/o/sLdjHYOB0mdgHjgr6zoFlvqxMplTsi/OdwS9i/ghuGf2R1Gr/QI4Q84jZuZu2exftr3+dA6gFiK8Ty+NWP06hsI6ujWcvtgvkvw7LRULUl3PUFFC1vdSqlLsl013U8FTaVe+3zgNC7hrcWhMtkjGH5weW8t+Y9tp7YSr1S9fi408dcW+naAnOtgqzOD9krzSnWNJoGOxdATF/o9jqEaeexCl4pFOUnVytutS+Fc6cgKrSaPbUgXCJjDCsPr2Ts+rGsPrKaykUr82rbV7mxxo3YRFviYmUr70d8CHvS/pqGQqkQ8IWrMz3DlniuptZioNVx8pQWhFwyxrD80HLGrR/HmqNrKF+oPMNih9Gzbs8CNfdQjtxuBtt/5Imwaew1V8CAn6BCE6tTKZVnNpharHPXotnvYzxHvvbQ+RgNnXfiZ8YYfjv4G2PXj2V90nquKHwFz7Z4ltvr3E6kvWBcbxVynhZ5z2s3eYaUfj+Qp8Pn8ZOrJcMd/dmUqRjolMoqVIx13sK45Hdhyw/Q5E6r4+QZLQgX4XA7mLV7Fp9t/oxtyduoUKQCI1qM4LY6t+kRQWY75sGPg+FsMiMcfZjs6kT2cxsqFfzmuK+GsnVh6XvQ+A4Ikf5CLQg5OJ1xmm+3f8vkrZM5knaEmiVq8vK1L9O9ZnctBJlEkc7wsK/gy7lQrgHcM53Jo/dffEOlgpjBBq0f9XwJSpjnudRmCNCCkMW+U/uYum0q3+74ljOOM8RWiOWFVi/QpnIb7SzOorHsYnT4GGrZDkHLwdDxRQiP4u+Xz1AqRDX9Fyx61XOr3SkkjhK0IABOt5MliUuYtm0ayw4uI0zC6Fy9M/c3vp9GZQr4eQTZCMPJIPsMHg37lmOU4O6MZ/mq2zNWx1Iqf9nDof0wz+Vet/wIjW61OtFlK9AF4djZY3y34zu+2f4Nh88cpnzh8gxuNpg76txB+cJ68lR2Gslu3ggfTyPbXuJcrRjh6MMpimqHsQpZF/zbvrI3LB/jOfmy/k2eIhHEClxByHBlsGj/In7c+SPLDizDZVy0qtiKYdcMo13VdoTZCtw/iW8cZ2HRa/wY8T4nKM6gjEeZ7Y61OpVS1rLZodNL8NW/YPX/IHaAxYEuT4H49DPGsOnYJn7c+SO/7P6FUxmnKF+4PH0a96FHrR5El4i2OmJg27MU4v4PTuzkW1d7Rjnv5hRFrU6lVGCo0wWi28KCV6Bhj6CemiXkC0KaI427Z97NzpSdRNoj6VitIz1q9aBFxRYF88I0uXH6MMx9ATZMhZLV4d8/8MyENKtTKRVYROCmd+Dj1jBrGNz5qdWJLlnIF4TC4YWJqRDDvQ3vpWt0V4pFBM+FLSzjzIAVH8Pi18GVAW2fhLaPQ0QRQPsKlPqHcnXhuqdg4SjPdRPqdrE60SUJ+YIAMKLlCKsj5Jt/XAs2V9vOoIstnqfDplLbdhDqdoOu/4UytfIkj1IhrfWjsOlb+On/4MGlUKSs1YlyTQfWK489y/g24iXGR7yLYKD3VLh76mUVA6UKlLAIuH0CpJ2A7waC2211olwrEEcIKgfGcNezbzMk7Aeus2+kkpTmaccAvnVdh2uSi/PNQ7k90lCqwKrYFLq9CjMfh6Vve5qRgogWhILI7Ybts2DpO0yNXEWSKc5/Hb35zNVVL3av1OWK6Qt7f/OMOipd0zPXUZDQglCQpCbBusmw+jNI3g0lqzPC0YdvXO20ECiVV0Sgx4dw6gB8/yAUKQc1rrM6lU+0IIQ6ZzokzIeN02DrDHA7oHpr6DACGt7K5OdmX/QltGNYqVwKLwS9v4ZPb4Cve0OvL6Fme6tTXZQWhBAUSQatbZvguzjY9jOkn4JCpeCa/hDTB8rVszqiUqGvUCm47wf44jb4sifcMREa3mJ1qgvSghACoofNoLYcoJ1tA5+HryfW9gdR4oDtJTx/gI1ugxrt/pxnRb/xK+U//xj6/cBMz9QW0+6D656E9sM9U14EIC0IwcjlgMMbYN8K2P87KyIXc4WcBCDBXYkvXZ1Y7G7K5yOe1IvaK2W1wqXhvjj4+SlY8ibsXQ63vB+QQ7q1IAQ6twuO7/QUgEPr4cAaOLAanGc9z5eoxu/uhvzubsBi15UcJNPJMFoMlAoMEYXh1jFQ/VrP9BYftYI2j0KrhyGquNXp/qQFwUJZm25KcpracoDpd5SGwxs9ReDIZnB45g9KN2H8Yaqx2t2Ovr3ugqotoERlhmoTkFL57pJmBWh+D9TuCLOf9UwNs2IctHgQrr4filfyU1Lf+VQQRKQbMBqwA58YY17L8nwk8DlwNXAcuMsYs8f73HCgH+AC/s8YM9uX1wxZLgek7IcTu+hr/4XacoBatoPUkoOUlVOedWYAkcWhQhO46j6o0BQqNqXxe7tweH9lfRvryWJKBaViFTwT4LV62NOEtPg1WPIG1OrgOWehVkcodoUl0cQYc+EVROzAdqAzkAisAnobY7ZkWmcw0NQY86CI9AJuM8bcJSINga+BWKASMA+o693sgq+ZnZiYGBMfH5/7d8nlzfHjM2PgXAqkHoXUI55xyMl7+WbeUqrakqgiSVSxnQDz1yntyaYoCaYSO92VSDCVPfdNJRJNOQy2v2XVzmClQsPfPoOO74R1X8H6rz2fGQDlG/H1wfJsNDXZ7K7Ojy/2ueSmJRFZbYyJ8WVdX44QYoEEY8wu74tPAXoAmT+8ewAvee9PBz4UEfEun2KMSQd2i0iC9/Xw4TXzzoE1NJcd2HBjw8Du4p4PZePy/HS7szz2/jx/c7s8z2WkeYZwpp+C9NN/3c4c+6sIuNL/sfu29lLsN+VY6a5PlXbXQqnqUKoGV328jxMUA4L/WqxKqUtUphZ0fB6ufw6ObIKEebB7Cd3sq+gtCz3rfDIZHl7l9yi+FITK/P2q6YlAi5zWMcY4RSQFKONd/nuWbSt771/sNfPOpBv5PvLsX48/u8zXs0dCZDHvrSgULgtlansujFH0Cu+tvKdNsERVWj4//89Nb+/w1zeDE5y8zCBKqZBhs3nmQqrYFNo+TvNhM6giSTSUvYzv0DxfIvhSELL7+pq1nSmndXJant0sq9m2XYnIQGCg92GqiGzLskpZ4Fh22/pXUm43KAsck9dzt1Fu179MFv1b5opmzBvBkBGCI2euM/r6/3ovsAyY8GquM2VW3dcVfSkIiUDVTI+rAAdzWCdRRMKAEsCJi2x7sdcEwBgzHhifUzgRife1fcxKwZBTM+YNzZh3giFnMGT0lS/XQ1gF1BGRGiISAfQC4rKsEwfc771/J7DAeHqr44BeIhIpIjWAOsBKH19TKaVUPrroEYK3T+BhYDaeIaKfGmM2i8hIIN4YEwdMBL7wdhqfwPMBj3e9aXg6i53AEGOMCyC718z7t6eUUspXPp2HYIz5Gfg5y7IXMt0/B/TMYdtRwChfXvMS5dicFGCCIadmzBuaMe8EQ85gyOiTi56HoJRSqmDQayorpZQCgqwgiMh/RGSDiKwTkTkiUsm7XETkfRFJ8D5/VaZt7heRHd7b/Tm/ep5lfFNE/vDm+F5ESmZ6brg34zYR6ZppeTfvsgQRGZYPGXuKyGYRcYtITJbnAiJjNpkt3X+WLJ+KyFER2ZRpWWkRmev9O5srIqW8y3P82/RzxqoislBEtnp/10MDLaeIRInIShFZ7834snd5DRFZ4c041TvwBO/glKnejCtEJNrfGTNltYvIWhGZEagZ84QxJmhuQPFM9/8P+Nh7/0bgFzznPbQEVniXlwZ2eX+W8t4v5eeMXYAw7/3Xgde99xsC64FIoAawE0+Hut17vyYQ4V2noZ8zNgDqAYuAmEzLAyZjlryW7j+bPNcBVwGbMi17AxjmvT8s0+8927/NfMhYEbjKe78YnqliGgZSTu++inrvhwMrvPueBvTyLv8YeMh7f3Cm//O9gKn5+Dt/HPgKmOF9HHAZ8+IWVEcIxphTmR4W4a+T2XoAnxuP34GSIlIR6ArMNcacMMYkA3OBbn7OOMcY4/Q+/B3PORbnM04xxqQbY3YD56fx+HNqEGNMBnB+Gg9/ZtxqjMl6gl9AZczC6v3/jTFmCZ7RdJn14K9z4D8Dbs20PLu/TX9nPGSMWeO9fxrYimeWgIDJ6d1XqvdhuPdmgA54psDJLuP57NOBjiLi93lfRKQKcBPwifexBFrGvBJUBQFAREaJyH7gHuD8SKfspteofIHl+aUvnm9dXCCL1RkzC9SMVu/fF1cYYw6B58MYKO9dbnl2b7NFczzfwAMqp7cpZh1wFM8Xtp3AyUxfqjLn+NsUOcD5KXL87T3gaeD8rJRlAjBjngi4giAi80RkUza3HgDGmOeMMVWBL4GHz2+WzUtdaOoMv2b0rvMcnnMvvgzUjNltlp8Zc8Hq/V8OS7OLSFHgW+DRLEfY/1g1m2V+z2mMcRljmuE5ko7F05yZU458zygi3YGjxpjVmRdfIEcw/60G3gVyjDGdfFz1K2Am8CI5T5GRCLTPsnyRvzOKp/O6O9DReBsTL5CRCyz3W8Yc5GvGXPBl+hSrHRGRisaYQ96mlqPe5ZZlF5FwPMXgS2PMd4GaE8AYc1JEFuHpQygpImHeb9iZc+Q0RY4/tQZuEZEbgSigOJ4jhkDKmGcC7gjhQkSkTqaHtwB/eO/HAfd5R0q0BFK8h8OzgS4iUso7mqKLd5k/M3YDngFuMcakZXoqGKbxCNSMVu/fF5mnb7kf+DHT8uz+Nv3K2249EdhqjHknEHOKSDnxjsITkUJAJzx9HQvxTIGTXcbspsjxG2PMcGNMFWNMNJ6/uwXGmHsCKWOesrpXOzc3PN92NgEbgJ+Ayuav0Qpj8LQ/buTvI2f64ukcTQD65EPGBDxtiOu8t48zPfecN+M24IZMy2/EMwpkJ/BcPmS8Dc83mXTgCDA70DJmk9nS/WfJ8jVwCHB4/x374Wknng/s8P4sfbG/TT9nbIOnqWJDpr/FGwMpJ9AUWOvNuAl4wbu8Jp4vIgnAN0Ckd3mU93GC9/ma+fx7b89fo4wCMuPl3vRMZaWUUkCQNRkppZTyHy0ISimlAC0ISimlvLQgKKWUArQgKKWU8tKCoJRSCtCCoJRSyksLglJKKQD+H0CrJsdQfovoAAAAAElFTkSuQmCC\n", | |
| "text/plain": [ | |
| "<Figure size 432x288 with 1 Axes>" | |
| ] | |
| }, | |
| "metadata": { | |
| "needs_background": "light" | |
| }, | |
| "output_type": "display_data" | |
| } | |
| ], | |
| "source": [ | |
| "# acceptance-rejection method\n", | |
| "b = plt.hist(X, bins=100, density=True)\n", | |
| "a = plt.plot(xvals, f_n(xvals))\n", | |
| "m = 45\n", | |
| "g_mu = 250\n", | |
| "g_std = 240\n", | |
| "c = plt.plot(xvals, m*K(xvals, g_mu, g_std))" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 11, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "def sample():\n", | |
| " candidate = np.random.normal(g_mu, g_std)\n", | |
| " u = np.random.uniform()\n", | |
| " if u <= f_n(candidate) / (m * K(candidate, g_mu, g_std)):\n", | |
| " return candidate\n", | |
| " return sample()" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 12, | |
| "metadata": {}, | |
| "outputs": [], | |
| "source": [ | |
| "# draw 2000 new samples from estimated density\n", | |
| "s = []\n", | |
| "for k in range(2000):\n", | |
| " s.append(sample())" | |
| ] | |
| }, | |
| { | |
| "cell_type": "code", | |
| "execution_count": 14, | |
| "metadata": {}, | |
| "outputs": [ | |
| { | |
| "data": { | |
| "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYYAAAD8CAYAAABzTgP2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl8VOXd///XJyELshP2RUAElCighE0B1wouFWm1otbbVqq97+qv7W03l9a2trS1tbZ3b+2tfMWlVgWqUhARN1SwZQubEgQJm2DYtxACWa/fH+cEZmKWCZnJmUnez8djHnPmzFneQ0I+c67rnOuYcw4REZEKSUEHEBGR+KLCICIiYVQYREQkjAqDiIiEUWEQEZEwKgwiIhJGhUFERMKoMIiISBgVBhERCdMs6AB10aFDB9e7d++gY4iIJIwVK1bsc851rMs6CVUYevfuTXZ2dtAxREQShpltq+s6akoSEZEwKgwiIhJGhUFERMKoMIiISBgVBhERCaPCICIiYVQYREQkjAqDiIiEUWEQEZEwCXXls4g0HYs37T8xPapvRoBJmh4dMYiISBgdMYhI3Ag9SpDg6IhBRETCqDCIiEgYFQYREQmjwiAiImFUGEREJIwKg4iIhFFhEBGRMCoMIiISRoVBRETCqDCIiEgYFQYREQmjwiAiImFUGEREJIwKg4iIhNGw2yISOA23HV90xCAiImEiKgxmNt7MNphZrpndW8X7aWY2w39/qZn1DnnvPn/+BjMbFzJ/q5l9bGarzSw7Gh9GRETqr9amJDNLBh4HvgTsAJab2Rzn3LqQxSYDB51zZ5rZJOBh4EYzGwhMAjKBbsA7ZtbfOVfmr3eJc25fFD+PiIjUUyRHDMOBXOfcZudcMTAdmFBpmQnAc/70y8BlZmb+/OnOuSLn3BYg19+eiIjEqUgKQ3dge8jrHf68KpdxzpUCh4GMWtZ1wFtmtsLM7qx7dBERiYVIzkqyKua5CJepad0LnXN5ZtYJeNvM1jvnFn5h517RuBPg9NNPjyCuiIjURyRHDDuAniGvewB51S1jZs2ANsCBmtZ1zlU87wFmUU0Tk3NuqnMuyzmX1bFjxwjiiohIfURSGJYD/cysj5ml4nUmz6m0zBzgNn/6emCBc8758yf5Zy31AfoBy8yshZm1AjCzFsAVwNr6fxwREamvWpuSnHOlZnY38CaQDDztnMsxs4eAbOfcHGAa8LyZ5eIdKUzy180xs5nAOqAUuMs5V2ZmnYFZXv80zYAXnXPzY/D5RESkjsz7Yp8YsrKyXHa2LnkQaWxqu/J5VN+MBkrS+JjZCudcVl3W0ZXPIiISRoVBRETCqDCIiEgYja4qIg1vy6KT033GBJdDqqTCICJxL7RzWh3RsaemJBERCaPCICIiYdSUJCIxVV0zUE5ePgD55bp7W7zREYOIiIRRYRARkTAqDCIiEkaFQUREwqgwiIhIGBUGEREJo8IgIiJhVBhERCSMCoOIiIRRYRARkTAqDCIiEkZjJYlI1NV2D+cKVlRMtyf/QsdXp5O6eyfH+vZn5ze/zf6rJ4JZjFNKdXTEICKBSD54mN53/4Ref3iIko6d2DvxRlLy99L/+3dy1n/fAuXlQUdssnTEICIN79hxTv/xL0nftI0Njz/LgXHXAHAg71o6/b+/0/H5mfTq+CDbHvh1wEGbJhUGEWl4D/6J5hty2T7l/hNFAYCkJPZ8+z9IKiyk2zNPcPiCsRy65IrgcjZRakoSkYb1+uvwyhvsu+V6jowZWeUiu++azNGzz6Hv/d8n6WhBAwcUFQYRiYrFm/afeFTr+HG46y7o34e937yp2sVcagqbH3qE1L176P7//jcGaaUmKgwi0nAefxy2bYOffReXmlLjogXnZbHv6uvo+tRfabZ/bwMFFFBhEJEGknS0AKZMgXHj4MKhEa2z/Xs/IanoOF2ffyrG6SSUCoOINIhO/3gBDh6EX/wi4nWOn9GPg5eNp/Pfnyap8GjswkmYiAqDmY03sw1mlmtm91bxfpqZzfDfX2pmvUPeu8+fv8HMxlVaL9nMVpnZ3Pp+EBGJY6WldH3mCRg9GkZW3eFcnbzJ3yHl0EE6vPZqjMJJZbUWBjNLBh4HrgQGAjeZ2cBKi00GDjrnzgT+BDzsrzsQmARkAuOBv/rbq/A94JP6fggRiW8Zb8wm/fPt8KMf1XndI1kjKTxzgHfEIQ0ikiOG4UCuc26zc64YmA5MqLTMBOA5f/pl4DIzM3/+dOdckXNuC5Drbw8z6wFcDajxUKSR6/LiMxzr1Qeuuab2hSszY88Nt9BqdTbNP9X3yIYQSWHoDmwPeb3Dn1flMs65UuAwkFHLun8GfgzouneRRix980ZaL1/CnhtvhaRT69bce93XKE9JodOr06OcTqoSyU+pqpGsXITLVDnfzK4B9jjnVtS6c7M7zSzbzLL37tUpayKJptM/XqC8WTP2TrzxlLdRmtGBwxdeRPs35rA4d1/t10tIvURSGHYAPUNe9wDyqlvGzJoBbYADNax7IXCtmW3Fa5q61Mz+XtXOnXNTnXNZzrmsjh07RhBXROKFFRfT8dUZHLx0HCUdO9drW/uvuo70z7fT8qNVUUon1YmkMCwH+plZHzNLxetMnlNpmTnAbf709cAC55zz50/yz1rqA/QDljnn7nPO9XDO9fa3t8A59/UofB4RiSNtPnyP1P172XP9zfXe1oHLr6Q8JYWMef+MQjKpSa2Fwe8zuBt4E+8MopnOuRwze8jMrvUXmwZkmFkucA9wr79uDjATWAfMB+5yzpVF/2OISDzqMHcWJW3bcXj0JfXeVlnrNhwafQkZ82aDq9yaLdEU0eiqzrl5wLxK8x4MmT4O3FDNulOAKTVs+33g/UhyiEjiSDp+jPbvvsG+a76CS02Nyjb3XzWB9u+9RYucNRw9Z0hUtilfpGG3RSQm2r7/DslHj7LvmokAtN61BJJan9K2Wu9aAkDZWW1xZnR57Rn2drgJ+l4dtbxykobEEJGY6DB3FsUdOpE//MKobbOsXRuODexPy38vj9o25Yt0xCAiUZdUcIR2773F7q99HZJPDnaQk5df720XjMqi47QXST54qN7bkqrpiEFEoq7twgUkFR1n/5XX1r5wHR0ZNQxzjpZLar0MSk6RCoOIRF37d+dT0q49R4aOiPq2j/fvS0lGe1qpOSlmVBhEJKqspIR277/NwUuuCGtGit4OjILh59Fi5UdQrhF1YkGFQUSiqtWKpTQ7fIiDl42P2T6ODh1Ms8P5sGZNzPbRlKkwiEhUtXvnDcpT0zg0+uKY7ePo0EHexIIFMdtHU6bCICLR4xzt353P4QvGUN6iZcx2U9qxA0Wn94DX/gFbFnkPiRoVBhGJmuYb15O+fRsHYtiMVOHo0MGw/CMoLon5vpoaFQYRiZr278wH4OCl42pZsv4KsgZD4TFYo5v3RJsKg4hETdsP3qEgcxAlnbvGfF+F550LZvBvXc8QbSoMIhIVyfmHabU6m0MXXd4g+ytr3Qoy+8ES7/4MFTfv0Q186k9DYohI5EI7efuMCXurzb8XYmVlHBpT/yG2IzZsMLw4G4qKIaXhdtvY6YhBRKKi7aIFlLZsRcGQrIbb6bBBXlHI+bTh9tkEqDCISP05R9uFCzh8wVhcSgN+dR96rvec/XHD7bMJUGEQkfpbv560nZ9zaOylDbvfju2hdw9Yriugo0l9DCJSf/O901QPjbn0xE116iuS7eTk5dMt82xaLVrijZuUpO+60aB/RRGpvzffpLBvP4q792zwXRcOGkiz/CM0z1U/Q7SoMIhI/Rw7Bh98wOExDdyM5CscnAlAqxXROVIRNSWJSB1V3IUtv3w/o/pmwMKFcPw4h06xMNS36am4e1dK2reldfZS9tz0jXptSzw6YhCR+pk/H9LTyR9xQTD7N6Nw0EBaZeuIIVpUGETklLTetcS74G3uLBg7lvL05oFlKRyUSfrn20ndmRdYhsZEhUFETt3nuyB3G4yP/WiqNSk892wAWq7ODjRHY6HCICKnbuEy73lc7EdTrUnRmX0oT02jlQpDVKjzWURO3cJl0LUTnH02bD4QWAyXksLxAWfQZtn7tN61hMWMPPHeqL4ZgeVKVDpiEJFTU1oG/1oBY4d7w18HrHDgAJpvyMVKdOOe+lJhEJFTclrOeig4CmOGBx0FgGOZZ5FUXEzapq1BR0l4ERUGMxtvZhvMLNfM7q3i/TQzm+G/v9TMeoe8d58/f4OZjfPnpZvZMjNbY2Y5ZvbLaH0gEWkYLZet9IaguHBo0FEAKMwcAPgFS+ql1sJgZsnA48CVwEDgJjMbWGmxycBB59yZwJ+Ah/11BwKTgExgPPBXf3tFwKXOucHAEGC8mY1ERBJGi2UrYcjZ0KZV0FEAKO3UgZIO7WmesyHoKAkvkiOG4UCuc26zc64YmA5MqLTMBOA5f/pl4DIzM3/+dOdckXNuC5ALDHeeAn/5FP/h6vlZRKSBJB/Op/n63LhpRgLAjGOZZ6kwREEkhaE7sD3k9Q5/XpXLOOdKgcNARk3rmlmyma0G9gBvO+eWVrVzM7vTzLLNLHvv3r0RxBWRWGuRvRpzzut4jiOFmQNI+3wnzfbvCzpKQoukMFR1ukHlb/fVLVPtus65MufcEKAHMNzMzqlq5865qc65LOdcVseOHSOIKyKx1nLpSkpbtYRBZwUdJcwxv5+h5ZoVASdJbJEUhh1A6Fi6PYDK152fWMbMmgFtgAORrOucOwS8j9cHISLxzjlaLlvJ0awhkJwcdJowxwaciUtOotVqFYb6iKQwLAf6mVkfM0vF60yeU2mZOcBt/vT1wALnnPPnT/LPWuoD9AOWmVlHM2sLYGbNgcsBnUogkgDStmwjZd8BCkacH3SUL3Dp6Rzv20dDY9RTrVc+O+dKzexu4E0gGXjaOZdjZg8B2c65OcA04Hkzy8U7Upjkr5tjZjOBdUApcJdzrszMugLP+WcoJQEznXNzY/EBRSS6Wi5dCUDBsPPChuCOF8cyz6L1Wx9AWVncHdEkioiGxHDOzQPmVZr3YMj0ceCGatadAkypNO8j4Ly6hhWR4LVctorjvXtS2jk++/wKMwfQftbrNN/0Kcf6nx10nISksZJEJHLHjnPamrUcvO6qsNnRus9zNFR0QLdauVyF4RRpSAwRidyyNSQVl1AwPP76FyoU9+hGSdt2OjOpHlQYRCRyi5ZTnprC0SGZQSepnhkFg4fSas3KoJMkLBUGEYncomUUDj4Hl54edJIaFQw+n+Yb15NUcCToKAlJhUFEIvPZZ7BxKwXD4/+8kSNDsjDnaPnxqqCjJCQVBhGJzOuvA3DkgmEBB6ldwWCvD6TlajUnnQoVBhGJzNy50Ks7xaf3CDpJrcratOVYn7661ecpUmEQkdodPQrvvguXXhAXd2uLxJEhWbRcsxKcBm6uKxUGEandggVQVOQVhgRRMHgoqfv2eH0jUicqDCJSu7lzoVUrGDYo6CQRKxjiX2uxJH4uvksUuvJZmrYti05O9xkTXI5Yq8/ndM4rDOPGQWoKcCyq0WKlcEAm5WnpJC1dCjfeGHSchKIjBhGp2apVkJcH11wTdJI6cSkpFJwzWEcMp0BHDCJVacgjiXg/anntNa/D+cor4Whi3DazYuym4n7dYNY8KC6G1NSAUyUOHTGISM1eeQUuvBA6dQo6SZ0dGzjA6zT/6KOgoyQUFQYRqd7GjfDxx/DVrwad5JQcG9jfm1BzUp2oMIjUxZZFJx9NIcMrr3jPX/lK7PcVAyWdO0LXrrB0adBREooKg4hU75VXYNgwOP30oJOcGjMYMUJHDHWkwiAiX7RlEXz4D8jOTthmpBNGjIDcXNgfP7cfjXcqDCJStfkLvedELwwjR3rPy5YFmyOBqDCISNVmvwXn9Iczzww6Sf1kZUFSkpqT6kCFQUS+aONWWPspXHdF0Enqr2VLOOccdUDXgQqDiHzRrDchORm+fHnQSaJjxAivKam8POgkCUGFQSSBLN60n5y8fHLy8mO3k/JymP02jM6Cju1jt5+GNGIEHDzoXZchtVJhEJFwCxfCzj0wcVzQSaKnogNazUkRUWEQkXBPPw0tT4MvjfZex8NFffV11lnesOHqgI6ICoOInLRvH8yc6R0tNE8POk30JCd7F+rpiCEiKgwictIzz3iDzt1yXdBJom/kSG8wvcLCoJPEPQ27LVJPizedvKJ2VN+MiJavGBY6s1vrmOWqs/JyeOIJGDsW+vcJOk30jRgBpaWwciWMHh10mrgW0RGDmY03sw1mlmtm91bxfpqZzfDfX2pmvUPeu8+fv8HMxvnzeprZe2b2iZnlmNn3ovWBROQUzZkDmzfDd74TdJLYGDHCe1ZzUq1qLQxmlgw8DlwJDARuMrOBlRabDBx0zp0J/Al42F93IDAJyATGA3/1t1cK/MA5dzYwEririm2KSCxU1ZnsHPz2t3DGGYk/BEZ1OneG3r3VAR2BSI4YhgO5zrnNzrliYDowodIyE4Dn/OmXgcvMzPz5051zRc65LUAuMNw5t9M5txLAOXcE+AToXv+PI1KzxZv2n3hIiPfe8y4A+9GPoFkjbmEeMUJHDBGIpDB0B7aHvN7BF/+In1jGOVcKHAYyIlnXb3Y6D6jyp2Vmd5pZtpll7927N4K4IlInzsGvfgVdusA3vhF0mtgaORK2b/fuYS3ViuSrgVUxz0W4TI3rmllL4BXg+865Ki/ldM5NBaYCZGVlVd6vSFTk5OWTX+4dRVTXgex1Gp/8Na1Y/gv8JpqcvHzyu4yscZuhVzBX1xEd86ObefPg/ffhf/4H0hvRKapVqbjQbfHixttkFgWRHDHsAHqGvO4BVC63J5Yxs2ZAG+BATeuaWQpeUXjBOffqqYQXkXoqLfWaj/r1g//8z6DTxN7553vF78MPg04S1yI5YlgO9DOzPsDneJ3JN1daZg5wG7AYuB5Y4JxzZjYHeNHMHgW6Af2AZX7/wzTgE+fco9H5KCJxzj+SCD3qCNxjj8Enn8CsWfB5I217r3zF9ogRsCiBr+JuALUeMfh9BncDb+J1Es90zuWY2UNmdq2/2DQgw8xygXuAe/11c4CZwDpgPnCXc64MuBC4FbjUzFb7j6ui/NlE6i0nL7/xdlZv3g733w9XXw0TKp9P0oiNGQOrVsGRI0EniVsRnX7gnJsHzKs078GQ6ePADdWsOwWYUmneh1Td/yAiDaGkFH78G0hLg6lTvXsjNxWjR3sX8y1dCpc3kmHFo6wRn5cmEh0VVykD4BzN9u0nZecekj8+QMr+fSQdK8TKyrCyMmiXDgU7IDWVtsfLSe24nfL05tC7IxzZDK1aklJYTlmLFpS3aO6N4eOLpAM8an71F1iZAy+9xOJjaeB3rId2gEfSMZ6QRo3y7ui2aJEKQzVUGERqkHSkgBar19Ji1cc0z1lP2tbtJB+NbKyd6i7M6R8yXdY8nfIWp1HWsgWlGe0o7J1JUdfucP5AOPdcyMys70f4oseeg7//E+6YBJMmQWNsJqtJ69YweLA6oGugwiDiO3FkUJIGb3wAcxcwYOEykkpLKU9N5djA/hwadwlFvU+nuHsX8vuOpqR9B8patMAlJ0NyMiPO7AgbP4DjRXz62X6OthlIUlERgzNS4dN/Q8FRPt+xn6SjhSQfLSSp4CjJhYUkHTlKyv4DtPnXB6Tu2XXyTmNJSZx3eg8Oj7iI/OEXkD/8wlNvhC0thd8/CU/N8EZP/fG3o/LvlpBGj4Zp06CkBFJSgk4Td1QYRHzN9h+k/atz4bX5cOAwdO3Egeu/zJHRIzg2cAAuNfwPyNEug7+4kfR0714GLU+jpDiZ413O9Ob3zYDOJQAcquHua/ldRmIlJYxMPgJvT4f1myhZsY6MebPpPON5AIq6d6XwvHPhshFw45nQtWvtA/mt2wg//xOsWAu3ToQHvxvWjNVUVDTXte83hAGFhV4n9PDhQceKOyoM0riEnprYZ8wX5oedKprkt5vnH6HTk8+RMXMOVlwMl18It30VRp7H7l0FDRA6nEtJgb79oNnFcOXFfJaXT37HYZy2PofWy/5NxsK5tH7/XzD3LfjvX8GAAfQZOoqjmYM53vsMKDrTG9Zi1y5vJNHXXoMFC6B9G3j0p3DdFQ3+meJJ611L4HS/yH/44cnCUN3vThOkwiCNXuUrlk8oL4cXZsOjT9Hx8BEOfeki9k6+heIe3bz3AygK1UpOpjBzkPe4cgiUlZF5aA98egjee48Os1+my4vPVr1u//7wg29591ho24g6keuhtEN76NXd64C+556g48QdFQZpktK2fAbf+6vXtDLqfDbdcRvH+/cNOlbkkpPh3AFw7Rj44Q9Z/uke0vJ2kL51MwObl3n9CZ06wYABcPrpiX1bzljJOhc++NAbK6opna4bARUGabyqaj4qLydjxj/pNPVv0OI0eOR+mDiO4zvrd7HTiY7rpIC+kScnU9SzF0U9e3n9GRWFoGwbbNkWTKZ4lzUIXpkP69bF5uyvBKbCIE1G8qHDdJ/yJ1otziZ/7ChaP3IvdGgXdKwqVdv8JdFzwfne83vvqTBUosIgTULzj9bR88HfkXw4n7x7/ouDE68is0OboGN9QdjFdBJbPbt5N+5ZsADuvjvoNHFFhUEajdBv2aFX6rad9w5df/8YJV068tkfHuV4vzNiliGnhlNRY7Gvaof+DlmmQlVXNde2fqN36aXeAIJlZUEniSsR3fNZJCGVlcHv/o/uv/kzhUPOYfPUP8W0KEgCuvRSOHgQ1qwJOklcUWGQxqnwGPzXT2HqSxy47iq2PfILylu3DDqVxJu+p3nPLz8VbI44o6YkaXSSCo7C96fAyrXwi++x8/IvBR0p4TVkE1mD6twB+p4O/14Jd9wUdJq4ocIgjUryocP0+sHPYdMW+MvP4apLIEZ/1Bqyo1id0jE06nx49U1vKPIU/UkENSVJI5KyZxe9776PtC2fse23PyVnyNDG+01XomfU+V7T40efBJ0kbqgwSOOwdSvnTPoyKbv3su2RX1AwaljQiSRRjDzPu/L5XyuCThI3dNwkiW/DBrj8cprlH2Hbn3/NscwBQSdqMFU2MdXh6ms1UQHt2njDi3ywFL77jaDTxAUdMUhiW7PGu4dvcTE5L85pUkVBouiikbB6HRw8HHSSuKDCIPFry6KTj6osWQIXX+zdA2HRIgrP0rAGcoouHuENprdoedBJ4oKakiQwtd5cpgY5f5/FWd++leSuXeDdd6FXrzrforKxdkw31s8VbWFXhQ86y2tS+mApXKv7QKswSOJ5/XXOnnwTx3v14bTnfwPln8GWzzTonJy65GQYPQwWLj15W9UmTE1JklhmzoTrrqNwwNnkvDgHOnUIOpE0FpeMhP2HYO2nQScJnI4YJK6FDvY26oPZcMcdcOGFrPvL3yhr1YqcPP0njpSamGoxZrh32uoHS2HCHUGnCZSOGCQYWxbReteSiE+X7PLcVJg8GS6/HObPp6xVqxgHlCYno6132ur7OoVXhUHim3N0fPYl+vzqfpg4EebMgdNOCzqVNFaXjPJOW929O+gkgVJhkPjlHJ0fe4pOT73Anok3ev0LaWl1PtoI0omcuudyYhg31jttdfbsoJMEKqLCYGbjzWyDmeWa2b1VvJ9mZjP895eaWe+Q9+7z528ws3Eh8582sz1mtjYaH0QambIyuPdhOsyYzf7rv8ymh/8XmqlLTGJswBnQqzu8+mrQSQJV6/80M0sGHge+BOwAlpvZHOfcupDFJgMHnXNnmtkk4GHgRjMbCEwCMoFuwDtm1t85VwY8CzwG/C2aH0gagaIiuOUWeGUee755E3tvvxmSEvvgNlE6fhMlZ8yYeUcNz7wMhw5B27ZBJwpEJP/bhgO5zrnNzrliYDowodIyE4Dn/OmXgcvMzPz5051zRc65LUCuvz2ccwuBA1H4DNKYHD0K114Lr7wCD9zN3sm3gNnJ5hg1yUisjRsLJSXwzB+b7O9bJIWhO7A95PUOf16VyzjnSoHDQEaE64qweNN+lq/IJX/spfDOOzBtGkz+WtCxpCkafLZ3A5+3mmZRgMgKg1Uxz0W4TCTr1rxzszvNLNvMsvfu3VuXVSWBpH6+ncwbr6blx6thxgy4/fagI0lTlZQEV4zxTlstKKx9zK5GKJLCsAPoGfK6B5BX3TJm1gxog9dMFMm6NXLOTXXOZTnnsjp27FiXVSVBnLbuY869fhype3fzybP/gOuvDzqSNHVfvhyOF8HbTacYhIqkMCwH+plZHzNLxetMnlNpmTnAbf709cAC55zz50/yz1rqA/QDlkUnujQGLZat5JybrsE1S2Ht9NfJH3Fh0JGkCcrJyz/xAGDoOdCjC8x+O9hgAam1MPh9BncDbwKfADOdczlm9pCZXesvNg3IMLNc4B7gXn/dHGAmsA6YD9zln5GEmb0ELAYGmNkOM5sc3Y8m8a7NG+/S60e/5HjP3qz9x3yO9T8r6EgiHjO49kvwYTbsbXrnyJj3xT4xZGVluezs7KBjSH05Bz+8Ax6dRsHQwax7ahZlrb5417FEuIBNGp/Mbv7v4satMO4/4Gf/H3zzBm9enzGB5TpVZrbCOZdVl3US++RwSTzHjnnXKDw6jUPjLuGzR35RZVEQCVy/3jCwH/yz6TUnqTBIw8nLg4sugpdegh/eyec/vQeXkhJ0KpHqTRwHH6+HTzYFnaRBqTBIw8jOhmHDYN06+Oc/4Ttf99pxReLY+gsuoDw1hQNPvRx0lAalwiCxN2MGjBnjjXX0r3/BhPAL5ysGmlOfgsSbsjatyb9kNG3eXABHC4OO02BUGCQmFm/az5L1u8j71l0waRIMHQrLl8PgwUFHE6mTg9eOJ7nwGMxdEHSUBqPhKiX6tiyi/dqt9Pj5w7RYkwNfvw4euBs6dap1VR01SLwpHDSQ431OJ/2F2fCT3zWJJlAdMUj0LVvNGbd/l+Ybctnx4A/goXsgLTXoVCKnxowDE6+GtRu8ptAmQEcMUmeLN+2vcv6oPu3gj3+E++6jvFsXtv351xSd0ZseFQtsWaRhnSUhHbrqMro98yL8/vcwenTQcWJORwwSFSl7d8OVV8KPfwxXjGHzU3+i6IzeQcdAPbPGAAANAklEQVQSiQqXng7/8RV47TXvzLpGToVB6q3te28x+OqLYNEiePJJeOyXlLfQfZmlkbl1IjRvDo88EnSSmFNhkFNmRcfp/ct7OfuOmynu3AVWrIA772wSnXPSBLVvCzdcCX9/Hj6YEXSamFJhkFNy2rqPGTTxS3R9/inyvvFtPn55Ppx9dtCxRGLrv74Oycnw56eDThJT6nyWuikupseff0f3J/5Mabv2fDJtOocuupzWu5aQ86/VQacTqbcaT5Do3AG+cT08+SJ89BEMGtRwwRqQjhgkcqtWwbBh9HzsEfZf8xXWvPEvDl10edCpRBrWt2+GVi3g/vuDThIzKgxSu2PH4Gc/88Y62rOH9U/+ndxH/kpp23ZBJxNpeG1awXduhddf9x6NkAqDnFTVvW3feAPOOQd+/Wu4+WbIyeHgZeODyygSD75xPQwcCHfdBUePBp0m6lQYpGrbt8NXvwpXXQUpKfDOO/C3v0H79kEnEwleago88QRs2wYPPRR0mqhT53MjFHpl8qi+GXVad92mPXSYMYuMF2aRjIPf/AbuuQfS0k45j652lkZpzBiYPNm72n/CBLjgAm9+6BF3At7xDVQYpEJpKUx/jX5/nEbK/gPsH3cNGU/8L/TuHXQykfj1xz/CggVeM+vq1dC2bdCJokJNSU1daSkbH32CwrMy4f4/UNKtM5v/7w98+vizLC5rxeJN+6sdG0mkyWvTBqZPh88/h299y7ufeSOgI4ZEEXp4CvU/RC0q8voMfvc7+m3eTGH/s/lsyv0cGTuq5iuXtyyi9S6vaSi/y8j6ZRBpDIYPhx/dAb/9P/je7fDftwedqN5UGJqabdtg6lR46inYsweGDWP9j37BwcvG03rPsqDTiSSmb02C3G3wv89C987wtauDTlQvKgxBiVUHVVXf6AsLYd48eO457xng0lFw60/g1v/m4OYDtW62ojmpYtsiEsIMfv1D2LUXHngEUprBxHFBpzplKgyNlBUV0e7d+fDzN2D2bCgogK5dvas1xw32vtWABrwTiZaUZvDXX8Gd98MPpsDRY/AznZUkDanSEcfiTftJ37qJLnNf4PQl2bRYtZakoiLvuoObb/buuzx2rDcAWMgNc/LLq+5Y1i02RSKzeNP+E0fSmd1aw9MPw10/hwcfhd1F3plLaWkJdRqrCkOiKi+HjVthxcewfirnvfcB6Z9vB6CoRzcOfvkKdl91KwO/fp13gZqINIy0NPi/X8MjU+Hxx2H5cnjmGWgedLDIqTAkAue8tsucT2Htp/DReq8g5Bd473fuzNEhw9k5+S5KMjMo6d4V8PsYVBREGl5KM7jvOzDueu8eJUOGwJ2T4E5/AL44p8IQsJy8/LDmHCsqIn3bZoYcyfMumFm1Clau9M4gApwZxad3p/CiC2g3dihkDYKLbuRTvwM5tAmo9a4lkNS6xv3Xp8lIzU3SlIRewZ/ZrfWJpqEaT8g4vxO8+SxMeQwefx7+/k+Y/DX4UT/o0sVbJg6bmCIqDGY2HvgfIBl4yjn3u0rvpwF/A4YC+4EbnXNb/ffuAyYDZcB3nXNvRrLNRu/wEfgsjzbZ62mz9x1O2/QpzXM3kP7ZVqyszFsmORkyM73xinq1gcz+rG/bifLTvGPSdt38P/rqQBaJXxlt4dGfwjdvgL88A49Og/95Fq65xuv7G9AK2tb8Ba6hmavlSj0zSwY+Bb4E7ACWAzc559aFLPMdYJBz7j/NbBIw0Tl3o5kNBF4ChgPdgHeA/v5qNW6zKllZWS47O7vunxJOVOWcvPwTp3HWdRyhiJWVwb59sHOn99i6lc9XraPlxpWk5u2m+a7dJ5uBAJecTFHPbhT36klR754U+c97R3wVl5bu5QzJLyKJJbNbyB/+TdvgrY+8C0x374akJDh3APvP7MuBkVdSOPBcBl82HPJCriuqx5GEma1wzmXVZZ1IjhiGA7nOuc3+TqYDE4DQP+ITgF/40y8Dj5mZ+fOnO+eKgC1mlutvjwi2GX1bttN8407K85K8b+U7Wnl/xMvKvM7cSJ7LyuDIEcjPh8OHTz4fPuz9kHft8pp9Kr71+7qmpFLStRPF3TrDiEHQsyv07MbGVhkU9+gKzb74o3Bp6TH95xCRAPTtBb//Ovz2t17H9ItPwpKVtJ33DhmvzPWWue02+PnkwCJGUhi6A9tDXu8ARlS3jHOu1MwOAxn+/CWV1u3uT9e2zej7/ZOc8ebC6G0vPd0bK6V1a++5a1c4/3yv7bBrV++5Sxfo3ZulhaknriwO/fZQrCMAkaYpORlGjoTOJcDtrN9+kJKjGZy2IYf+Q4O9f3okhaGqBuzK7U/VLVPd/KoG76uyTcvM7gTu9F8WmNmGanICdAD21fB+dB0/7j127450jYbNV3fxng/iP2O854P4z6h89ReasVddV46kMOwAeoa87gHkVbPMDjNrBrQBDtSybm3bBMA5NxWYGkFOzCy7rm1pDUn56i/eM8Z7Poj/jMpXf/XNGMmw28uBfmbWx8xSgUnAnErLzAFu86evBxY4r1d7DjDJzNLMrA/QD1gW4TZFRCQAtR4x+H0GdwNv4p1a+rRzLsfMHgKynXNzgGnA837n8gG8P/T4y83E61QuBe5yzpUBVLXN6H88ERGpq4iuY3DOzQPmVZr3YMj0ceCGatadAkyJZJtREFGTU4CUr/7iPWO854P4z6h89VevjLVexyAiIk2Lbu0pIiJhErIwmNmvzOwjM1ttZm+ZWTd/vpnZX8ws13///JB1bjOzjf7jtuq3HrWMfzCz9X6OWWbWNuS9+/yMG8xsXMj88f68XDO7N8b5bjCzHDMrN7OsSu8Fnq+KvIHtu1KOp81sj5mtDZnX3sze9n+33jazdv78an8fY5ivp5m9Z2af+D/f78VTRjNLN7NlZrbGz/dLf34fM1vq55vhn5SCf+LKDD/fUjPrHct8ITmTzWyVmc2N03xbzexj/29gtj8vej9j51zCPYDWIdPfBZ7wp68C3sC7fmIksNSf3x7Y7D+386fbxTjjFUAzf/ph4GF/eiCwBkgD+gCb8Drgk/3pM4BUf5mBMcx3NjAAeB/ICpkfF/kqZQ1s31VkGQucD6wNmfd74F5/+t6Qn3WVv48xztcVON+fboU39MzAeMno76elP50CLPX3OxOY5M9/Avgvf/o7If+/JwEzGujnfA/wIjDXfx1v+bYCHSrNi9rPOCGPGJxzoZcLt+DkxXETgL85zxKgrZl1BcYBbzvnDjjnDgJvA+NjnPEt51yp/3IJ3rUaFRmnO+eKnHNbgIphQk4MPeKcKwYqhgmJVb5PnHNVXSwYF/kqCXLfYZxzC/HOvAs1AXjOn34OuC5kflW/j7HMt9M5t9KfPgJ8gjfaQFxk9PdTMVBYiv9wwKV4w+lUla8i98vAZWaxHTXSzHoAVwNP+a8tnvLVIGo/44QsDABmNsXMtgO3ABVnSFU1fEf3GuY3lNvxKjY1ZAk6Y4V4zBcv/zbV6eyc2wneH2agkz8/0Nx+s8Z5eN/K4yaj30yzGtiD9yVtE3Ao5ItUaIaw4XaAiuF2YunPwI+Bcv91RpzlA6+YvmVmK8wbHQKi+DOO2/sxmNk7QJcq3nrAOTfbOfcA8IB5w3rfDfycug/NEdOM/jIP4F3D8ULFatVkiXiYkGjmq2q1anJEPV8dxOTn1wACy21mLYFXgO875/Jr+BLb4Bmddy3TEPP63WbhNWtWl6FB85nZNcAe59wKM7s4ggxB/YwvdM7lmVkn4G0zW1/DsnXOGLeFwTl3eYSLvgi8jlcYqhuCYwdwcaX578c6o3md3NcAlzm/sa+GjNQwPyb5qtFg+aKUKR7sNrOuzrmd/iH6Hn9+ILnNLAWvKLzgnHs1HjMCOOcOmdn7eO3ebc2smf+tOzRDdcPtxMqFwLVmdhWQDrTGO4KIl3wAOOfy/Oc9ZjYLr7k1aj/jhGxKMrN+IS+vBSqq5RzgP/xe+JHAYf+Q6k3gCjNr5/fUX+HPi2XG8cBPgGudc4Uhb8X7MCHxmC9e/m2qEzokzG3A7JD5Vf0+xozfvj0N+MQ592i8ZTSzjv6RAmbWHLgcrx/kPbzhdKrKV9VwOzHhnLvPOdfDOdcb7/dsgXPulnjJB2BmLcysVcU03t+ztUTzZxzLnvNYPfC+Da0FPgJeA7q7k2c8PI7XZvkx4Wfb3I7XkZoLfLMBMubiteut9h9PhLz3gJ9xA3BlyPyr8M4i2YTX3BPLfBPxvkkUAbuBN+MpXxV5A9t3pRwvATuBEv/fbzJem/K7wEb/uX1tv48xzDcar5ngo5DfvaviJSMwCFjl51sLPOjPPwPvC0gu8A8gzZ+f7r/O9d8/owF/1hdz8qykuMnnZ1njP3Iq/j9E82esK59FRCRMQjYliYhI7KgwiIhIGBUGEREJo8IgIiJhVBhERCSMCoOIiIRRYRARkTAqDCIiEub/B5FXL4LjTogSAAAAAElFTkSuQmCC\n", | |
| "text/plain": [ | |
| "<Figure size 432x288 with 1 Axes>" | |
| ] | |
| }, | |
| "metadata": { | |
| "needs_background": "light" | |
| }, | |
| "output_type": "display_data" | |
| } | |
| ], | |
| "source": [ | |
| "plt.hist(X, bins=100, alpha=0.25, density=True)\n", | |
| "plt.hist(s, bins=100, alpha=0.25, density=True)\n", | |
| "plt.plot(xvals, f_n(xvals), color='red');" | |
| ] | |
| } | |
| ], | |
| "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.7.3" | |
| } | |
| }, | |
| "nbformat": 4, | |
| "nbformat_minor": 2 | |
| } |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment