{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Initial Conditions (Experiences, Observables, Lagged Choices)\n", "\n", "While working with **respy**, you often want to simulate the effects of policies in counterfactual environments. At the start of the simulation, you need to sample the characteristics of individuals. You normally want to do at least one of the following points:\n", "\n", "- Individuals should start with nonzero years of experience for some choice.\n", "- The previous (lagged) choice in the first period can only be a subset of all choices in the model.\n", "- An observed characteristic is not evenly distributed in the population.\n", "\n", "Taken together these assumptions are called the initial conditions of a model. An initial condition is also called a *seed value* and determines the value of a variable in the first period of a dynamic system. \n", "\n", "In the following, we describe how to set the initial condition for each of the the three points for a small Robinson Crusoe Economy. A more thorough presentation of a similar model can be found [here](../tutorials/robinson_crusoe.ipynb)." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "\n", "import io\n", "import pandas as pd\n", "import respy as rp\n", "import matplotlib.pyplot as plt\n", "import seaborn as sns" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Experiences\n", "\n", "In a nutshell, Robinson can choose between fishing and staying in the hammock every period. He can accumulate experience in fishing which makes him more productive. To describe such a simple model, we write the following parameterization and simulate data with ten periods." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "params = \"\"\"\n", "category,name,value\n", "delta,delta,0.95\n", "wage_fishing,exp_fishing,0.01\n", "nonpec_hammock,constant,1\n", "shocks_sdcorr,sd_fishing,1\n", "shocks_sdcorr,sd_hammock,1\n", "shocks_sdcorr,corr_hammock_fishing,0\n", "\"\"\"\n", "\n", "options = {\n", " \"n_periods\": 10,\n", " \"simulation_agents\": 1_000,\n", " \"covariates\": {\"constant\": \"1\"}\n", "}" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
value
categoryname
deltadelta0.95
wage_fishingexp_fishing0.01
nonpec_hammockconstant1.00
shocks_sdcorrsd_fishing1.00
sd_hammock1.00
corr_hammock_fishing0.00
\n", "
" ], "text/plain": [ " value\n", "category name \n", "delta delta 0.95\n", "wage_fishing exp_fishing 0.01\n", "nonpec_hammock constant 1.00\n", "shocks_sdcorr sd_fishing 1.00\n", " sd_hammock 1.00\n", " corr_hammock_fishing 0.00" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "params = pd.read_csv(io.StringIO(params), index_col=[\"category\", \"name\"])\n", "params" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "simulate = rp.get_simulate_func(params, options)\n", "df = simulate(params)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAu8AAAFwCAYAAAAIW5VaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deXxV9ZnH8e+ThH0ViCiBEJYEkoBQQXBDGKwOaN0Qt7YqWgcdsXWZTms7tcVqW3VqF62ttdRB24o61LZanSpqXbDFIlYqS1AUgbAoCLJvSZ7545zQ6zXLDbnJuSf5vF+v+yL3nuX3nJPwe577u79zrrm7AAAAAGS+rKgDAAAAAJAaincAAAAgJijeAQAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmKN5xSMxsppn9+hC3vdfMbkp3TOliZgVm5maWc4jbu5kNrmXZ58zsmZrWre+8mNnXzWzWocQEILOY2TQzmx91HDUxs95m9pKZ7TCzO6OOpy5mttPMBkYdR1Qy+e+ooczs/8zs0kPc9j0z+3S6Y8pUFO+olZl91sxeCzvHDeF/rBMbu193v8rdb0lHjNXCNxMHwlg/MrO/mNlx6WwjHdz9N+5+ai3LDp4XM5tgZuVJy7/r7lc0R5wAGs/MTgz7om1mtsXMXjGzY6KOKwXTJW2W1NXd/yN5oZnNNrP9YX9b/Vjc/GFK7t7Z3d+Nou2atKYiMsxTVeHvf4eZrTCzyw51f+4+2d0fSGeMLRXFO2pkZjdI+pGk70rqLSlf0k8lnRVlXPV4xN07S8qVNF/SY2ZmySsd6og6AKTKzLpK+qOkuyX1kJQn6WZJ+5qgrXT3af0lLfO6v8XxjrBwrn6MSHMMdaIfb151nO/1Yd7tKumrkn5hZiUN3LeZGfVoA3Cy8Alm1k3StyXNcPfH3H2Xux9w9yfc/T8TVm1rZg+G77iXmtnohH0Um9kL4Sj4UjM7M2HZbDO7NeH5WWb2hpltN7N3zGxSdRxm9stw1H+dmd1qZtn1xe/uByQ9IOkIST3DjxVfMbMfmtkWSTPNLMvMvmFmq83sg/A4uiXt6nIzWx+2f3D0yczGmNlfw2PbYGY/MbO2SdueZmbvmtlmM/vv6o6pro84q8+LmXWS9H+S+iSMavWxpKlKZnZsOKr3kZktNrMJCcumhe3vMLNVZva5+s4bgLQqkiR3n+Pule6+x92fcfd/JK5kZt83s63h/9PJCa9fZmbLw//D75rZlQnLJphZuZl91cw2Svqf8PXPhH1p9aePR9UWnJkdb2YLw08FFprZ8eHrsyVdKukrYd/ToFFkM7sgjLdr+HyymW00s9zwuZvZl2rqH8Pll4fHvdXMnjaz/gnL3MxmmNnbkt5OeK166mG78HyuMbP3LZiK2CHpnP1H2OdvsIRRYjPrYGZ3hjlhm5nNT9i21r62geemxt+Pmd1oZnOT1v2xmd0V/nxIudD+OQV0ei25LCts+x0z+9DMHjWzHknbfsHM1kh6vq62PPB7SVsllYT7qCtHvWBm3zGzVyTtljQwfO2KhNhqzdFmdnG47EMz+6/6zkWL4+48eHzsIWmSpApJOXWsM1PSXkmnScqW9D1JC8JlbSStlPR1SW0lTZS0Q9KQcPlsSbeGP4+RtE3SKQreTOZJGhou+72kn0vqJOlwSX+TdGUd8fw6/LmdpP+WtDZ8Pi08ni9KypHUQdLlYYwDJXWW9JikX4XrF0hySXPCtodL2iTp0+HyUZKODfdVIGm5pOsSYnFJf1Yw2pYv6S1JVyTEMj9p3cE1nJcJksrrOMY8SR+G5z8rPH8fKvjUoZOk7Qnn+0hJpVH/XfHg0ZoeCkYiP1QwkDBZ0mFJy6dJOiDp38I+9N8lrZdk4fLTJQ2SZJLGKyhwjg6XTQj7tNvD/q6DpKMlfSBpbLi/SyW9J6ldDbH1UFBkXRz2YxeFz3uGyw/2RbUcW33LfxOu0zM8ps8kLKurfzxbQb9cHMb1DUl/Sdp2Xrhth4TXqvvQH0l6PFzeRdITkr6XdM6+rSBHnRae08PC5fdIeiHsW7MlHR+e21r72lqO/T2FuSLp9Vp/Pwo+6ditYJqSwuUbJB0bPq81FyoppyS1WaC6c9l1khZI6hvG8XNJc5K2fTDctkMN+5+gME+F5+YcBX/TQ+o7b+G5XiOpNPxdtwlfq/5bqCtHl0jaKemkMO4fhL/bT5z3lvqIPAAemfeQ9DlJG+tZZ6akZxOel0jaE/48TtJGSVkJy+dImhn+PFv/LFJ/LumHNey/t4KPlzskvHaRpD/XEc9+SR+FHeTzkkaFy6ZJWpO0/nOSrk54PiTsdKoLclf4JiJcfoekX9bS9nWSfpfw3CVNSnh+taTnEmJJR/H+1eqOLGH50woSQqfwPJyrGjpcHjx4NM9DQRE6W1J5WFw8Lql3uGyapJUJ63YM+4MjatnX7yVdG/48Iezv2ics/5mkW5K2WSFpfA37uljS35Je+6ukaeHPB/uiWmKZrWDw5qOExwMJy7srKMzelPTzpG3r6h//T9IXEpZlKShq+ydsO7GG/Q1W8CZnl6RBCcuOk7Qq4ZztUcKglIJccWzYzh5JI2o41lr72lrOzXuquXiv8/ejYKrnJeHPp0h6J/y5zlyo1Ir3GnOZgoGnkxOWHalP5sGBdfwdTJBUFf7+t0h6Q9KFqZw3BYX6t5OWv6B/Fu915ehvSno4YVknBf8fWk3xzpwx1ORDSb3MLMfdK+pYb2PCz7sltbdgXlwfBaPeVQnLVyt4J56sn6Snani9v4J34hvsn9PWsyStrSOeR93987UsS96uTxhTYnw5CjrKmrZZrWDUQmZWpOCd/mgFCTdH0qI62lsdtpdO/SWdZ2ZnJLzWRkGHvsvMLpD0ZUm/DD+W/A93L0tzDADq4O7LFRRXMrOhkn6tYHT4onCVjQnr7g77us7h+pMlfUvB9JssBX3Nmwm73+TuexOe95d0qZl9MeG1tqq570nu/6Ta++jafN/dv1HTAnf/yMz+V9INCgYRktXWP/aX9GP7+B1uLIxrdQ3bJspVcI4WJeQMUzCKXe3DpJy2W8H57iWpvaR3athvrX1tLXHUpr7fz0MK/i4elPTZ8Hn1dg3NhclqzGXhvn9nZom5ulK158GarHf3vjW8nsp5q2vfdeXoPonbhjnvw3ribFGY846a/FXBqMrZh7j9ekn97OMXoORLWlfDumsVfDRc0+v7JPVy9+7ho6u7lx5iTF5DjP2T4quQ9H7Ca/2Slq8Pf/6ZpDJJhe7eVcH0oOQLY2vb9lDjTbZWwahG94RHJ3e/TZLc/Wl3P0XBSEqZpF80sH0AaRS+eZ4taVh965pZO0m/lfR9BSP13RUMciT2M8l9xFpJ30nqEzq6+5wamkju/6Ta++gGM7ORCqY9zJF0Vw2r1NY/rlUwHSTxGDq4+18S1q+tb9ysYPS8NGHbbh5cTFmfzQpyXm25qNa+tgHq+/38r6QJZtZXwfSThxK2a2wurOt8T06Kqb27J/4d1JeLapPKeatr33Xl6A2Jx2RmHRVM0Wo1KN7xCe6+TcHHUveY2dlm1tHM2oQXHt2Rwi5eVfDx5VfC7SZIOkPSwzWs+0tJl5nZyeEFKnlmNtTdN0h6RtKdZtY1XDbIzMan5yg1R9L1ZjbAzDoruKvOI0mjMjeFx14q6TJJj4Svd1Ewp3xnOJr27zXs/z/N7DAz6yfp2oRtU/W+gottky+irfZrSWeY2b+aWbaZtQ8vyOprwT2az7Tgwtd9CuYGVjawfQCNYGZDw4sj+4bP+ykYWV2QwuZtFczl3SSpIhyFr/EWswl+IekqMxtrgU5mdrqZdalh3ackFVlwO+Cc8JO6EgV3x2kUM2uvoH/6uoJ+M8/Mrk5arbb+8V5JXwv73OoLNc9Lpd3wk95fSPqhmR0ebp9nZv+a4rb3S/qBBTcHyDaz48I3UbX2tXXssk24XvUjR/X8ftx9k4JpI/+jYKrP8vD1dOTC2nLZvZK+Y+FFwWaWa2bpuqPcoZy3RHXl6LmSPmPBrVjbKriOoVXVs63qYJE6d/+Bgo88v6EggayVdI2CeZf1bbtf0pkKLtLarOAWk5fUNG3D3f+moDP5oYILV1/UP99tX6IgiS1TcDHVXAUjyelwv6RfSXpJ0ioFoy5fTFrnRQUXzDyn4CPi6i9X+rKCjzV3KOiQayrM/6BgKs0bkp5U8CYlZeG5miPpXQuu1O+TtHytgtt2fl3//P38p4L/01mS/kPByMUWBRe7JSdPAE1rh4KLE181s10KivYlCv5v1sndd0j6kqRHFfR9n1UwX76ubV5TcPHrT8JtViqcslPDuh9K+kwYy4eSvqLgotLNKRxXteq70VQ/qrf9noLrdX7m7vskfV7SrWZWmLBtjf2ju/9OwUW4D5vZdgXna7JS91UFx70g3P5ZBXOlU/FlBdOSFiroN29XcN1WXX1tbZ5S8ClA9WNmir+fhyR9Wv8cda/W2FxYWy77sYK/q2fMbIeCv9GxDdhvrQ7xvCWqNUe7+1JJMxScpw0Kzkl5zbtpmaqvagcAAGhSZuYKphyujDqWls7MChQUvm3quX4NMcPIOwAAABATFO8AAABATDBtBgAAAIgJRt4BAACAmKB4BwAAAGIism9Y7dWrlxcUFETVPABklEWLFm1299yo44gSeQEAAnXlhMiK94KCAr322mtRNQ8AGcXMkr+uvtUhLwBAoK6cwLQZAAAAICYo3gEAAICYoHgHAAAAYoLiHQAAAIgJincAAAAgJijeAQAAgJigeAcAAABiot7i3czuN7MPzGxJLcvNzO4ys5Vm9g8zOzr9YQIAMgV5AQCik8rI+2xJk+pYPllSYfiYLulnjQ8LAJDBZou8AACRqLd4d/eXJG2pY5WzJD3ogQWSupvZkekKEACQWcgLABCdnDTsI0/S2oTn5eFrG5JXNLPpCkZhlJ+fX/deZ3ZrfGQztzVy+0bGEPf2MyGGqNvPhBiibj8TYoi6/UyJIT6aJi/sbuQ57Ni436GXL29c+5Ksb3Gjtq985bFGbZ99wpTGtf/ALY3aXpKyL72pUdtX3HRxo7bPueVXjdp+z/njGrV9h0dfbtT2klQ+cmijtu/7Rlmjtn+5d99GbT/u/fJGbS9J93fNbdT2l2/f1Kjtr7Kujdr+Xt9+yNum44JVq+E1r2lFd7/P3Ue7++jc3MaddABAxiIvAEATScfIe7mkfgnP+0pan4b9AgDiibwAoEUb06VdZG2nY+T9cUmXhHcXOFbSNnf/xEejAIBWg7wAAE2k3pF3M5sjaYKkXmZWLulbktpIkrvfK+kpSadJWilpt6TLmipYAED0yAsAolbYIbqR76jVW7y7+0X1LHdJM9IWEQAgo5EXACA66ZjzDgAAgFZiwJGdow6hVUvHnHcAAAAAzYDiHQAAAIgJincAAAAgJijeAQAAgJjgglUAAIAY6VnEtxG3Zoy8AwAAADFB8Q4AAADEBMU7AAAAEBMU7wAAAEBMcMEqAABAitoMKYg6BGSAwg7tImubkXcAAAAgJijeAQAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmKN4BAACAmKB4BwAAAGKC4h0AAACICb6kCQAAALEy4MjOUYcQGUbeAQAAgJigeAcAAABiguIdAAAAiAnmvAMAgNiwwUVRhwBEipF3AAAAICYo3gEAAICYoHgHAAAAYoLiHQAAAIgJincAAAAgJijeAQAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmKN4BAACAmMiJOgAAAADER8+i3KhDaNUYeQcAAABiguIdAAAAiAmKdwAAACAmKN4BAACAmKB4BwAAAGIipeLdzCaZ2QozW2lmN9awvJuZPWFmi81sqZldlv5QAQCZgJwAANGpt3g3s2xJ90iaLKlE0kVmVpK02gxJy9x9hKQJku40s7ZpjhUAEDFyAgBEK5WR9zGSVrr7u+6+X9LDks5KWscldTEzk9RZ0hZJFWmNFACQCcgJABChVIr3PElrE56Xh68l+omkYknrJb0p6Vp3r0pLhACATEJOAIAIpfINq1bDa570/F8lvSFpoqRBkuaZ2cvuvv1jOzKbLmm6JOXn5zc8WgBA1NKWEyTyQuwMLo06AqDVS2XkvVxSv4TnfRWMpiS6TNJjHlgpaZWkock7cvf73H20u4/OzeWrdQEghtKWEyTyAgA0VCrF+0JJhWY2ILzg6EJJjyets0bSyZJkZr0lDZH0bjoDBQBkBHICAESo3mkz7l5hZtdIelpStqT73X2pmV0VLr9X0i2SZpvZmwo+Uv2qu29uwrgBABEgJwBAtFKZ8y53f0rSU0mv3Zvw83pJp6Y3NABAJiInAEB0+IZVAAAAICZSGnkHAAAAEBhwZOfI2mbkHQAAAIgJRt4BAABipM2QgqhDQIQYeQcAAABiguIdAAAAiAmKdwAAACAmKN4BAACAmKB4BwAAAGKC4h0AAACICYp3AAAAICYo3gEAAICYoHgHAAAAYoLiHQAAAIgJincAAAAgJijeAQAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmcqIOAAAAAGiInkW5UYcQGUbeAQAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmKN4BAACAmKB4BwAAAGKCW0UCAACkyAYXRR0CWjlG3gEAAICYoHgHAAAAYoLiHQAAAIgJincAAAAgJijeAQAAgJigeAcAAABigltFAgAQE1n9i6MOAUDEGHkHAAAAYoLiHQAAAIgJincAAAAgJijeAQAAgJigeAcAAABiguIdAAAAiAluFQnUo2DvQ43a/r30hAEAAJDayLuZTTKzFWa20sxurGWdCWb2hpktNbMX0xsmACBTkBMAIDr1jrybWbakeySdIqlc0kIze9zdlyWs013STyVNcvc1ZnZ4UwWM1oVRbyCzkBMAIFqpjLyPkbTS3d919/2SHpZ0VtI6n5X0mLuvkSR3/yC9YQIAMgQ5AQAilErxnidpbcLz8vC1REWSDjOzF8xskZldkq4AAQAZhZwAABFK5YJVq+E1r2E/oySdLKmDpL+a2QJ3f+tjOzKbLmm6JOXn5zc8WgBA1NKWEyTyAhBHbYYURB1Cq5ZK8V4uqV/C876S1tewzmZ33yVpl5m9JGmEpI911O5+n6T7JGn06NHJnf3HNHaus8R8ZwBoAmnLCVLD8gIAILXifaGkQjMbIGmdpAsVzGdM9AdJPzGzHEltJY2V9MN0Bgqg9eLC5YxCTgCACNVbvLt7hZldI+lpSdmS7nf3pWZ2Vbj8XndfbmZ/kvQPSVWSZrn7kqYMHADQ/MgJABCtlL6kyd2fkvRU0mv3Jj3/b0n/nb7QosdoHzIBf4fINK01JwBAJuAbVgEgBbyJAgBkAor3DEaxAAAAgESp3OcdAAAAQAageAcAAABigmkzqBNTdwAAADIHI+8AAABATDDyDgAA4mNwadQRAOpZlBtZ2xTvAOrF9CkAADID02YAAACAmKB4BwAAAGKC4h0AAACICYp3AAAAICYo3gEAAICYoHgHAAAAYoLiHQAAAIgJincAAAAgJijeAQAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmKN4BAACAmKB4BwAAAGKC4h0AAACICYp3AAAAICYo3gEAAICYoHgHAAAAYiIn6gAAAACQOhtcFHUIiBAj7wAAAEBMULwDAAAAMUHxDgAAAMQExTsAAAAQExTvAAAAQExQvAMAAAAxQfEOAAAAxAT3eQcAAECstBlSEHUIkWHkHQAAAIgJincAAAAgJijeAQAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmUirezWySma0ws5VmdmMd6x1jZpVmNjV9IQIAMgk5AQCiU2/xbmbZku6RNFlSiaSLzKyklvVul/R0uoMEAGQGcgIARCuVL2kaI2mlu78rSWb2sKSzJC1LWu+Lkn4r6Zi0RggAyCStNidYjz5RhwAAKU2byZO0NuF5efjaQWaWJ+kcSffWtSMzm25mr5nZa5s2bWporACA6KUtJ4TrkhcAoAFSKd6thtc86fmPJH3V3Svr2pG73+fuo919dG5ubqoxAgAyR9pygkReAICGSmXaTLmkfgnP+0pan7TOaEkPm5kk9ZJ0mplVuPvv0xIlACBTkBMAIEKpFO8LJRWa2QBJ6yRdKOmziSu4+4Dqn81stqQ/0kkDQItETgCACNVbvLt7hZldo+COAdmS7nf3pWZ2Vbi83jmNAICWgZwAANFKZeRd7v6UpKeSXquxg3b3aY0PCwCQqcgJABAdvmEVAAAAiImURt4BAACy+hdHHQLQ6jHyDgAAAMQExTsAAAAQExTvAAAAQExQvAMAAAAxQfEOAAAAxATFOwAAABAT3CoSAAAgVYNLo44ArRwj7wAAAEBMULwDAAAAMUHxDgAAAMQExTsAAAAQExTvAAAAQExQvAMAAAAxwa0iAQAAkDIbXBR1CK0axTsAAADQAG2GFETWNtNmAAAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmKN4BAACAmKB4BwAAAGKC4h0AAACICYp3AAAAICYo3gEAAICYoHgHAAAAYoLiHQAAAIgJincAAAAgJijeAQAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmKN4BAACAmKB4BwAAAGKC4h0AAACICYp3AAAAICYo3gEAAICYoHgHAAAAYoLiHQAAAIiJnFRWMrNJkn4sKVvSLHe/LWn55yR9NXy6U9K/u/vidAYKAMgM5AREKat/cdQhAJGqt3g3s2xJ90g6RVK5pIVm9ri7L0tYbZWk8e6+1cwmS7pP0timCBgAEB1yAoBMYIOLog4hMqlMmxkjaaW7v+vu+yU9LOmsxBXc/S/uvjV8ukBS3/SGCQDIEOQEAIhQKtNm8iStTXherrpHUL4g6f8aExQAIGORE4CoDS6NOgJEKJXi3Wp4zWtc0exfFHTUJ9ayfLqk6ZKUn5+fYogAgAyStpwQrkNeAIAGSGXaTLmkfgnP+0pan7ySmR0laZaks9z9w5p25O73uftodx+dm5t7KPECAKKVtpwgkRcAoKFSKd4XSio0swFm1lbShZIeT1zBzPIlPSbpYnd/K/1hAgAyBDkBACJU77QZd68ws2skPa3gtmD3u/tSM7sqXH6vpG9K6inpp2YmSRXuPrrpwgYARIGcAADRSuk+7+7+lKSnkl67N+HnKyRdkd7QAACZiJwAANFJqXgHAADRsx59og4BQMRSmfMOAAAAIANQvAMAAAAxQfEOAAAAxATFOwAAABATFO8AAABATFC8AwAAADFB8Q4AAADEBMU7AAAAEBMU7wAAAEBMULwDAAAAMUHxDgAAAMQExTsAAAAQEzlRBwAAAADEiQ0uiqxtRt4BAACAmKB4BwAAAGKCaTMAAABI3eDSqCNo1Rh5BwAAAGKCkXcAAIAUZfUvjjoEtHKMvAMAAAAxQfEOAAAAxATFOwAAABATFO8AAABATFC8AwAAADFB8Q4AAADEBLeKBAAAKbEefaIOAWj1GHkHAAAAYoLiHQAAAIgJincAAAAgJpjzDgAAgHgZXBp1BJFh5B0AAACICYp3AAAAICaYNgMAABAjWf2Low4BEWLkHQAAAIgJincAAAAgJijeAQAAgJhgzjsAAIgN69En6hCASG9Vycg7AAAAEBOMvAMAACBl3O0mWoy8AwAAADHByDsAAABipTWP/qdUvJvZJEk/lpQtaZa735a03MLlp0naLWmau7+e5lgBABmAnIDWjAtmEbV6i3czy5Z0j6RTJJVLWmhmj7v7soTVJksqDB9jJf0s/BcA0IKQE4Do8QaidUtl5H2MpJXu/q4kmdnDks6SlNhRnyXpQXd3SQvMrLuZHenuG9IeMQAgSuQEoJXjzUO003ZSKd7zJK1NeF6uT46g1LROniQ6agBoWcgJACLXmt9ApFK8Ww2v+SGsIzObLml6+HSnma1Iof3a9JK0ua4V7PZG7D0NMbSC9jMhhqjbz4QYom4/E2KIuv10xNA/ncE0obTlBKn580ITi7r9TIgh6vYzIYao28+EGKJuPxNiaGz7teaEVIr3ckn9Ep73lbT+ENaRu98n6b4U2qyXmb3m7qPTsa+4xhB1+5kQQ9TtZ0IMUbefCTFE3X6mxNBM0pYTpJaVF6JuPxNiiLr9TIgh6vYzIYao28+EGJqy/VTu875QUqGZDTCztpIulPR40jqPS7rEAsdK2sbcRgBokcgJABChekfe3b3CzK6R9LSC24Ld7+5LzeyqcPm9kp5ScEuwlQpuC3ZZ04UMAIgKOQEAopXSfd7d/SkFnXHia/cm/OySZqQ3tHql5WPWRoo6hqjbl6KPIer2pehjiLp9KfoYom5fyowYmkWG5gQp+t9B1O1L0ccQdftS9DFE3b4UfQxRty9FH0OTtW9BHwsAAAAg06Uy5x0AAABABohl8W5mk8xshZmtNLMbI2j/fjP7wMyWNHfbYfv9zOzPZrbczJaa2bXN3H57M/ubmS0O27+5OdtPiiXbzP5uZn+MoO33zOxNM3vDzF5r7vbDGLqb2VwzKwv/Ho5rxraHhMde/dhuZtc1V/sJcVwf/h0uMbM5Zta+mdu/Nmx7aRTHD3JCGAN5QdHmhLD9SPNClDkhbD/yvBB1TghjaNq84O6xeii4QOodSQMltZW0WFJJM8dwkqSjJS2J6BwcKeno8Ocukt5qznOg4B7OncOf20h6VdKxEZ2LGyQ9JOmPEbT9nqReURx3QgwPSLoi/LmtpO4RxZEtaaOk/s3cbp6kVZI6hM8flTStGdsfJmmJpI4KriF6VlJhlH8Tre1BTjgYA3nBo80JYfuR5oVMyQlh+82eF6LOCWGbTZ4X4jjyfvCrud19v6Tqr+ZuNu7+kqQtzdlmUvsb3P318OcdkpYr+INtrvbd3XeGT9uEj2a/eMLM+ko6XdKs5m47E5hZVwVFwy8lyd33u/tHEYVzsqR33H11BG3nSOpgZjkKOssa7yfeRIolLXD33e5eIelFSec0Y/sgJ1TH0OrzAjkho3KCFF1eiDInSM2QF+JYvNf2tdutkpkVSPqUglGO5mw328zekPSBpHnu3qzth34k6SuSqiJoWwoS0zNmtscalZIAABSdSURBVMiCb4lsbgMlbZL0P+HHxLPMrFMEcUjBvb7nNHej7r5O0vclrZG0QcH9xJ9pxhCWSDrJzHqaWUcFt0fsV882SC9yQpJWnBeizglStHkhk3KCFEFeyICcIDVDXohj8Z7y1263dGbWWdJvJV3n7tubs213r3T3kQq+OXGMmQ1rzvbN7DOSPnD3Rc3ZbpIT3P1oSZMlzTCzk5q5/RwFH9X/zN0/JWmXpCjm+7aVdKak/42g7cMUjLIOkNRHUicz+3xzte/uyyXdLmmepD8pmLJR0VztQxI54WNaa17IkJwgRZsXMiInSNHlhahzgtQ8eSGOxXvKX7vdkplZGwUd9G/c/bGo4gg/kntB0qRmbvoESWea2XsKPiafaGa/bs4A3H19+O8Hkn6n4OP75lQuqTxhdGuugo67uU2W9Lq7vx9B25+WtMrdN7n7AUmPSTq+OQNw91+6+9HufpKCqRNvN2f7ICdUa+V5IfKcIEWeFzIlJ0jR5YXIc4LU9HkhjsV7Kl/N3aKZmSmY07bc3X8QQfu5ZtY9/LmDgv8sZc0Zg7t/zd37unuBgr+B59292d5dm1knM+tS/bOkUxV8VNZs3H2jpLVmNiR86WRJy5ozhtBFimDKTGiNpGPNrGP4/+JkBXN9m42ZHR7+my9piqI7F61Vq88JEnkh6pwgRZ8XMignSNHlhchzgtT0eSGlb1jNJF7LV3M3ZwxmNkfSBEm9zKxc0rfc/ZfNGMIJki6W9GY4v1CSvu7Btx42hyMlPWBm2QreAD7q7pHclitCvSX9LugblCPpIXf/UwRxfFHSb8Ki5V0189fQh/P5TpF0ZXO2W83dXzWzuZJeV/Cx5N/V/N+q91sz6ynpgKQZ7r61mdtv1cgJB5EXopcJeSHSnCBFmxcyJCdITZwX+IZVAAAAICbiOG0GAAAAaJUo3gEAAICYoHgHAAAAYoLiHQAAAIgJincAAAAgJije0SKYWaWZvWFmS8zsf8NbVTVk+1lmVtKA9aeZ2U8aHikAoKmRE9CSUbyjpdjj7iPdfZik/ZKuSnVDM8t29yvcPaovswAApBc5AS0WxTtaopclDZYkM/u8mf0tHIH5efgFIjKznWb2bTN7VdJxZvaCmY0Ol11kZm+GIza3V+/UzC4zs7fM7EUFX4gCAMh85AS0KBTvaFHMLEfSZAXfMlgs6QJJJ7j7SEmVkj4XrtpJ0hJ3H+vu8xO27yPpdkkTJY2UdIyZnW1mR0q6WUEHfYqklD9OBQBEg5yAlign6gCANOmQ8JXgL0v6paTpkkZJWhh+XXUHSR+E61RK+m0N+zlG0gvuvkmSzOw3kk4KlyW+/oikoiY4DgBA45ET0GJRvKOl2BOOpBxkQe/8gLt/rYb197p7ZQ2vWx1teGMCBAA0G3ICWiymzaAle07SVDM7XJLMrIeZ9a9nm1cljTezXuFcyIskvRi+PsHMeppZG0nnNWXgAIC0IyegRWDkHS2Wuy8zs29IesbMsiQdkDRD0uo6ttlgZl+T9GcFIy5PufsfJMnMZkr6q6QNkl6XlN20RwAASBdyAloKc+dTHwAAACAOmDYDAAAAxATFOwAAABATFO8AAABATFC8AwAAADFB8Q4AAADEBMU7AAAAEBMU7wAAAEBMULwDAAAAMUHxDgAAAMQExTsAAAAQExTvAAAAQExQvAMAAAAxQfEOAAAAxATFOwAAABATOVEHgIZZtGjR4Tk5ObMkDRNvvpB5qiQtqaiouGLUqFEfRB0M0BD0rwAUgzxG8R4zOTk5s4444oji3NzcrVlZWR51PECiqqoq27RpU8nGjRtnSToz6niAhqB/BRCHPMbIQvwMy83N3U5iQSbKysry3NzcbQpGLoG4oX8FWrk45DGK9/jJIrEgk4V/n/QtiCP6VwAZn8cyNjBkruzs7FFDhw4tqX6sWLGi7UsvvdRx2rRp/Wrb5o9//GOXf/mXfxlc07ILLrig/6JFi9o3XcSZq2PHjp9KfH7XXXf1vOSSS/Kjiqc+Y8aMGfLSSy91jDoOoCWbO3du14KCgmH5+fnDvv71rx8RdTzpsHLlyjZjx44tGjhwYOngwYNLb7nllsOjjimdKioqVFxcXFJbnourzZs3Z0+aNGnggAEDSgcOHFj67LPPdoo6pnS4+eabDx88eHBpYWFh6RlnnDFg9+7dFnVMDcGc95gruPHJUenc33u3nb6ovnXatWtXVVZWtizxtSFDhuw/6aSTdh9Km4888sjqQ9ku7WZ2S+u51Mxt9Z5LABls97b09gkdu9XbJ1RUVOj666/Pf/rpp98aOHDggREjRhSfe+65H40aNWpvOkOpfOWxtB5b9glT6jy2Nm3a6M477yw/8cQTd2/dujXrU5/6VMlpp522Pd3HVXHTxWk9rpxbfpVSP37rrbf2Hjx48J6dO3dmp7P9auUjh6b1uPq+UZbScU2fPr3fqaeeuv1Pf/rTu3v37rWdO3emfdD3/q65aT22y7dvqvPYVq1a1ea+++7rvWLFiiWdO3f20047beCsWbN6fOlLX/ownXE0JUbekRaJI+tPPvlk5+pR+eLi4pKtW7dmSdKuXbsOvoM/88wzB1RVVUn6+Ghux44dP/XFL34xb8iQISUjRowYunbt2hxJWrp0absRI0YMHTZsWPF1113XJ3nEuiV66KGHuh111FFDi4uLS44//vii6nNxww039JkyZUrBCSecUJiXlzf8gQce6H7VVVf1LSoqKhk3blzhvn37TJLy8vKGX3PNNXkjR44cOmzYsOL58+d3PPHEEwv79es37I477siVpKqqKl155ZV9CwsLS4uKikp+8YtfHFbd/je+8Y3eRUVFJUOGDCm5+uqr8xJjq6ys1JQpUwq+9KUv9WnOcwK0dC+88EKn/v377yspKdnfvn17nzJlypa5c+d2jzquxurfv/+BE088cbckHXbYYVWDBg3as2bNmrZRx5UO77zzTpunn36627/9279tjjqWdNqyZUvWq6++2uW6667bLEnt27f3Xr16VUYdVzpUVlbarl27sg4cOKA9e/Zk9e3b90DUMTUExTsabN++fVnVxfkpp5wyKHn5nXfeecRdd921uqysbNmCBQvKOnfuXCVJy5cv73DPPfesXbly5dI1a9a0mzdvXufkbffs2ZN13HHH7VyxYsWy4447bufdd9+dK0nXXHNNv6uvvvqDJUuWLO/Tp0+s/pPVJfFcDh06tOR73/vewWL4lFNO2fnGG2+ULV++fNnUqVO3fPvb3z748fnq1avbPf/88yvnzp278qqrrhowceLE7W+99day9u3bVz366KPdqtfr16/f/jfeeKNs7NixOy+//PKCJ5544p1XX3217LbbbusjSQ8++GD3N998s8Py5cuXPvfcc29985vf7Lt69eo2jz76aNcnn3zysEWLFpWtWLFi2be+9a2N1fs8cOCAnX322QMKCwv33nXXXeub61wBrcHatWvb5uXl7a9+3rdv3/3r1q1rEUVutRUrVrRdtmxZx/Hjx++MOpZ0mDFjRr877rijPCurZZVUZWVl7Xr06FFx3nnnFRQXF5dccMEF/bdv3x77gxwwYMCBGTNmbBwwYMBRhx9++IguXbpUTpkyZXvUcTVE7H8JaH7V02bKysqWzZs3753k5ccee+zOL3/5y/1uvfXWwzdv3pzdpk0bSdLw4cN3DRo06EB2drZKS0t3v/POO59ISG3atPELL7xwmySNGjVq1+rVq9tK0t///vfOl19++RZJuuKKK2Lz0VZ9Es9lWVnZsq997WsHi+FVq1a1HTduXGFRUVHJXXfddURZWVmH6mWf/vSnt7Vr187HjBmzp7Ky0qZOnbpdkkpLS/esWrXq4Hk9//zzP5Kk4cOH7z766KN3HXbYYVV9+vSpaNeuXdXmzZuzX3755S7nn3/+lpycHPXr169i7NixO+fPn99x3rx5XT//+c9v7tKlS5Uk9e7d++Boy9VXX92/pKRkz+23336woAeQHu6fvF7WzFrMRbTbtm3LmjJlyqDbbrttbY8ePaqijqex5syZ061Xr14V48aNO6Rpo5msoqLCli9f3nHGjBmbli9fvqxjx45VN910U+yvwdi0aVP2k08+2X3lypVvbty48R+7d+/O+ulPf9oj6rgaguIdaffd735346xZs1bv2bMn6/jjjy/++9//3l6S2rVrdzABZWdnq6Ki4hMXiOTk5Hj16EVOTk6N67QW11xzTf7VV1/9wVtvvbXsJz/5yep9+/Yd/P9afS6zs7M/ds6ysrI+ds7at2/v1a+3bdv24PnPysrSgQMHrKZCQQoKCLOaT/3o0aN3vvzyy13jdoEPEAf5+fkfG2kvLy9v21I+bdy3b5+dfvrpg84777wtl1566UdRx5MO8+fP7zxv3rzueXl5w6dNmzZwwYIFXc4666wBUceVDgUFBft79+69f+LEibsk6YILLti6ePHi2N+w4Iknnuian5+/LxzI8rPPPvujv/zlL5+YCZDJKN6RdkuXLm03ZsyYPd/5znc2Dh8+fNeSJUsafSeZkSNH7pw9e/ZhknT//ffH6h3yodqxY0d2fn7+AUmaPXt2z6ZoY/z48Tvmzp3bo6KiQuvXr8/529/+1nncuHG7Jk2atP1Xv/pVrx07dmRJ0vvvv3/wIqwrr7xy86mnnrrtM5/5zKADB1pETQFkjPHjx+9677332peVlbXdu3evPfbYYz3OPffc2Be6VVVVuvDCC/sXFRXtnTlz5vtRx5Mu99xzz7r333//H+vWrXtz9uzZ7x577LE7/vCHP6yKOq50yM/PrzjiiCP2L168uJ0kPfPMM12HDBmS1guMo1BQULD/9ddf77xjx46sqqoqPf/8812Ki4tjdVwU70i7O+644/DCwsLSIUOGlHTo0KFq6tSp2xq7z7vvvnvt3Xff3Xv48OHFGzZsaNO5c+cWcdFMXf7rv/5r/UUXXTRo1KhRQ3r27FnRFG1cfPHFH5WWlu4pLi4unTBhQtHNN99cnp+fXzF16tTtkydP/mjkyJHFQ4cOLbnllls+9lHpzJkz3x8xYsTuKVOmDKisbPG/CqDZhHdlWTNp0qSiwsLC0rPPPnvL6NGjY1VY1GTevHmdf//73/ecP39+l+prfB555JFu9W+JKN19991rPve5zw0sKioq+cc//tHh1ltv3RB1TI01ceLEXWecccbWo446qnjIkCGlVVVVdsMNN2yKOq6GqPVjc2SmxYsXvzdixIgWdUV7Knbs2JHVqVOnqqysLN13332HPfLIIz2ee+65T8y3R2ZYvHhxrxEjRhREHQfQEK21fwXwSZmcx7jPO2LhlVde6Xjttdfmu7u6du1aOXv27PeijgkAAKC5UbwjFiZNmrRzxYoVy+pfEwAAoOVizjsAAAAQExTv8VNVVVXFLfqQscK/z9jfvxkAgExE8R4/SzZt2tSNAh6ZqKqqyjZt2tRN0pKoYwEAoCViznvMVFRUXLFx48ZZGzduHCbefCHzVElaUlFRcUXUgQAA0BJRvMfMqFGjPpB0ZtRxAADS77zzzit47rnnuvXs2bPi7bffXhp1POmye/duGzt27ND9+/dbZWWlnXHGGVt/+MMfro86rnTJy8sb3qlTp8qsrCzl5OT4kiVLlkcdU2MtXry43QUXXDCo+nl5eXm7r3zlK+u++c1vfhBlXOlwyy23HP7ggw/mursuueSSTXE7Jop3AABq4OXLR6Vzf9a3eFF961x++eWbr7322g8uu+yyAelsO1nlA7ek9diyL72pzmNr3769z58/f0W3bt2q9u3bZ8ccc8yQ5557btvJJ5+8K51x7Dl/XFqPq8OjL9f7O6v24osvvnXkkUc2yRfqvdy7b1qPa9z75fUe14gRI/aVlZUtk6SKigodccQRIy688MK0f9vvVdY1rcd2r2+v89gWLlzY/sEHH8x9/fXXl7dv375q/PjxReecc8624cOH70tnHE2JaRcAAGSIyZMn78zNzW2SAjBKWVlZ6tatW5Uk7d+/3yoqKsyMS7fi4vHHH++an5+/r6ioaH/UsTTWm2++2eHoo4/e2aVLl6o2bdrohBNO2PHII490jzquhqB4BwAATa6iokJDhw4t6d2794jx48dvnzhxYlpH3aN28sknF5aWlhZ///vf7xV1LOk2Z86cHlOnTv0w6jjSYeTIkXteffXVLhs3bszesWNH1rx587qtXbu2bdRxNQTTZgAAQJPLyclRWVnZss2bN2effvrpgxYuXNj+mGOO2Rt1XOnwyiuvlBUUFBxYt25dzsSJE4tKS0v3Tp48eWfUcaXD3r177dlnn+32gx/8oDzqWNLh6KOP3nvttddunDhxYlHHjh2rSkpKdufkxKscZuQdAAA0m169elWeeOKJO5544oluUceSLgUFBQckKS8vr+L000//6K9//WunqGNKl7lz53YrKSnZ3a9fvxYznev666/fvGzZsuWvvfbaih49elQWFhbG6k0kxTsAAGhS69evz9m8eXO2JO3cudNeeOGFrsXFxbEqmGqzffv2rK1bt2ZV//znP/+561FHHbUn6rjS5eGHH+5x/vnnb4k6jnRat25djiS9/fbbbZ988snuX/jCF2J1fPH6nAAAgBbsjDPOGLBgwYIuW7duzendu/dRN9544/rrr79+c9RxNdbatWvbTJs2bUBlZaXc3c4666wtF1100bao40qH8vLynHPOOWewJFVWVtq555774dSpU7dHHVc67NixI2v+/PldH3jggdVRx5JOZ5555qCPPvooJycnx3/0ox+tyc3NrYw6poYwd486BgAAIrd48eL3RowYEftCGUDjLV68uNeIESMKoo6jJkybAQAAAGKC4h0AAACICYp3AAAAICYo3gEACFRVVVXxtZ9AKxf2A1VRx1EbincAAAJLNm3a1I0CHmi9qqqqbNOmTd0kLYk6ltpwq0gAACRVVFRcsXHjxlkbN24cJga3gNaqStKSioqKK6IOpDbcKhIAAACICUYWAAAAgJigeAcAAABiguIdAAAAiAmKdwAAACAmKN4BAACAmPh/UfeeCV6vXxAAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, axs = plt.subplots(1, 2, figsize=(12.8, 4.8))\n", "\n", "(df.groupby(\"Period\").Choice.value_counts(normalize=True).unstack()\n", " .plot.bar(ax=axs[0], stacked=True, rot=0, title=\"Choice Probabilities\"))\n", "(df.groupby(\"Period\").Experience_Fishing.value_counts(normalize=True).unstack().plot\n", " .bar(ax=axs[1], stacked=True, rot=0, title=\"Share of Experience Level per Period\", cmap=\"Reds\"))\n", "\n", "axs[0].legend([\"Fishing\", \"Hammock\"], loc=\"upper center\", bbox_to_anchor=(0.5, -0.15), ncol=2)\n", "axs[1].legend(loc=\"upper center\", bbox_to_anchor=(0.5, -0.15), ncol=5)\n", "\n", "plt.show()\n", "plt.close()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The figure on the left-hand-side shows the choice probabilities for Robinson in each period. On the right-hand-side, one can see the average experience in fishing. By default Robinson starts with zero experience in fishing.\n", "\n", "What if Robinson has an equal probability start with zero, one or two periods of experience in fishing? The first way to add more complex initial conditions to the model is via **probability mass functions**. This type of distributions is handy if the probabilities do not depend on any information. To feed the information to **respy** use the keyword ``initial_exp_fishing_*`` in the category-level of the index. Replace ``*`` with the experience level. In the name-level, use ``probability`` to signal that the float in ``value`` is a probability. The new parameter specification is below.\n", "\n", "Note that one probability is set to 0.34 such that all probabilities sum to one. If that is not the case, respy will emit a warning and normalize probabilities." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "params = \"\"\"\n", "category,name,value\n", "delta,delta,0.95\n", "wage_fishing,exp_fishing,0.01\n", "nonpec_hammock,constant,1\n", "shocks_sdcorr,sd_fishing,1\n", "shocks_sdcorr,sd_hammock,1\n", "shocks_sdcorr,corr_hammock_fishing,0\n", "initial_exp_fishing_0,probability,0.33\n", "initial_exp_fishing_1,probability,0.33\n", "initial_exp_fishing_2,probability,0.34\n", "\"\"\"" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
value
categoryname
deltadelta0.95
wage_fishingexp_fishing0.01
nonpec_hammockconstant1.00
shocks_sdcorrsd_fishing1.00
sd_hammock1.00
corr_hammock_fishing0.00
initial_exp_fishing_0probability0.33
initial_exp_fishing_1probability0.33
initial_exp_fishing_2probability0.34
\n", "
" ], "text/plain": [ " value\n", "category name \n", "delta delta 0.95\n", "wage_fishing exp_fishing 0.01\n", "nonpec_hammock constant 1.00\n", "shocks_sdcorr sd_fishing 1.00\n", " sd_hammock 1.00\n", " corr_hammock_fishing 0.00\n", "initial_exp_fishing_0 probability 0.33\n", "initial_exp_fishing_1 probability 0.33\n", "initial_exp_fishing_2 probability 0.34" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "params = pd.read_csv(io.StringIO(params), index_col=[\"category\", \"name\"])\n", "params" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "C:\\Users\\tobia\\git\\respy\\respy\\state_space.py:84: UserWarning: Some choices in the model are not admissible all the time. Thus, respy applies a penalty to the utility for these choices which is -400000 by default. For the full solution, the penalty only needs to be larger than all other value functions to be effective. Choose a milder penalty for the interpolation which does not dominate the linear interpolation model.\n", " \"Some choices in the model are not admissible all the time. Thus, respy\"\n" ] } ], "source": [ "simulate = rp.get_simulate_func(params, options)\n", "df = simulate(params)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAu8AAAFwCAYAAAAIW5VaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deXxU9b3/8fcnCfsqqxCyAAZIAEMBxQWFi9WCFhfErYtFa9Gf2LrcLtbbulRr1Vtrr0tVSr1oW1Evta3WFVFUrLhgRVmCogiERRNA9i3J9/fHOaHjOEkGMsl3TvJ6Ph7zIHOW+X7OSfh+PvOd7zljzjkBAAAASH8ZvgMAAAAAkByKdwAAACAiKN4BAACAiKB4BwAAACKC4h0AAACICIp3AAAAICIo3nFQzOx6M/vTQe57n5n9PNUxpYqZ5ZuZM7Osg9zfmdlhNaz7ppk9n2jbus6LmV1jZjMOJiYA6cXMppjZfN9xJGJmPc3sFTPbZma3+46nNma23cz6+Y7Dl3T+OzpQZvaMmX3nIPf9xMy+muqY0hXFO2pkZt8ws7fDznF9+B9rdH1f1zl3iXPuxlTEWC18M7EvjPVzM/unmR2dyjZSwTn3Z+fcSTWs239ezGysmZXGrb/ZOXdRY8QJoP7MbHTYF20xs01m9pqZHeE7riRMlVQuqaNz7j/jV5rZTDPbG/a31Y9FjR+m5Jxr75z72EfbiTSnIjLMU1Xh73+bmS03swsO9vWccxOccw+mMsamiuIdCZnZVZJ+K+lmST0l5Ur6naTTfMZVh0edc+0ldZc0X9LjZmbxGx3siDoAJMvMOkr6h6S7JHWRlC3pBkl7GqCtVPdpeZKWutq/xfG2sHCufhSnOIZa0Y83rlrO97ow73aU9BNJvzezogN8bTMz6tEDwMnCl5hZJ0m/kDTNOfe4c26Hc26fc+5J59yPYjZtaWYPhe+4l5jZyJjXKDSzeeEo+BIzOzVm3Uwzuynm+Wlm9q6ZbTWzj8xsfHUcZvaHcNR/rZndZGaZdcXvnNsn6UFJh0rqGn6s+JqZ3WFmmyRdb2YZZvYzM1tlZp+Fx9Ep7qUuNLN1Yfv7R5/M7Egzez08tvVmdreZtYzb92Qz+9jMys3sv6s7pto+4qw+L2bWTtIzknrHjGr1tripSmZ2VDiq97mZLTKzsTHrpoTtbzOzlWb2zbrOG4CUGiBJzrlZzrlK59wu59zzzrn3Yjcys1+b2ebw/+mEmOUXmNmy8P/wx2Z2ccy6sWZWamY/MbMNkv43XP71sC+t/vTx8JqCM7NjzOyt8FOBt8zsmHD5TEnfkfTjsO85oFFkMzsnjLdj+HyCmW0ws+7hc2dmP0jUP4brLwyPe7OZPWdmeTHrnJlNM7MPJX0Ys6x66mGr8HyuNrNPLZiK2CbunP1n2Oevt5hRYjNrY2a3hzlhi5nNj9m3xr72AM9Nwt+PmV1tZrPjtv0fM7sz/PmgcqH9ewro1BpyWUbY9kdmttHMHjOzLnH7ftfMVkt6sba2XOBvkjZLKgpfo7YcNc/Mfmlmr0naKalfuOyimNhqzNFm9u1w3UYz+6+6zkWT45zjweMLD0njJVVIyqplm+sl7ZZ0sqRMSb+StCBc10LSCknXSGopaZykbZIGhutnSrop/PlISVsknajgzWS2pEHhur9Jul9SO0k9JL0p6eJa4vlT+HMrSf8taU34fEp4PN+XlCWpjaQLwxj7SWov6XFJfwy3z5fkJM0K2x4qqUzSV8P1IyQdFb5WvqRlkq6IicVJeknBaFuupA8kXRQTy/y4bQ9LcF7GSiqt5RizJW0Mz39GeP42KvjUoZ2krTHnu5ekwb7/rnjwaE4PBSORGxUMJEyQdEjc+imS9kn6XtiH/j9J6yRZuP4USf0lmaQxCgqc4eG6sWGfdmvY37WRNFzSZ5JGha/3HUmfSGqVILYuCoqsb4f92Hnh867h+v19UQ3HVtf6P4fbdA2P6esx62rrH09X0C8XhnH9TNI/4/adE+7bJmZZdR/6W0lPhOs7SHpS0q/iztkvFOSok8Nzeki4/h5J88K+NVPSMeG5rbGvreHYP1GYK+KW1/j7UfBJx04F05QUrl8v6ajweY25UHE5Ja7NfNWey66QtEBSnzCO+yXNitv3oXDfNglef6zCPBWemzMU/E0PrOu8hed6taTB4e+6Rbis+m+hthxdJGm7pOPDuH8T/m6/dN6b6sN7ADzS7yHpm5I21LHN9ZJeiHleJGlX+PNxkjZIyohZP0vS9eHPM/XvIvV+SXckeP2eCj5ebhOz7DxJL9USz15Jn4cd5IuSRoTrpkhaHbf9XEmXxjwfGHY61QW5U/gmIlx/m6Q/1ND2FZL+GvPcSRof8/xSSXNjYklF8f6T6o4sZv1zChJCu/A8nKkEHS4PHjwa56GgCJ0pqTQsLp6Q1DNcN0XSipht24b9waE1vNbfJF0e/jw27O9ax6y/V9KNcfsslzQmwWt9W9KbcctelzQl/Hl/X1RDLDMVDN58HvN4MGZ9ZwWF2fuS7o/bt7b+8RlJ341Zl6GgqM2L2Xdcgtc7TMGbnB2S+sesO1rSyphztksxg1IKcsVRYTu7JBUnONYa+9oazs0nSly81/r7UTDV8/zw5xMlfRT+XGsuVHLFe8JcpmDg6YSYdb305TzYr5a/g7GSqsLf/yZJ70o6N5nzpqBQ/0Xc+nn6d/FeW46+VtIjMevaKfj/0GyKd+aMIZGNkrqZWZZzrqKW7TbE/LxTUmsL5sX1VjDqXRWzfpWCd+LxciQ9nWB5noJ34uvt39PWMyStqSWex5xz36phXfx+vcOYYuPLUtBRJtpnlYJRC5nZAAXv9EcqSLhZkhbW0t6qsL1UypN0lplNjFnWQkGHvsPMzpH0Q0l/CD+W/E/nXEmKYwBQC+fcMgXFlcxskKQ/KRgdPi/cZEPMtjvDvq59uP0ESdcpmH6ToaCveT/m5cucc7tjnudJ+o6ZfT9mWUsl7nvi+z+p5j66Jr92zv0s0Qrn3Odm9n+SrlIwiBCvpv4xT9L/2BfvcGNhXKsS7Buru4JztDAmZ5iCUexqG+Ny2k4F57ubpNaSPkrwujX2tTXEUZO6fj8PK/i7eEjSN8Ln1fsdaC6MlzCXha/9VzOLzdWVqjkPJrLOOdcnwfJkzlttr11bju4du2+Y8zbWEWeTwpx3JPK6glGV0w9y/3WScuyLF6DkSlqbYNs1Cj4aTrR8j6RuzrnO4aOjc27wQcbkEsSYFxdfhaRPY5blxK1fF/58r6QSSQXOuY4KpgfFXxhb074HG2+8NQpGNTrHPNo5526RJOfcc865ExWMpJRI+v0Btg8ghcI3zzMlDalrWzNrJekvkn6tYKS+s4JBjth+Jr6PWCPpl3F9Qlvn3KwETcT3f1LNffQBM7NhCqY9zJJ0Z4JNauof1yiYDhJ7DG2cc/+M2b6mvrFcwej54Jh9O7ngYsq6lCvIeTXlohr72gNQ1+/n/ySNNbM+CqafPByzX31zYW3ne0JcTK2dc7F/B3Xlopokc95qe+3acvT62GMys7YKpmg1GxTv+BLn3BYFH0vdY2anm1lbM2sRXnh0WxIv8YaCjy9/HO43VtJESY8k2PYPki4wsxPCC1SyzWyQc269pOcl3W5mHcN1/c1sTGqOUrMkXWlmfc2svYK76jwaNyrz8/DYB0u6QNKj4fIOCuaUbw9H0/5fgtf/kZkdYmY5ki6P2TdZnyq42Db+Itpqf5I00cy+ZmaZZtY6vCCrjwX3aD7Vggtf9yiYG1h5gO0DqAczGxReHNknfJ6jYGR1QRK7t1Qwl7dMUkU4Cp/wFrMxfi/pEjMbZYF2ZnaKmXVIsO3TkgZYcDvgrPCTuiIFd8epFzNrraB/ukZBv5ltZpfGbVZT/3ifpJ+GfW71hZpnJdNu+Env7yXdYWY9wv2zzexrSe77gKTfWHBzgEwzOzp8E1VjX1vLS7YIt6t+ZKmO349zrkzBtJH/VTDVZ1m4PBW5sKZcdp+kX1p4UbCZdTezVN1R7mDOW6zacvRsSV+34FasLRVcx9Cs6tlmdbBInnPuNwo+8vyZggSyRtJlCuZd1rXvXkmnKrhIq1zBLSbPTzRtwzn3poLO5A4FF66+rH+/2z5fQRJbquBiqtkKRpJT4QFJf5T0iqSVCkZdvh+3zcsKLpiZq+Aj4uovV/qhgo81tynokBMV5n9XMJXmXUlPKXiTkrTwXM2S9LEFV+r3jlu/RsFtO6/Rv38/P1LwfzpD0n8qGLnYpOBit/jkCaBhbVNwceIbZrZDQdG+WMH/zVo557ZJ+oGkxxT0fd9QMF++tn3eVnDx693hPisUTtlJsO1GSV8PY9ko6ccKLiotT+K4qlXfjab6Ub3vrxRcr3Ovc26PpG9JusnMCmL2Tdg/Ouf+quAi3EfMbKuC8zVByfuJguNeEO7/goK50sn4oYJpSW8p6DdvVXDdVm19bU2eVvApQPXj+iR/Pw9L+qr+Peperb65sKZc9j8K/q6eN7NtCv5GRx3A69boIM9brBpztHNuiaRpCs7TegXnpDTxyzRN1Ve1AwAANCgzcwqmHK7wHUtTZ2b5CgrfFnVcv4aIYeQdAAAAiAiKdwAAACAimDYDAAAARAQj7wAAAEBEULwDAAAAEeHtG1a7devm8vPzfTUPAGll4cKF5c657r7j8Im8AACB2nKCt+I9Pz9fb7/9tq/mASCtmFn819U3O+QFAAjUlhOYNgMAAABEBMU7AAAAEBEU7wAAAEBEULwDAAAAEUHxDgAAAEQExTsAAAAQERTvAAAAQETUWbyb2QNm9pmZLa5hvZnZnWa2wszeM7PhqQ8TAJAuyAsA4E8yI+8zJY2vZf0ESQXhY6qke+sfFgAgjc0UeQEAvKizeHfOvSJpUy2bnCbpIRdYIKmzmfVKVYAAgPRCXgAAf7JS8BrZktbEPC8Nl62P39DMpioYhVFubm7tr3p9p/pHdv2Weu5fzxii3n46xOC7/XSIwXf76RCD7/bTJYboaJC8ULVobr2Cyig+oV77Vz41vV77S1LmKVPrF8PvflK/9i+9tV77V/z0vHrtL0lZv5pVr/13f2NsvfZv/fC8eu2/8egh9dq/6+sJZ5sdkPfy8uu1/+GrPqnX/s90612v/SeUr6vX/pJ0T/tu9dp/2vbyyLafiuLdEixziTZ0zk2XNF2SRo4cmXAbAEDkkRcANKh+rVv6DsGbVBTvpZJyYp73kVT/t1QAgKgiLwBNWHa75ls4p4NUFO9PSLrMzB6RNErSFufclz4aBQA0G+QFoAFl927vOwR4VGfxbmazJI2V1M3MSiVdJ6mFJDnn7pP0tKSTJa2QtFPSBQ0VLADAP/ICAPhTZ/HunKv16hTnnJM0LWURAQDSGnkBAPxJxbQZAACAZqFd3/rdZQRNg88LZpP5kiYAAAAAaYDiHQAAAIgIincAAAAgIijeAQAAgIjgglUAABAZWXm9fIeANNCcvyiKkXcAAAAgIijeAQAAgIigeAcAAAAiguIdAAAAiAiKdwAAACAiuNsMAAAAkpbdu73vEJo1Rt4BAACAiGDkHQAAJMVycn2HADR7jLwDAAAAEUHxDgAAAEQExTsAAAAQERTvAAAAQERQvAMAAAARwd1mAAAAIqRd326+Q4BHjLwDAAAAEUHxDgAAAEQE02YAAACAA5DdrqW3thl5BwAAACKCkXcAAKIir7/vCAB4xsg7AAAAEBEU7wAAAEBEULwDAAAAEcGcdwAAgCRl5fXyHQKaOUbeAQAAgIigeAcAAAAigmkzAAAAiJTs3u19h+ANI+8AAABARFC8AwAAABFB8Q4AAABEBMU7AAAAEBEU7wAAAEBEULwDAAAAEUHxDgAAAEQExTsAAAAQEUkV72Y23syWm9kKM7s6wfpOZvakmS0ysyVmdkHqQwUApANyAgD4U2fxbmaZku6RNEFSkaTzzKwobrNpkpY654oljZV0u5m1THGsAADPyAkA4FdWEtscKWmFc+5jSTKzRySdJmlpzDZOUgczM0ntJW2SVJHiWAEA/pET4JXl5PoOAfAqmWkz2ZLWxDwvDZfFultSoaR1kt6XdLlzriolEQIA0gk5AQA8SqZ4twTLXNzzr0l6V1JvScMk3W1mHb/0QmZTzextM3u7rKzsgIMFAHiXspwgkRcA4EAlM22mVFJOzPM+CkZTYl0g6RbnnJO0wsxWShok6c3YjZxz0yVNl6SRI0fGd/YAgPSXspwgkReAKGrXt5vvEJq1ZIr3tyQVmFlfSWslnSvpG3HbrJZ0gqRXzaynpIGSPk5loACAtNBsc4L16e87BACou3h3zlWY2WWSnpOUKekB59wSM7skXH+fpBslzTSz9xV8pPoT51x5A8YNAPCAnAAAfiUz8i7n3NOSno5bdl/Mz+sknZTa0AAA6YicAAD+8A2rAAAAQERQvAMAAAARQfEOAAAARATFOwAAABARFO8AAABARFC8AwAAABFB8Q4AAABEBMU7AAAAEBEU7wAAAEBEULwDAAAAEZHlOwAAAAAgSrJ7t/fWNiPvAAAAQERQvAMAAAARQfEOAAAARATFOwAAABARXLAKAAAQIVl5vXyHAI8o3gEAQHLy+vuOAGj2mDYDAAAARATFOwAAABARFO8AAABARFC8AwAAABFB8Q4AAABEBMU7AAAAEBEU7wAAAEBEULwDAAAAEUHxDgAAAEQExTsAAAAQERTvAAAAQERQvAMAAAARkeU7AAAAAOBAtOvbzXcI3jDyDgAAAEQExTsAAAAQEUybAQBEgvXI9R0CAHjHyDsAAAAQERTvAAAAQERQvAMAAAARQfEOAAAARAQXrAIAACTJcrhwGn4x8g4AAABEBMU7AAAAEBFMmwHqkL/74Xrt/0lqwgAAAEhu5N3MxpvZcjNbYWZX17DNWDN718yWmNnLqQ0TAJAuyAkA4E+dI+9mlinpHkknSiqV9JaZPeGcWxqzTWdJv5M03jm32sx6NFTAaF4Y9QbSCzkBAPxKZuT9SEkrnHMfO+f2SnpE0mlx23xD0uPOudWS5Jz7LLVhAgDSBDkBADxKpnjPlrQm5nlpuCzWAEmHmNk8M1toZuenKkAAQFohJwCAR8lcsGoJlrkErzNC0gmS2kh63cwWOOc++MILmU2VNFWScnO5TyoARFDKcoJEXgCAA5VM8V4qKSfmeR9J6xJsU+6c2yFph5m9IqlY0hc6aufcdEnTJWnkyJHxnf0X1Heus8R8ZwBoACnLCdKB5QUAQHLF+1uSCsysr6S1ks5VMJ8x1t8l3W1mWZJaShol6Y5UBgqg+eLC5bRCTvDI+vT3HQIAz+os3p1zFWZ2maTnJGVKesA5t8TMLgnX3+ecW2Zmz0p6T1KVpBnOucUNGXhjoGAAgC9qzjkBANJBUl/S5Jx7WtLTccvui3v+35L+O3WhAZB4E4n0Q04AAH/4htU0RtEGpA/+PwIA0kFS37AKAAAAwD+KdwAAACAimDaDWjFVAAAAIH1QvAMAACBpWXm9fIfQrDFtBgAAAIgIRt4B1InpUwAApAdG3gEAAICIoHgHAAAAIoLiHQAAAIgIincAAAAgIijeAQAAgIjgbjMAACA68vr7jgBQu77dvLXNyDsAAAAQERTvAAAAQERQvAMAAAARQfEOAAAARATFOwAAABARFO8AAABARFC8AwAAABFB8Q4AAABEBMU7AAAAEBEU7wAAAEBEULwDAAAAEUHxDgAAAEQExTsAAAAQERTvAAAAQERQvAMAAAARQfEOAAAARATFOwAAABARWb4DAAAAQPIsJ9d3CPCIkXcAAAAgIijeAQAAgIigeAcAAAAigjnvAAAkwXowzxiAfxTvAAAAiJSsvF6+Q/CGaTMAAABARFC8AwAAABFB8Q4AAABEBMU7AAAAEBEU7wAAAEBEJFW8m9l4M1tuZivM7OpatjvCzCrNbHLqQgQApBNyAgD4U2fxbmaZku6RNEFSkaTzzKyohu1ulfRcqoMEAKQHcgIA+JXMyPuRklY45z52zu2V9Iik0xJs931Jf5H0WQrjAwCkF3ICAHiUTPGeLWlNzPPScNl+ZpYt6QxJ99X2QmY21czeNrO3y8rKDjRWAIB/KcsJ4bbkBQA4AMkU75ZgmYt7/ltJP3HOVdb2Qs656c65kc65kd27d082RgBA+khZTpDICwBwoLKS2KZUUk7M8z6S1sVtM1LSI2YmSd0knWxmFc65v6UkSgBAuiAnNGPWp7/vEIBmL5ni/S1JBWbWV9JaSedK+kbsBs65vtU/m9lMSf+gkwaAJomcAAAe1Vm8O+cqzOwyBXcMyJT0gHNuiZldEq6vc04jAKBpICcAgF/JjLzLOfe0pKfjliXsoJ1zU+ofFgAgXZETAMAfvmEVAAAAiAiKdwAAACAiKN4BAACAiKB4BwAAACIiqQtWAQAAICmPe93DL0beAQAAgIigeAcAAAAiguIdAAAAiAiKdwAAACAiKN4BAACAiKB4BwAAACKC4h0AAACICO7zDgAAAByArLxe/tr21jIAAAAix3JyfYfQrDFtBgAAAIgIincAAAAgIijeAQAAgIigeAcAAAAiggtWAQCRYJ16+A4BALxj5B0AAACICIp3AAAAICIo3gEAAICIoHgHAAAAIoLiHQAAAIgIincAAAAgIijeAQAAgIigeAcAAAAiguIdAAAAiAiKdwAAACAiKN4BAACAiKB4BwAAACKC4h0AAACIiCzfAQAAgORYj1zfIQDwjOIdAABEhvXp7zsEwCumzQAAAAARQfEOAAAARATFOwAAABARzHkHAABApFhO8714m+IdAAAgSvK4aLc5Y9oMAAAAEBFJFe9mNt7MlpvZCjO7OsH6b5rZe+Hjn2ZWnPpQAQDpgJwAAP7UWbybWaakeyRNkFQk6TwzK4rbbKWkMc65wyXdKGl6qgMFAPhHTgAAv5IZeT9S0grn3MfOub2SHpF0WuwGzrl/Ouc2h08XSOqT2jABAGmCnAAAHiVTvGdLWhPzvDRcVpPvSnqmPkEBANIWOQEAPErmbjOWYJlLuKHZfyjoqEfXsH6qpKmSlJvbfG/xAwARlrKcEG5DXgCAA5DMyHuppJyY530krYvfyMwOlzRD0mnOuY2JXsg5N905N9I5N7J79+4HEy8AwK+U5QSJvAAAByqZ4v0tSQVm1tfMWko6V9ITsRuYWa6kxyV92zn3QerDBACkCXICAHhU57QZ51yFmV0m6TlJmZIecM4tMbNLwvX3SbpWUldJvzMzSapwzo1suLABAD6QEwDAr6S+YdU597Skp+OW3Rfz80WSLkptaACAdEROAAB/+IZVAAAAICIo3gEAAICIoHgHAAAAIoLiHQAAAIgIincAAAAgIijeAQAAgIigeAcAAAAiguIdAAAAiIikvqQJAAAAQMBycr21TfEOAACSYj38FSwAAhTvAAAASbI+/X2HgGaO4h0AAADJy+MNjE9csAoAAABEBCPvAAAkwTr18B0CADDyDgAAAEQFxTsAAAAQERTvAAAAQERQvAMAAAARQfEOAAAARATFOwAAABARFO8AAABARFC8AwAAABFB8Q4AAABEBN+wCgAAECHWp7/vEPzLa77ngJF3AAAAICIYeQcAAJFhPXJ9hwB4xcg7AAAAEBEU7wAAAEBEULwDAAAAEUHxDgAAAEQEF6wCABAR1qmH7xAAeEbxDgAAgKRxn3l5vc8802YAAACAiKB4BwAAACKCaTMAAABJ4kui4Bsj7wAAAEBEMPIOAACASGnOF80y8g4AAABEBMU7AAAAEBFMmwEAAEnhS6IA/yjeAQAAIoQ73vjnc859UsW7mY2X9D+SMiXNcM7dErfewvUnS9opaYpz7p0UxwoASAPkBKB5482DX3UW72aWKekeSSdKKpX0lpk94ZxbGrPZBEkF4WOUpHvDfwEATQg5Ab4xdQdS834DkczI+5GSVjjnPpYkM3tE0mmSYjvq0yQ95JxzkhaYWWcz6+WcW5/yiAEAPpET0Kylw5uHdIgB/iRTvGdLWhPzvFRfHkFJtE22JDpqAGhayAlAM5cObx58x+Bz5D+Z4t0SLHMHsY3MbKqkqeHT7Wa2PIn2a9JNUnltG9it9Xj1FMTQDNpPhxh8t58OMfhuPx1i8N1+KmLIS2UwDShlOUFq/LzQwHy3nw4x+G4/HWLw3X46xOC7/XSIob7t15gTkineSyXlxDzvI2ndQWwj59x0SdOTaLNOZva2c25kKl4rqjH4bj8dYvDdfjrE4Lv9dIjBd/vpEkMjSVlOkJpWXvDdfjrE4Lv9dIjBd/vpEIPv9tMhhoZsP5kvaXpLUoGZ9TWzlpLOlfRE3DZPSDrfAkdJ2sLcRgBoksgJAOBRnSPvzrkKM7tM0nMKbgv2gHNuiZldEq6/T9LTCm4JtkLBbcEuaLiQAQC+kBMAwK+k7vPunHtaQWccu+y+mJ+dpGmpDa1OKfmYtZ58x+C7fcl/DL7bl/zH4Lt9yX8MvtuX0iOGRpGmOUHy/zvw3b7kPwbf7Uv+Y/DdvuQ/Bt/tS/5jaLD2LehjAQAAAKS7ZOa8AwAAAEgDkSzezWy8mS03sxVmdrWH9h8ws8/MbHFjtx22n2NmL5nZMjNbYmaXN3L7rc3sTTNbFLZ/Q2O2HxdLppn9y8z+4aHtT8zsfTN718zebuz2wxg6m9lsMysJ/x6ObsS2B4bHXv3YamZXNFb7MXFcGf4dLjazWWbWupHbvzxse4mP4wc5IYyBvCC/OSFs32te8JkTwva95wXfOSGMoWHzgnMuUg8FF0h9JKmfpJaSFkkqauQYjpc0XNJiT+egl6Th4c8dJH3QmOdAwT2c24c/t5D0hqSjPJ2LqyQ9LOkfHtr+RFI3H8cdE8ODki4Kf24pqbOnODIlbZCU18jtZktaKalN+PwxSVMasf0hkhZLaqvgGqIXJBX4/Jtobg9ywv4YyAvOb04I2/eaF9IlJ4TtN3pe8J0TwjYbPC9EceR9/1dzO+f2Sqr+au5G45x7RdKmxmwzrv31zrl3wp+3SVqm4A+2sdp3zrnt4dMW4aPRL54wsz6STpE0o7HbTgdm1lFB0fAHSXLO7XXOfe4pnBMkfeScW+Wh7SxJbcwsS0FnmfB+4g2kUNIC58TsQl8AABZLSURBVNxO51yFpJclndGI7YOcUB1Ds88L5IS0ygmSv7zgMydIjZAXoli81/S1282SmeVL+oqCUY7GbDfTzN6V9JmkOc65Rm0/9FtJP5ZU5aFtKUhMz5vZQgu+JbKx9ZNUJul/w4+JZ5hZOw9xSMG9vmc1dqPOubWSfi1ptaT1Cu4n/nwjhrBY0vFm1tXM2iq4PWJOHfsgtcgJcZpxXvCdEyS/eSGdcoLkIS+kQU6QGiEvRLF4T/prt5s6M2sv6S+SrnDObW3Mtp1zlc65YQq+OfFIMxvSmO2b2dclfeacW9iY7cY51jk3XNIESdPM7PhGbj9LwUf19zrnviJphyQf831bSjpV0v95aPsQBaOsfSX1ltTOzL7VWO0755ZJulXSHEnPKpiyUdFY7UMSOeELmmteSJOcIPnNC2mREyR/ecF3TpAaJy9EsXhP+mu3mzIza6Ggg/6zc+5xX3GEH8nNkzS+kZs+VtKpZvaJgo/Jx5nZnxozAOfcuvDfzyT9VcHH942pVFJpzOjWbAUdd2ObIOkd59ynHtr+qqSVzrky59w+SY9LOqYxA3DO/cE5N9w5d7yCqRMfNmb7ICdUa+Z5wXtOkLznhXTJCZK/vOA9J0gNnxeiWLwn89XcTZqZmYI5bcucc7/x0H53M+sc/txGwX+WksaMwTn3U+dcH+dcvoK/gRedc4327trM2plZh+qfJZ2k4KOyRuOc2yBpjZkNDBedIGlpY8YQOk8epsyEVks6yszahv8vTlAw17fRmFmP8N9cSZPk71w0V80+J0jkBd85QfKfF9IoJ0j+8oL3nCA1fF5I6htW04mr4au5GzMGM5slaaykbmZWKuk659wfGjGEYyV9W9L74fxCSbrGBd962Bh6SXrQzDIVvAF8zDnn5bZcHvWU9Negb1CWpIedc896iOP7kv4cFi0fq5G/hj6cz3eipIsbs91qzrk3zGy2pHcUfCz5LzX+t+r9xcy6StonaZpzbnMjt9+skRP2Iy/4lw55wWtOkPzmhTTJCVID5wW+YRUAAACIiChOmwEAAACaJYp3AAAAICIo3gEAAICIoHgHAAAAIoLiHQAAAIgIinc0CWZWaWbvmtliM/u/8FZVB7L/DDMrOoDtp5jZ3QceKQCgoZET0JRRvKOp2OWcG+acGyJpr6RLkt3RzDKdcxc553x9mQUAILXICWiyKN7RFL0q6TBJMrNvmdmb4QjM/eEXiMjMtpvZL8zsDUlHm9k8MxsZrjvPzN4PR2xurX5RM7vAzD4ws5cVfCEKACD9kRPQpFC8o0kxsyxJExR8y2ChpHMkHeucGyapUtI3w03bSVrsnBvlnJsfs39vSbdKGidpmKQjzOx0M+sl6QYFHfSJkpL+OBUA4Ac5AU1Rlu8AgBRpE/OV4K9K+oOkqZJGSHor/LrqNpI+C7eplPSXBK9zhKR5zrkySTKzP0s6PlwXu/xRSQMa4DgAAPVHTkCTRfGOpmJXOJKynwW984POuZ8m2H63c64ywXKrpQ1XnwABAI2GnIAmi2kzaMrmSppsZj0kycy6mFleHfu8IWmMmXUL50KeJ+nlcPlYM+tqZi0kndWQgQMAUo6cgCaBkXc0Wc65pWb2M0nPm1mGpH2SpklaVcs+683sp5JeUjDi8rRz7u+SZGbXS3pd0npJ70jKbNgjAACkCjkBTYU5x6c+AAAAQBQwbQYAAACICIp3AAAAICIo3gEAAICIoHgHAAAAIoLiHQAAAIgIincAAAAgIijeAQAAgIigeAcAAAAiguIdAAAAiAiKdwAAACAiKN4BAACAiKB4BwAAACKC4h0AAACICIp3AAAAICKyfAeAA7Nw4cIeWVlZMyQNEW++kH6qJC2uqKi4aMSIEZ/5DgbRRD8HoBlKOn9SvEdMVlbWjEMPPbSwe/fumzMyMpzveIBYVVVVVlZWVrRhw4YZkk71HQ+iiX4OQHNzIPmTEY3oGdK9e/etJDSko4yMDNe9e/ctCkZMgYNFPwegWTmQ/EnxHj0ZJDSks/Dvk74F9UE/B6DZSTZ/kmBxwDIzM0cMGjSoqPqxfPnylq+88krbKVOm5NS0zz/+8Y8O//Ef/3FYonXnnHNO3sKFC1s3XMTpq23btl+JfX7nnXd2Pf/883N9xVOXI488cuArr7zS1nccQEObPXt2x/z8/CG5ublDrrnmmkN9x5MqK1asaDFq1KgB/fr1G3zYYYcNvvHGG3v4jimVKioqVFhYWFRTvomq8vLyzPHjx/fr27fv4H79+g1+4YUX2vmOKRVuuOGGHocddtjggoKCwRMnTuy7c+dO8x1TKpx11ln5Xbp0KS4oKBhcvezTTz/NPOaYYwry8vKGHHPMMQVlZWWZB/v6zHmPuPyrnxqRytf75JZTFta1TatWrapKSkqWxi4bOHDg3uOPP37nwbT56KOPrjqY/VLu+k4pPZe6fkud5xJAEnZuSe3/zbadav2/WVFRoSuvvDL3ueee+6Bfv377iouLC88888zPR4wYsTuVYVQtmpvS48ooPqHOPqdFixa6/fbbS0ePHr1z8+bNGV/5yleKTj755K2pPrbK3/0kpceWeemtSfWnN910U8/DDjts1/bt2w+6MKrN7m+MTelxtX54XlLHNXXq1JyTTjpp67PPPvvx7t27bfv27SkdfH0vLz+lx3X4qk/qPK6VK1e2mD59es/ly5cvbt++vTv55JP7zZgxo8sPfvCDjamM5Z723VJ6bNO2l9d5bBdeeGH55Zdf/tkFF1zQt3rZdddd12vs2LHbbr755g+vueaaQ6+99tpD77333rUHEwMj70iJ2JH1p556qn31qHxhYWHR5s2bMyRpx44d+0cOTj311L5VVVWSvjia27Zt2698//vfzx44cGBRcXHxoDVr1mRJ0pIlS1oVFxcPGjJkSOEVV1zRO37Euil6+OGHOx1++OGDCgsLi4455pgB1efiqquu6j1p0qT8Y489tiA7O3vogw8+2PmSSy7pM2DAgKLjjjuuYM+ePSZJ2dnZQy+77LLsYcOGDRoyZEjh/Pnz244ePbogJydnyG233dZdkqqqqnTxxRf3KSgoGDxgwICi3//+94dUt/+zn/2s54ABA4oGDhxYdOmll2bHxlZZWalJkybl/+AHP+jdmOcEaAzz5s1rl5eXt6eoqGhv69at3aRJkzbNnj27s++4UiEvL2/f6NGjd0rSIYccUtW/f/9dq1evbuk7rlT46KOPWjz33HOdvve975X7jiWVNm3alPHGG290uOKKK8olqXXr1q5bt26VvuNKhcrKStuxY0fGvn37tGvXrow+ffrs8x1TKkyYMGF79+7dK2KXPfvss50vvvjijZJ08cUXb3zmmWcOSbx33SjeccD27NmTUV2cn3jiif3j199+++2H3nnnnatKSkqWLliwoKR9+/ZVkrRs2bI299xzz5oVK1YsWb16das5c+a0j993165dGUcfffT25cuXLz366KO333XXXd0l6bLLLsu59NJLP1u8ePGy3r17N4n/3NIXz+WgQYOKfvWrX+0vhk888cTt7777bsmyZcuWTp48edMvfvGL/R/dr1q1qtWLL764Yvbs2SsuueSSvuPGjdv6wQcfLG3dunXVY4891ql6u5ycnL3vvvtuyahRo7ZfeOGF+U8++eRHb7zxRsktt9zSW5Ieeuihzu+//36bZcuWLZk7d+4H1157bZ9Vq1a1eOyxxzo+9dRThyxcuLBk+fLlS6+77roN1a+5b98+O/300/sWFBTsvvPOO9c11rkCGsuaNWtaZmdn761+3qdPn71r165tEgVurOXLl7dcunRp2zFjxmz3HUsqTJs2Lee2224rzchoWqVNSUlJqy5dulScddZZ+YWFhUXnnHNO3tatWyN/kH379t03bdq0DX379j28R48exR06dKicNGnSVt9xNZSNGzdm5eXl7ZOCN9GbNm066Nkvkf/lo/FVT5spKSlZOmfOnI/i1x911FHbf/jDH+bcdNNNPcrLyzNbtGghSRo6dOiO/v3778vMzNTgwYN3fvTRR19Khi1atHDnnnvuFkkaMWLEjlWrVrWUpH/961/tL7zwwk2SdNFFF6X0IzWfYs9lSUnJ0p/+9Kf7i+GVK1e2PO644woGDBhQdOeddx5aUlLSpnrdV7/61S2tWrVyRx555K7KykqbPHnyVkkaPHjwrpUrV+4/r2efffbnkjR06NCdw4cP33HIIYdU9e7du6JVq1ZV5eXlma+++mqHs88+e1NWVpZycnIqRo0atX3+/Plt58yZ0/Fb3/pWeYcOHaokqWfPnvtHeS699NK8oqKiXbfeeuv+gh5oSpz78rWyZtakLqDdsmVLxqRJk/rfcssta7p06VLlO576mjVrVqdu3bpVHHfccQc1fTOdVVRU2LJly9pOmzatbNmyZUvbtm1b9fOf/zzy12GUlZVlPvXUU51XrFjx/oYNG97buXNnxu9+97suvuOKAop3pNzNN9+8YcaMGat27dqVccwxxxT+61//ai1JrVq12p/8MjMzVVFR8aULU7Kyslz1qElWVlbCbZqLyy67LPfSSy/97IMPPlh69913r9qzZ8/+/6/V5zIzM/ML5ywjI+ML56x169auennLli33n/+MjAzt27fPEhUpUlC8mCU+9SNHjtz+6quvdmwqFxYB8XJzc78w0l5aWtqyiX3iZ6ecckr/s846a9N3vvOdz33Hkwrz589vP2fOnM7Z2dlDp0yZ0m/BggUdTjvttL5175n+8vPz9/bs2XPvuHHjdkjSOeecs3nRokWRv3HAk08+2TE3N3dPOKDkTj/99M//+c9/fukT+aaia9euFatWrWohSatWrWrRpUuXirr2qQnFO1JuyZIlrY488shdv/zlLzcMHTp0x+LFi+t9J5lhw4Ztnzlz5iGS9MADDzSLd+bbtm3LzM3N3SdJM2fO7NoQbYwZM2bb7Nmzu1RUVGjdunVZb775Zvvjjjtux/jx47f+8Y9/7LZt27YMKbhKvnqfiy++uPykk07a8vWvf73/vn1Npp4B9hszZsyOTz75pHVJSUnL3bt32+OPP97lzDPPbBJFblVVlc4999y8AQMG7L7++us/9R1Pqtxzzz1rP/300/fWrl37/syZMz8+6qijtv39739f6TuuVMjNza049NBD9y5atKiVJD3//PMdBw4cmNILjH3Iz8/f+84777Tftm1bRlVVlV588cUOhYWFkT+umnzta1/7/P777+8qSffff3/X8ePHH3SfQvGOlLvtttt6FBQUDB44cGBRmzZtqiZPnrylvq951113rbnrrrt6Dh06tHD9+vUt2rdv3yQu1qnNf/3Xf60777zz+o8YMWJg165dD/odem2+/e1vfz548OBdhYWFg8eOHTvghhtuKM3Nza2YPHny1gkTJnw+bNiwwkGDBhXdeOONX/iI9vrrr/+0uLh456RJk/pWVjb5XwWamfCOLKvHjx8/oKCgYPDpp5++aeTIkU2iqJgzZ077v/3tb13nz5/fofpam0cffbRT3XvCp7vuumv1N7/5zX4DBgwoeu+999rcdNNN633HVF/jxo3bMXHixM2HH3544cCBAwdXVVXZVVddVeY7rlSYOHFi39GjRw9auXJlq549ex5+xx13dLvhhhvWv/TSSx3z8vKGvPTSSx1vuOGGg/4d1vixOdLTokWLPikuLm5SV9InY9u2bRnt2rWrysjI0PTp0w959NFHu8ydO/dL8+2RHhYtWtStuLg433cciKbm2s8BQDL5k/u8IxJee+21tpdffnmuc04dO3asnDlz5ie+YwIAAGhsFO+IhPHjx29fvnz50rq3BAAAaLqY8w4AAABEBMV79FRVVVVxiz6krfDvM/L3jQYAIB1RvEfP4rKysk4U8EhHVVVVVlZW1knSYt+xAADQFDHnPWIqKiou2rBhw4wNGzYMEW++kH6qJC2uqKi4yHcgAAA0RRTvETNixIjPJJ3qOw4AaMrOOuus/Llz53bq2rVrxYcffrjEdzypsnPnThs1atSgvXv3WmVlpU2cOHHzHXfcsc53XKmQnZ09tF27dpUZGRnKyspyixcvXuY7plRYtGhRq3POOad/9fPS0tJWP/7xj9dee+21n/mMKxVuvPHGHg899FB355zOP//8sqZwTFLi/uOBBx445Oabb+798ccft543b96y448/fufBvj7FOwAgrbn1H45I5etZr4KFdW1z4YUXll9++eWfXXDBBX1T2Xasyqemp/S4Mk+ZWudxtW7d2s2fP395p06dqvbs2WNHHHHEwLlz52454YQTdqQyloqfnpfSY8v61aw6j02SXn755Q969erVIF9qJ0kbjx6S0uPq+vriOo+ruLh4T0lJyVJJqqio0KGHHlp87rnnpvQbf5/p1julxzWhfF2dx/XWW2+1fuihh7q/8847y1q3bl01ZsyYAWecccaWoUOH7kllLJdYx5Qe231u60H1H8OGDdv1l7/8ZcX3vve9/PrGwLQLAADiTJgwYXv37t0brAj0JSMjQ506daqSpL1791pFRYWZcQlVVDzxxBMdc3Nz9wwYMGCv71jq6/33328zfPjw7R06dKhq0aKFjj322G2PPvpoZ99xpUKi/mP48OG7i4uLU/LGhOIdAIBmpKKiQoMGDSrq2bNn8ZgxY7aOGzcupaPuPp1wwgkFgwcPLvz1r3/dzXcsDWHWrFldJk+evNF3HKkwbNiwXW+88UaHDRs2ZG7bti1jzpw5ndasWdPSd1xRwLQZAACakaysLJWUlCwtLy/PPOWUU/q/9dZbrY844ojdvuOqr9dee60kPz9/39q1a7PGjRs3YPDgwbsnTJiw3XdcqbJ792574YUXOv3mN78p9R1LKgwfPnz35ZdfvmHcuHED2rZtW1VUVLQzK4uyNBmMvAMA0Ax169atcvTo0duefPLJTr5jSYX8/Px9kpSdnV1xyimnfP7666+38x1TKs2ePbtTUVHRzpycnCYznevKK68sX7p06bK33357eZcuXSoLCgoi/yayMVC8AwDQTKxbty6rvLw8U5K2b99u8+bN61hYWBj5gmnr1q0Zmzdvzqj++aWXXup4+OGH7/IdVyo98sgjXc4+++xNvuNIpbVr12ZJ0ocfftjyqaee6vzd7363SR1fQ+HzCQAA4kycOLHvggULOmzevDmrZ8+eh1999dXrrrzyynLfcdXXmjVrWkyZMqVvZWWlnHN22mmnbTrvvPO2+I6rvkpLS7POOOOMwySpsrLSzjzzzI2TJ0/e6juuVNm2bVvG/PnzOz744IOrfMeSSqeeemr/zz//PCsrK8v99re/Xd29e/dK3zGlQqL+o2vXrhU/+tGPcjdv3px1xhlnFBQWFu6cP3/+hwfz+uacS3XMAAActEWLFn1SXFwc+UIZAA7UokWLuhUXF+fXtg3TZgAAAICIoHgHAAAAIoLiHQAAAIgIincAQLqpqqqq4ms/ATQrYb9XVdd2FO8AgHSzuKysrBMFPIDmoqqqysrKyjpJWlzXttwqEgCQVioqKi7asGHDjA0bNgwRg0wAmocqSYsrKiouqmtDbhUJAAAARAQjGgAAAEBEULwDAAAAEUHxDgAAAEQExTsAAAAQERTvAAAAQET8f1OkBRW7kv1TAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, axs = plt.subplots(1, 2, figsize=(12.8, 4.8))\n", "\n", "(df.groupby(\"Period\").Choice.value_counts(normalize=True).unstack()\n", " .plot.bar(ax=axs[0], stacked=True, rot=0, title=\"Choice Probabilities\"))\n", "(df.groupby(\"Period\").Experience_Fishing.value_counts(normalize=True).unstack().plot\n", " .bar(ax=axs[1], stacked=True, rot=0, title=\"Share of Experience Level per Period\", cmap=\"Reds\"))\n", "\n", "axs[0].legend([\"Fishing\", \"Hammock\"], loc=\"upper center\", bbox_to_anchor=(0.5, -0.15), ncol=2)\n", "axs[1].legend(loc=\"upper center\", bbox_to_anchor=(0.5, -0.15), ncol=6)\n", "\n", "plt.show()\n", "plt.close()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "One can clearly see the proportions of experiences in the first period.\n", "\n", "## Lagged Choices\n", "\n", "In the next step, we want to enrich the model with information on the choice in the previous period. Let us assume that Robinson has troubles storing food on a tropical island so that he incurs a penalty if he enjoys his life in the hammock two times in a row. We add the parameter ``\"not_fishing_last_period\"`` to the nonpecuniary reward and add a corresponding covariate to the options." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "params = \"\"\"\n", "category,name,value\n", "delta,delta,0.95\n", "wage_fishing,exp_fishing,0.01\n", "nonpec_hammock,constant,1\n", "nonpec_hammock,not_fishing_last_period,-0.5\n", "shocks_sdcorr,sd_fishing,1\n", "shocks_sdcorr,sd_hammock,1\n", "shocks_sdcorr,corr_hammock_fishing,0\n", "initial_exp_fishing_0,probability,0.33\n", "initial_exp_fishing_1,probability,0.33\n", "initial_exp_fishing_2,probability,0.34\n", "\"\"\"\n", "\n", "options = {\n", " \"n_periods\": 10,\n", " \"simulation_agents\": 1_000,\n", " \"covariates\": {\n", " \"constant\": \"1\",\n", " \"not_fishing_last_period\": \"lagged_choice_1 != 'fishing'\"\n", " }\n", "}" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
value
categoryname
deltadelta0.95
wage_fishingexp_fishing0.01
nonpec_hammockconstant1.00
not_fishing_last_period-0.50
shocks_sdcorrsd_fishing1.00
sd_hammock1.00
corr_hammock_fishing0.00
initial_exp_fishing_0probability0.33
initial_exp_fishing_1probability0.33
initial_exp_fishing_2probability0.34
\n", "
" ], "text/plain": [ " value\n", "category name \n", "delta delta 0.95\n", "wage_fishing exp_fishing 0.01\n", "nonpec_hammock constant 1.00\n", " not_fishing_last_period -0.50\n", "shocks_sdcorr sd_fishing 1.00\n", " sd_hammock 1.00\n", " corr_hammock_fishing 0.00\n", "initial_exp_fishing_0 probability 0.33\n", "initial_exp_fishing_1 probability 0.33\n", "initial_exp_fishing_2 probability 0.34" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "params = pd.read_csv(io.StringIO(params), index_col=[\"category\", \"name\"])\n", "params" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "C:\\Users\\tobia\\git\\respy\\respy\\pre_processing\\model_processing.py:470: UserWarning: The distribution of initial lagged choices is insufficiently specified in the parameters. Covariates require 1 lagged choices and parameters define 0. Missing lags have equiprobable choices.\n", " category=UserWarning,\n", "C:\\Users\\tobia\\git\\respy\\respy\\pre_processing\\model_processing.py:470: UserWarning: The distribution of initial lagged choices is insufficiently specified in the parameters. Covariates require 1 lagged choices and parameters define 0. Missing lags have equiprobable choices.\n", " category=UserWarning,\n" ] } ], "source": [ "simulate = rp.get_simulate_func(params, options)\n", "df = simulate(params)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The warning is raised because we forgot to specify what individuals in period 0 have been doing in the previous period. As a default, **respy** assumes that all choices have the same probability for being the previous choice. Note that, it might lead to inconsistent states where an individual should have accumulated experience in the previous period, but still starts with zero experience.\n", "\n", "Below we see the choice probabilities on the left-hand-side and the shares of lagged choices on the right-hand-side. Without further information, **respy** makes all previous choices in the first period equiprobable.\n", "\n", "If we had set the covariate with the lagged choice but not added a parameter using the covariate, **respy** would have discarded the covariate and created a model without lagged choices." ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAu8AAAFhCAYAAADA0URkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deXxU9b3/8fcnCTuyCHEhgKAmgQCCgriBW9Vi25/2WtyXqhcRFbd6a63X22q1rdrS24vaohWL2uJStV735WqLWnfqRthEBUFAwQXZRJJ8fn+cEzqOk2SAyZz5Jq/n4zEPZuYs3885JN/PJ9/zPTPm7gIAAABQ+IqSDgAAAABAdijeAQAAgEBQvAMAAACBoHgHAAAAAkHxDgAAAASC4h0AAAAIBMU7toiZXW5mf9rCbaeY2X/lOqZcMbN+ZuZmVrKF27uZ7drAshPN7IlM6zZ1XszsUjO7eUtiAhAOMzvVzJ5LOo5MzGx7M3vGzFab2aSk49kSZjbNzK7aiu0b7OOb2K6vma0xs+ItbTt0W5P/zezvZjYu1zGFiOIdDTKzE8zs1bizWWZmj5rZqK3dr7tPcPcrcxFjvfiPiY1xrJ+Z2fNmtk8u28gFd/+zux/WwLJN58XMDjSzJWnLf+HudFxAC2Bmo+J+apWZfWJm/zCzPZOOKwvjJa2U1MXdL0pfuLWFcSEwsx3NbGqc91ab2Vwzu8LMOm3Nft39fXfv7O61uYq1uaUMZq2JHwvN7JIt3V9z5P/WiOIdGZnZDyT9VtIvJG0vqa+k30k6Msm4mnCXu3eWVCrpOUn3mZmlr7SlI+oAkAtm1kXSQ5Kuk7StpDJJV0ja0Axt5bq/20nSbG+h3/BoZttKekFSB0n7uPs2kg6V1E3SLknG1tya+FnpFufX4yX9xMzGbMH+W+0Vh1yjeMfXmFlXST+TdI673+fua919o7s/6O4/TFm1rZndFo9MVJvZiJR9DIwvcX0WLzsiZdlXRmbM7Egze93MPjezd+o7BTPrmjL68YGZXZXNL7+7b5R0q6QdJPWIL0H/w8z+28w+kXS5mRWZ2WVmtsjMPoqPo2vark43s6Vx+5tGmMxspJm9EB/bMjO73szapm37LTN718xWmtmvzKwo3rbBy+H15yUe3XlUUq+U0Y5eljZVycz2jkfuPjOzN8zswJRlp8btrzaz98zsxKbOG4C8qZAkd7/D3Wvdfb27P+Hub6auZGa/NrNP49/hw1PeP83M5sS/3++a2Zkpyw40syVm9iMzWy7pj/H734n72fork7s1FJyZ7Wtmr8RXBV4xs33j96dJ+r6ki+N+6ZDNOWgz+x8zWxz39TPNbHTKsg5mdmt8vHPM7GJLufpoZnuY2WvxMf/FzO5KyyMNHp+Z7W5m/4y3vUtS+0bC/IGk1ZJOcveFkuTui939/LT/n0PM7O043hvMooGixnKLpU3JNLNtzeyPcZ751Mzuz+Z4MpxXN7PzMuWcePnp8Tn91MweN7Od0rY9x8zelvR2I+dF8bl4QVK1pMHx9gPM7EmLrh7NM7NjUvY9zcx+b2aPmNlaSQfZ1/P/GWa2IN7+ATPrlbLsUIuueqwys+slfW0wrtVydx48vvKQNEZSjaSSRta5XNIXkr4lqVjSLyW9GC9rI2mBpEsltZV0sKLOsDJePk3SVfHzkZJWKRrZKFI0AjUgXna/pBsldZK0naSXJZ3ZSDx/ip+3k/QrSYvj16fGx3OupBJFIyqnxzHuLKmzpPsk3R6v30+SS7ojbnuIpBWSDomXD5e0d7yvfpLmSLogJRaX9DdFI2p9Jc2XNC4llufS1t01w3k5UNKSRo6xTNLH8fkvis/fx4quOnSS9HnK+d5R0qCkf6548OARPSR1iX9fb5V0uKTuactPlbRR0hlx/3qWpKWSLF7+bUWjwCbpAEnrJO0RLzsw7u+uifvCDpL2kPSRpL3i/X1f0kJJ7TLEtq2kTyWdHPdxx8eve8TLN/VTDRxbg8slnSSpR7zfiyQtl9Q+Xna1pBmSukvqLenN+j5QUR5ZJOl8RfnlKElfpvSXDR5fyrYXxtuOjc9tQzG+KOmKJv7/XNGVk26K+vgVksbEy7LJLSXx64cl3RUfcxtJBzR1PI3E01DO+W4cz8D4vF8m6fm0bZ+Mt+2QYd+bYlb087afop+3byjKNYslnRYv30PRlKpBKT8Lq+JtihT90bTp50NRbbAy3q6doitRz8TLeirKY2Pjc3Ohop/rcUn//hbCI/EAeBTeQ9KJkpY3sc7lkv4v5XWVpPXx89GKOuWilOV3SLo8fp76y3ujpP/OsP/tFV1C7pDy3vGS/tZIPF9K+izu9J6WNDxedqqk99PWf0rS2SmvKxV16PUFuSv+IyJefq2kqQ20fYGkv6a8dsUdefz6bElPpcSSi+L9R4oTQsryxxV18p3i8/A9ZeiMefDgkfxDUTE1TdKSuCh5QNL28bJTJS1IWbdj3Ffs0MC+7pd0fvz8wLgvbJ+y/PeSrkzbZp7iYjHt/ZMlvZz23guSTo2fb+qnGoil0eVp634qaWj8/F1J30xZNk7/Kt73l/SB4j9e4veeS+kvGzy+eNulads+31CMikafJzQRt0salfL6bkmXxM+zyS0ligZV6pT2h9vm/n+lxNNQznlU0r+nLCtSVHzvlLLtwY0ca33Mn8X/X3MknRcvO1bSs2nr3yjppyk/C7c19PMhaaqka1OWdY7PVT9JpygeEIyXmaLfFYp3d6bNIKOPJfW0pudKLk95vk5S+3ibXopGvetSli9SNFqcro+kdzK8v5Oiv7aXxZcNP1PUKWzXSDx3u3s3d9/O3Q9295kpyxanrdsrjik1vhJFfzRk2mZRvI3MrMLMHjKz5Wb2uaL7Anqm7T/jtjm0k6Sj689NfH5GSdrR3dcq6lQnKDp/D5vZgBy3D2AruPscdz/V3XsrmoLQS9F9RvWWp6y7Ln7aWZLM7HAzezGeavCZoitwqX3QCnf/IuX1TpIuSusv+ihzv5TeN0oN99+bxcwuiqdvrIpj6JoSdy99td9Mfd5L0gceV3EZljd2fJm2TT++VB8rKqybkp7/OqfE2lRuURzfJ+7+aYZ9b87/V72Gcs5Okv4nZT+fKCqEyxrYtiE93b27uw9098kp+94rLc4TFU1ZzWbfXzlX7r5G0fkvU9rPQ/z/l02crQLFOzJ5QdGUmO9u4fZLJfVJnXOn6FLeBxnWXazMNwEtVjTy3jMuyLu5exd3H7SFMXna66WKOp7U+GokfZjyXp+05Uvj57+XNFdSubt3UTQ9KH0uXkPbbmm86RYrGnnvlvLo5O5XS5K7P+7uhypKQnMl/WEz2weQJ+4+V9GI5OCm1jWzdpLulfRrRSP13SQ9oq/2Qen9x2JJP0/rLzq6+x0ZmkjvG6WG+++sxfPbfyTpGEWjzd0UTamoj3uZouky9VL70GWSyurnlWdY3tjxZdq2byOh/p+kf0vLX5sjm9xSH/O2ZtYtwz425/+rXkM5Z7Gi6aap++rg7s+nrN9UvmnIYkkz0vbd2d3PynLfXzlXFt3v1UPRz9qy1GOK///6pO+gtaJ4x9e4+ypJP5F0g5l918w6mlmbeLTn2ix28ZKktYpuampj0Y2U/0/SnRnWnSrpNDP7RnyjT5mZDXD3ZZKekDTJzLrEy3YxswNyc5S6Q9KFZtbfzDorGj2/y91rUtb5r/jYByma03dX/P42iubirYlHtFM7qno/NLPuZtZH0TzNuzKs05gPFd1sm34Tbb0/Sfp/ZvZNMys2s/YW3ajW26LPYT4i7gg3SFojKZiPJgNauvgmv4vMrHf8uo+iaYEvZrF5W0Xzg1dIqrHoRtaMHz+b4g+SJpjZXhbpZGbfNrNtMqz7iKQKiz4quMTMjlU0LfKhLA9Pkur7pPpHW0X9Zk0cd4mZ/UTR3P96d0v6cdxvlkmamLLsBUV92MQ4piMV3S+VzfG9ELd7XrztUWnbpvtNHNetFt/YGeel31gjN42myCa3KM5xj0r6XXzMbcxs/yyOpyEN5Zwpis7roPhYuprZ0VkcRzYeUvSzcnIcfxsz29PMBma5/XRF+X9Y/EfpLyS95NGNwg9LGmRmR1l0Rf88fXVEv1WjeEdG7v4bRXfdX6aos12sqDO9v7Ht4m2/lHSEohuxVir6iMlT4tGl9HVfVlQY/7eiUZgZ+tdf4qcoSlSzFc21u0fZXc7Mxi2Sbpf0jKT3FF1pODdtnRmKbvR5StKv3b3+y5X+Q9IJim7C/YMyF+b/K2mmpNcVdUJTNye4+FzdIend+HJkr7TlixV9bOel+tf/zw8V/U4XKboZbKmiS6QHKJoDCaAwrFZ0M+JLFn0Kx4uSZin6vW2Uu69WVMjcrahfPEHRfPnGtnlV0c2v18fbLFA0rz7Tuh9L+k4cy8eSLpb0HXdfmcVx1btE0vqUx9OK7sl5VNHNlIsU9bmp0yB+pmhO83uKRr/vUfzRmXFOOUrSvyuae32SosKxfnmDx5ey7anxsmMV3USakbt/ImlfRXOvXzKz1YpywKp4v03JJrfUOzluZ66ie7UuaOp4GpEx57j7XxXdvHynRdM8ZynKzVst/lk8TNJxivLNcv3rRulstn9K0n8pupK0TNFV+OPiZSslHa3oRuaPJZVL+kcu4m4J6u9cBwAAKAhmdpak49w949VWM3tJ0hR3/2N+Iys8ZuaKpnFm88cFWgBG3gEAQKIs+lbT/eIpkpWKRv7/mrL8ADPbIZ768n1Ju0l6LKl4gSTxTZMAACBpbRV9olh/RVNj7lQ05bJepaKpQp0VfULZ2HjeONDqMG0GAAAACATTZgAAAIBAULwDAAAAgUhsznvPnj29X79+STUPAAVl5syZK929NOk4kkReAIBIYzkhseK9X79+evXVV5NqHgAKipk19pXtrQJ5AQAijeUEps0AAAAAgaB4BwAAAAJB8Q4AAAAEguIdAAAACATFOwAAABAIincAAAAgEBTvAAAAQCCaLN7N7BYz+8jMZjWw3MxsspktMLM3zWyP3IcJACgU5AUASE42I+/TJI1pZPnhksrjx3hJv9/6sAAABWyayAsAkIgmi3d3f0bSJ42scqSk2zzyoqRuZrZjrgIEABQW8gIAJKckB/sok7Q45fWS+L1l6Sua2XhFozDq27dv43u9vOvWR3b5qq3cfitjCL39Qogh6fYLIYak2y+EGJJuv1BiCEdh5oWkfw4LIYak2y+EGEJvvxBiSLr9QoghwfZzccOqZXjPM63o7je5+wh3H1FaWpqDpgEABYi8AADNJBfF+xJJfVJe95a0NAf7BQCEibwAAM0kF8X7A5JOiT9dYG9Jq9z9a5dGAQCtBnkBAJpJk3PezewOSQdK6mlmSyT9VFIbSXL3KZIekfQtSQskrZN0WnMFCwBIHnkBAJLTZPHu7sc3sdwlnZOziAAABY28AADJ4RtWAQAAgEBQvAMAAACBoHgHAAAAAkHxDgAAAASC4h0AAAAIBMU7AAAAEAiKdwAAACAQFO8AAABAICjeAQAAgEBQvAMAAACBoHgHAAAAAkHxDgAAAASC4h0AAAAIBMU7AAAAEAiKdwAAACAQFO8AAABAICjeAQAAgEBQvAMAAACBoHgHAAAAAkHxDgAAAASC4h0AAAAIBMU7AAAAEAiKdwAAACAQFO8AAABAICjeAQAAgEBQvAMAAACBoHgHAAAAAkHxDgAAAASC4h0AAAAIBMU7AAAAEAiKdwAAACAQFO8AAABAICjeAQAAgEBQvAMAAACBoHgHAAAAAkHxDgAAAASC4h0AAAAIBMU7AAAAEAiKdwAAACAQFO8AAABAILIq3s1sjJnNM7MFZnZJhuVdzexBM3vDzKrN7LTchwoAKATkBABITpPFu5kVS7pB0uGSqiQdb2ZVaaudI2m2uw+VdKCkSWbWNsexAgASRk4AgGRlM/I+UtICd3/X3b+UdKekI9PWcUnbmJlJ6izpE0k1OY0UAFAIyAkAkKBsivcySYtTXi+J30t1vaSBkpZKekvS+e5el5MIAQCFhJwAAAnKpni3DO952utvSnpdUi9JwyRdb2ZdvrYjs/Fm9qqZvbpixYrNDhYAkLic5QSJvAAAmyub4n2JpD4pr3srGk1JdZqk+zyyQNJ7kgak78jdb3L3Ee4+orS0dEtjBgAkJ2c5QSIvAMDmyqZ4f0VSuZn1j284Ok7SA2nrvC/pG5JkZttLqpT0bi4DBQAUBHICACSopKkV3L3GzCZKelxSsaRb3L3azCbEy6dIulLSNDN7S9El1R+5+8pmjBsAkAByAgAkq8niXZLc/RFJj6S9NyXl+VJJh+U2NABAISInAEBy+IZVAAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQiJKkA2hIvy+mb/U+Fm59GAAASCIvASgMjLwDAAAAgaB4BwAAAAJB8Q4AAAAEguIdAAAACATFOwAAABCIgv20GUDa+k93WNhCYgAAFIZCyAmFEAOSw8g7AAAAEAiKdwAAACAQTJsBUPC4RAwAKCRJ5iWK9wJGwQIAKCR8yyyQPKbNAAAAAIFg5L0RjHyjEPBzCKAe/QEKAT+HyWLkHQAAAAgExTsAAAAQCKbNoFFcGgMAFBLyElo7Rt4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCD4qEkCT+Gg2AEAhac15iZF3AAAAIBCMvANAFlrzKA8AoHBkNfJuZmPMbJ6ZLTCzSxpY50Aze93Mqs1sRm7DBAAUCnICACSnyZF3MyuWdIOkQyUtkfSKmT3g7rNT1ukm6XeSxrj7+2a2XXMFDABIDjkBAJKVzcj7SEkL3P1dd/9S0p2Sjkxb5wRJ97n7+5Lk7h/lNkwAQIEgJwBAgrIp3sskLU55vSR+L1WFpO5m9nczm2lmp+QqQABAQSEnAECCsrlh1TK85xn2M1zSNyR1kPSCmb3o7vO/siOz8ZLGS1Lfvn03P1oAQNJylhMk8gIAbK5sRt6XSOqT8rq3pKUZ1nnM3de6+0pJz0gamr4jd7/J3Ue4+4jS0tItjRkAkJyc5QSJvAAAmyub4v0VSeVm1t/M2ko6TtIDaev8r6TRZlZiZh0l7SVpTm5DBQAUAHICACSoyWkz7l5jZhMlPS6pWNIt7l5tZhPi5VPcfY6ZPSbpTUl1km5291nNGTgAIP/ICQCQrKy+pMndH5H0SNp7U9Je/0rSr3IXGgCgEJETACA5WX1JEwAAAIDkUbwDAAAAgaB4BwAAAAJB8Q4AAAAEguIdAAAACATFOwAAABAIincAAAAgEBTvAAAAQCAo3gEAAIBAULwDAAAAgaB4BwAAAAJB8Q4AAAAEguIdAAAACATFOwAAABAIincAAAAgEBTvAAAAQCAo3gEAAIBAULwDAAAAgaB4BwAAAAJB8Q4AAAAEoiTpAAAAyEa/L6Zv1fYLcxMGACSKkXcAAAAgEBTvAAAAQCAo3gEAAIBAULwDAAAAgaB4BwAAAAJB8Q4AAAAEguIdAAAACATFOwAAABAIincAAAAgEBTvAAAAQCAo3gEAAIBAULwDAAAAgaB4BwAAAAJB8Q4AAAAEguIdAAAACATFOwAAABAIincAAAAgEBTvAAAAQCAo3gEAAIBAULwDAAAAgciqeDezMWY2z8wWmNkljay3p5nVmtnY3IUIACgk5AQASE6TxbuZFUu6QdLhkqokHW9mVQ2sd42kx3MdJACgMJATACBZ2Yy8j5S0wN3fdfcvJd0p6cgM650r6V5JH+UwPgBAYSEnAECCsineyyQtTnm9JH5vEzMrk/RvkqY0tiMzG29mr5rZqytWrNjcWAEAyctZTojXJS8AwGbIpni3DO952uvfSvqRu9c2tiN3v8ndR7j7iNLS0mxjBAAUjpzlBIm8AACbqySLdZZI6pPyurekpWnrjJB0p5lJUk9J3zKzGne/PydRAgAKBTkBABKUTfH+iqRyM+sv6QNJx0k6IXUFd+9f/9zMpkl6iE4aAFokcgIAJKjJ4t3da8xsoqJPDCiWdIu7V5vZhHh5k3MaAQAtAzkBAJKVzci73P0RSY+kvZexg3b3U7c+LABAoSInAEBy+IZVAAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQiKyKdzMbY2bzzGyBmV2SYfmJZvZm/HjezIbmPlQAQCEgJwBAcpos3s2sWNINkg6XVCXpeDOrSlvtPUkHuPtukq6UdFOuAwUAJI+cAADJymbkfaSkBe7+rrt/KelOSUemruDuz7v7p/HLFyX1zm2YAIACQU4AgARlU7yXSVqc8npJ/F5D/l3So1sTFACgYJETACBBJVmsYxne84wrmh2kqKMe1cDy8ZLGS1Lfvn2zDBEAUEBylhPidcgLALAZshl5XyKpT8rr3pKWpq9kZrtJulnSke7+caYduftN7j7C3UeUlpZuSbwAgGTlLCdI5AUA2FzZFO+vSCo3s/5m1lbScZIeSF3BzPpKuk/Sye4+P/dhAgAKBDkBABLU5LQZd68xs4mSHpdULOkWd682swnx8imSfiKph6TfmZkk1bj7iOYLGwCQBHICACQrmznvcvdHJD2S9t6UlOfjJI3LbWgAgEJETgCA5PANqwAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgKN4BAACAQFC8AwAAAIGgeAcAAAACQfEOAAAABILiHQAAAAgExTsAAAAQCIp3AAAAIBAU7wAAAEAgsirezWyMmc0zswVmdkmG5WZmk+Plb5rZHrkPFQBQCMgJAJCcJot3MyuWdIOkwyVVSTrezKrSVjtcUnn8GC/p9zmOEwBQAMgJAJCsbEbeR0pa4O7vuvuXku6UdGTaOkdKus0jL0rqZmY75jhWAEDyyAkAkKBsivcySYtTXi+J39vcdQAA4SMnAECCzN0bX8HsaEnfdPdx8euTJY1093NT1nlY0i/d/bn49VOSLnb3mWn7Gq/oEqokVUqatxWx95S0ciu2z4WkY0i6/UKIIen2CyGGpNsvhBiSbj8XMezk7qW5Cqa55DInxMtaUl5Iuv1CiCHp9gshhqTbL4QYkm6/EGJotpxQksXGSyT1SXndW9LSLVhH7n6TpJuyaLNJZvaqu4/Ixb5CjSHp9gshhqTbL4QYkm6/EGJIuv1CiSFPcpYTpJaVF5JuvxBiSLr9Qogh6fYLIYak2y+EGJqz/WymzbwiqdzM+ptZW0nHSXogbZ0HJJ0Sf8LA3pJWufuyHMcKAEgeOQEAEtTkyLu715jZREmPSyqWdIu7V5vZhHj5FEmPSPqWpAWS1kk6rflCBgAkhZwAAMnKZtqM3P0RRZ1x6ntTUp67pHNyG1qTcnKZdSslHUPS7UvJx5B0+1LyMSTdvpR8DEm3LxVGDHlRoDlBSv7/IOn2peRjSLp9KfkYkm5fSj6GpNuXko+h2dpv8oZVAAAAAIUhq29YBQAAAJC8IIv3pr6aOw/t32JmH5nZrHy3Hbffx8z+ZmZzzKzazM7Pc/vtzexlM3sjbv+KfLafFkuxmb1mZg8l0PZCM3vLzF43s1fz3X4cQzczu8fM5sY/D/vkse3K+NjrH5+b2QX5aj8ljgvjn8NZZnaHmbXPc/vnx21XJ3H8ICfEMZAXlGxOiNtPNC8kmRPi9hPPC0nnhDiG5s0L7h7UQ9ENUu9I2llSW0lvSKrKcwz7S9pD0qyEzsGOkvaIn28jaX4+z4Ekk9Q5ft5G0kuS9k7oXPxA0nRJDyXQ9kJJPZM47pQYbpU0Ln7eVlK3hOIolrRc0efS5rPdMknvSeoQv75b0ql5bH+wpFmSOiq6h+j/JJUn+TPR2h7khE0xkBc82ZwQt59oXiiUnBC3n/e8kHROiNts9rwQ4sh7Nl/N3azc/RlJn+SzzbT2l7n7P+PnqyXNUR6/vdAja+KXbeJH3m+eMLPekr4t6eZ8t10IzKyLoqJhqiS5+5fu/llC4XxD0jvuviiBtkskdTCzEkWdZcbPE28mAyW96O7r3L1G0gxJ/5bH9kFOqI+h1ecFckJB5QQpubyQZE6Q8pAXQize+drtFGbWT9LuikY58tlusZm9LukjSU+6e17bj/1W0sWS6hJoW4oS0xNmNtOib4nMt50lrZD0x/gy8c1m1imBOKTos77vyHej7v6BpF9Lel/SMkWfJ/5EHkOYJWl/M+thZh0VfTxinya2QW6RE9K04ryQdE6Qks0LhZQTpATyQgHkBCkPeSHE4t0yvNcqPzLHzDpLulfSBe7+eT7bdvdadx+m6JsTR5rZ4Hy2b2bfkfSRZ/i69Tzaz933kHS4pHPMbP88t1+i6FL97919d0lrJSUx37etpCMk/SWBtrsrGmXtL6mXpE5mdlK+2nf3OZKukfSkpMcUTdmoyVf7kERO+MH4fBsAABC5SURBVIrWmhcKJCdIyeaFgsgJUnJ5IemcIOUnL4RYvGf9tdstmZm1UdRB/9nd70sqjviS3N8ljclz0/tJOsLMFiq6TH6wmf0pnwG4+9L4348k/VXR5ft8WiJpScro1j2KOu58O1zSP939wwTaPkTSe+6+wt03SrpP0r75DMDdp7r7Hu6+v6KpE2/ns32QE+q18ryQeE6QEs8LhZITpOTyQuI5QWr+vBBi8Z7NV3O3aGZmiua0zXH33yTQfqmZdYufd1D0yzI3nzG4+4/dvbe791P0M/C0u+ftr2sz62Rm29Q/l3SYoktleePuyyUtNrPK+K1vSJqdzxhixyuBKTOx9yXtbWYd49+Lbyia65s3ZrZd/G9fSUcpuXPRWrX6nCCRF5LOCVLyeaGAcoKUXF5IPCdIzZ8XsvqG1ULiDXw1dz5jMLM7JB0oqaeZLZH0U3efmscQ9pN0sqS34vmFknSpR996mA87SrrVzIoV/QF4t7sn8rFcCdpe0l+jvkElkqa7+2MJxHGupD/HRcu7yvPX0Mfz+Q6VdGY+263n7i+Z2T2S/qnosuRryv+36t1rZj0kbZR0jrt/muf2WzVywibkheQVQl5INCdIyeaFAskJUjPnBb5hFQAAAAhEiNNmAAAAgFaJ4h0AAAAIBMU7AAAAEAiKdwAAACAQFO8AAABAICje0SKYWa2ZvW5ms8zsL/FHVW3O9jebWdVmrH+qmV2/+ZECAJobOQEtGcU7Wor17j7M3QdL+lLShGw3NLNidx/n7kl9mQUAILfICWixKN7REj0raVdJMrOTzOzleATmxvgLRGRma8zsZ2b2kqR9zOzvZjYiXna8mb0Vj9hcU79TMzvNzOab2QxFX4gCACh85AS0KBTvaFHMrETS4Yq+ZXCgpGMl7efuwyTVSjoxXrWTpFnuvpe7P5eyfS9J10g6WNIwSXua2XfNbEdJVyjqoA+VlPXlVABAMsgJaIlKkg4AyJEOKV8J/qykqZLGSxou6ZX466o7SPooXqdW0r0Z9rOnpL+7+wpJMrM/S9o/Xpb6/l2SKprhOAAAW4+cgBaL4h0txfp4JGUTi3rnW939xxnW/8LdazO8b4204VsTIAAgb8gJaLGYNoOW7ClJY81sO0kys23NbKcmtnlJ0gFm1jOeC3m8pBnx+weaWQ8zayPp6OYMHACQc+QEtAiMvKPFcvfZZnaZpCfMrEjSRknnSFrUyDbLzOzHkv6maMTlEXf/X0kys8slvSBpmaR/Sipu3iMAAOQKOQEthblz1QcAAAAIAdNmAAAAgEBQvAMAAACBoHgHAAAAAkHxDgAAAASC4h0AAAAIBMU7AAAAEAiKdwAAACAQFO8AAABAICjeAQAAgEBQvAMAAACBoHgHAAAAAkHxDgAAAASC4h0AAAAIBMU7AAAAEIiSpAPA5pk5c+Z2JSUlN0saLP74QuGpkzSrpqZm3PDhwz9KOhigtSA3oMCQC5oRxXtgSkpKbt5hhx0GlpaWflpUVORJxwOkqqursxUrVlQtX778ZklHJB0P0FqQG1BIyAXNi7/OwzO4tLT0czpnFKKioiIvLS1dpWj0D0D+kBtQMMgFzYviPTxFdM4oZPHPJ30LkF/kBhQUckHz4aRisxUXFw8fMGBAVf1j3rx5bZ955pmOp556ap+GtnnooYe2Oeigg3bNtOzYY4/daebMme2bL+LC1bFjx91TX0+ePLnHKaec0jepeJoycuTIymeeeaZj0nEAKDxXXXXVdjvvvPOgLl26DLv00kt3aGi9xvq5Aw44YNeVK1cWN1+UhW3evHlty8vLByUdRzbKysqGLFu2jOnXCeCkB67fJQ8Pz+X+Fl797ZlNrdOuXbu6uXPnzk59r7Ky8sv9999/3Za0eddddy3aku1y7vKuOT2XunxVk+cSAJpDErlh6tSppY8++ujbAwYM+HJL25kxY8aCLd0258gJKFCMvCMnUkfWH3744c71o/IDBw6s+vTTT4skae3atcVjxozZuX///oOOOOKI/nV1dZK+OprbsWPH3c8999yyysrKqqFDhw5YvHhxiSRVV1e3Gzp06IDBgwcPvOCCC3qlj1i3RNOnT++62267DRg4cGDVvvvuW1F/Ln7wgx/0Ouqoo/rtt99+5WVlZUNuvfXWbhMmTOhdUVFRNXr06PINGzaYFI2KTJw4sWzYsGEDBg8ePPC5557rOGrUqPI+ffoMvvbaa0slqa6uTmeeeWbv8vLyQRUVFVV/+MMfute3f9lll21fUVFRVVlZWXX22WeXpcZWW1uro446qt95553XK5/nBEBhOuGEE/ouWbKk3RFHHLHrFVdcsV39yPott9zSvby8fFBlZWXViBEjKuvXX758eZvRo0eX77TTToMnTJjQu/79+tHcefPmtd15550HHXfccTvtuuuug/bbb7/yNWvWmCTNmDGjY0VFRdWwYcMG1Pdf+T/i5lNbW6v04540aVLPwYMHD6ysrKz65je/ucvq1auLJOl73/tevxNPPLHvXnvtVdG7d+8hDz/8cOejjz6638477zzoe9/7Xr/6fXbs2HH3s846q2zQoEED991334q//e1vHUeOHFnZu3fvIX/+85+7StK6dets7Nix/SoqKqoGDhxY9eCDD24jSTU1NRo/fnzvioqKqoqKiqqf//zn26XGu2bNGhs9enT5pEmTeubxNLVqFO/YbBs2bCiqL84PPfTQXdKXT5o0aYfJkycvmjt37uwXX3xxbufOneskac6cOR1uuOGGxQsWLKh+//332z355JOd07ddv3590T777LNm3rx5s/fZZ5811113XakkTZw4sc/ZZ5/90axZs+b06tVrY/MfZX6knssBAwZU/fKXv9xUDB966KFrXn/99blz5syZPXbs2E9+9rOfbboMvWjRonZPP/30gnvuuWfBhAkT+h988MGfz58/f3b79u3r7r777q716/Xp0+fL119/fe5ee+215vTTT+/34IMPvvPSSy/Nvfrqq3tJ0m233dbtrbfe6jBnzpzqp556av5PfvKT3osWLWpz9913d3n44Ye7z5w5c+68efNm//SnP11ev8+NGzfad7/73f7l5eVfTJ48eWm+zhWAwjV9+vT3t9tuu40zZsyY371799r696+++uodn3jiifnz5s2b/dhjj20aVZ89e3bH+++//905c+ZUP/DAA90XLFjQJn2f77//fvvzzjvvowULFlR37dq19rbbbusuSePGjet/ww03LHr99dfnFhcXt7h5/pmO+8QTT/x01qxZc+bNmze7srJy/eTJkzcVyqtWrSp54YUX5l999dWLjz322PIf/vCHH7799tvVc+fO7fD88893kKLcetBBB62urq6e06lTp9rLLrus7Nlnn53/l7/8ZcGVV15ZJknXXHPNdpI0f/782dOnT393/Pjx/datW2eTJk0qXbRoUbvq6urZ8+fPnz1u3LiP69v+/PPPiw477LDyY4899pOLLrpoZb7PVWtF8Y7NVj9tZu7cubOffPLJd9KX77333mv+4z/+o89VV1213cqVK4vbtIn65CFDhqzdZZddNhYXF2vQoEHr3nnnnbbp27Zp08aPO+64VZI0fPjwtYsWLWorSa+99lrn008//RNJSu04Qpd6LufOnTv7xz/+8aZi+L333ms7evTo8oqKiqrJkyfvMHfu3A71yw455JBV7dq185EjR66vra21sWPHfi5JgwYNWv/ee+9tOq/HHHPMZ5I0ZMiQdXvsscfa7t271/Xq1aumXbt2dStXrix+9tlntznmmGM+KSkpUZ8+fWr22muvNc8991zHJ598sstJJ520cptttqmTpO23335TMj777LN3qqqqWn/NNddsKugBIJMRI0asOfHEE/tNmjSpZ01Nzab3R40a9XmPHj1qO3bs6LvuuusX77zzTrv0bcvKyjbsu+++6yVp9913X7dw4cJ2K1euLF67dm3RoYceulaSvv/973+St4PJk0zHPXPmzA7Dhw+vrKioqLr33nt7VFdXb7pP7Nvf/vZnRUVF2mOPPdb16NFj48iRI9cXFxeroqJiff15bdOmjafmiVGjRq2uzyEffPBBW0l6/vnnO59yyikfx+1+0atXry/feuut9k8//XSXCRMmrKjP5an54Igjjtj15JNPXjlx4sQWk5dDQPGOnPvFL36x/Oabb160fv36on333Xfga6+91l6S2rVrt2mEpLi4WDU1NZa+bUlJiRcVFdU/z7hOazFx4sS+Z5999kfz58+fff311y/asGHDpt/X+nNZXFz8lXNWVFT0lXPWvn17r3+/bdu2m85/UVGRNm7caO6ZB63cXWaZT/2IESPWPPvss13WrVvXav9vAGRn+vTp71911VVLFy9e3HbYsGGDli9fXizpK/1RcXGxb9y48Wv9Sfo6NTU1DfZZLUmm4x4/fnz/66+//v358+fP/tGPfrQ0NR/U9/PFxcVf6+fr80F6nkjNIbW1tSZF/X4mcT7IuHDPPfdc89hjj3WtnwaL/KB4R85VV1e3Gzly5Pqf//zny4cMGbJ21qxZW/1JMsOGDVszbdq07pJ0yy23bLv1URa+1atXF/ft23ejJE2bNq1Hc7RxwAEHrL7nnnu2ramp0dKlS0tefvnlzqNHj147ZsyYz2+//fae9fMqP/zww02f/nDmmWeuPOyww1Z95zvf2WXjxhYzgwlAM6iurm538MEHr/3tb3+7tHv37jXvvvvu1664bo7S0tLaTp061T311FOdJOn2229vFflg3bp1RX379t24YcMGu/POO5vlmEeNGrXmT3/607aS9Oabb7ZbtmxZ29122+2LQw455PMpU6aU1vf3qfngV7/61dJtt9225uSTTy7YT0lriSjekXPXXnvtdvU3KHXo0KFu7Nixq7Z2n9ddd93i6667bvshQ4YMXLZsWZvOnTvXNr1V2P7zP/9z6fHHH7/L8OHDK3v06FHT9Bab7+STT/5s0KBB6wcOHDjowAMPrLjiiiuW9O3bt2bs2LGfH3744Z8NGzZs4IABA6quvPLKr3zs2+WXX/7h0KFD1x111FH9a2tb/H8FgC104YUX9q6oqKgqLy8ftPfee6/ee++912/tPm+88caFZ5111k7Dhg0b4O7aZpttWnwndMkllywdOXLkwNGjR1eUl5d/0RxtXHzxxR/V1tZaRUVF1bHHHrvLjTfeuLBDhw5+4YUXrujdu/eXAwYMGFRZWVk1derUr/zxMHXq1MUbNmwoSr3xGM2rVVyCakneeOONhUOHDm11N4WsXr26qFOnTnVFRUW66aabut91113bPvXUU1+bb4/C8MYbb/QcOnRov6TjAFqL1pQbVq1aVdS1a9c6Sbr00kt3WLZsWZs//vGPi5OOC19HLmgefM47gvCPf/yj4/nnn9/X3dWlS5faadOmLUw6JgBA/t19991dJ02atGNtba2VlZVtmD59+sKkYwLyieIdQRgzZsyaefPmzW56TQBAS3bGGWd8esYZZ3yadBxAUpjzDgAAAASC4j08dXV1dXxEHwpW/PPJ54YB+UVuQEEhFzQfivfwzFqxYkVXOmkUorq6OluxYkVXSbOSjgVoZcgNKBjkgubFnPfA1NTUjFu+fPnNy5cvHyz++ELhqZM0q6amZlzSgQCtCbkBBYZc0Iz4qEgAAAAgEPx1DgAAAASC4h0AAAAIBMU7AAAAEAiKdwAAACAQFO8AAABAIP4/wbRKFz5/3NEAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, axs = plt.subplots(1, 2, figsize=(12.8, 4.8))\n", "\n", "(df.groupby(\"Period\").Choice.value_counts(normalize=True).unstack()\n", " .plot.bar(ax=axs[0], stacked=True, rot=0, title=\"Choice Probabilities\"))\n", "(df.groupby(\"Period\").Lagged_Choice_1.value_counts(normalize=True).unstack().plot\n", " .bar(ax=axs[1], stacked=True, rot=0, title=\"Share of Lagged Choice per Period\"))\n", "\n", "axs[0].legend([\"Fishing\", \"Hammock\"], loc=\"upper center\", bbox_to_anchor=(0.5, -0.15), ncol=2)\n", "axs[1].legend(loc=\"upper center\", bbox_to_anchor=(0.5, -0.15), ncol=2)\n", "\n", "plt.show()\n", "plt.close()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "What if we wanted a more complex distribution of lagged choices. Remember that Robinson starts with 0, 1 or 2 periods of experience in fishing. Thus, it is natural to assume that higher experience levels correspond to a higher probability for having fishing as a previous choice.\n", "\n", "This kind of distribution is more complex because it involves covariates. A flexible solution is to use a **multinomial logit** or **softmax function** to map parameters with covariates to probabilities. This allows to specify very complex distributions, ones which could have been estimated outside the structural model. Here is a short explanation.\n", "\n", "For each lag, we have a set of parameters, $\\beta^f$ and $\\beta^h$, and their corresponding covariates, $x^f$ and $x^h$. The probability $p^f$ for one individual having fishing as their first lagged choice is\n", "\n", "$$\n", " p^f = \\frac{e^{x^f \\beta^f}}{e^{x^f \\beta^f} + e^{x^h \\beta^h}}\n", "$$\n", "\n", "The probability for hammock as the lagged choice is defined analogously.\n", "\n", "In the parameters, you can use the ``lagged_choice_1_fishing`` keyword to define the parameters for fishing as the first lag. You can also define higher order lags using, for example, ``lagged_choice_2_*``.\n", "\n", "Let us assume that with each level of initial experience the probability for choosing fishing the previous period rises from 0.5 over 0.75 to 1. For zero experience in fishing, the resulting probabilities should be one halve for fishing and one halve for hammock. First, we create a covariate which is ``True`` if an individual has zero experience in fishing. This covariate is ``{\"zero_exp_fishing\": \"exp_fishing == 0\"}``. Then, the parameters for this covariate just have to be equal and can take any value. That is because the softmax function is shift-invariant. The resulting lines in the csv are\n", "\n", " lagged_choice_1_fishing,zero_exp_fishing,1\n", " lagged_choice_1_hammock,zero_exp_fishing,1\n", " \n", "For one period of experience, the weights are 0.75 for fishing and 0.25 for the hammock. First, define a covariate ``{\"one_exp_fishing\": \"exp_fishing == 1\"}``. To get the corresponding softmax coefficients for the probabilities, note that you can rewrite the softmax formula replacing $x_i$ with 1. Recognize that the sum in the denominator, $C$, applies to all probabilities and can be discarded. Thus, the coefficients are the logs of the probabilities.\n", "\n", "$$\\begin{align}\n", " p_i &= \\frac{e^{x_i \\beta_i}}{\\sum_j e^{x_j \\beta_j}} \\\\\n", " &= \\frac{e^{\\beta_i}}{\\sum_j e^{\\beta_j}} \\\\\n", " log(p_i) &= \\beta_i - \\log(\\sum_j e^{\\beta_j}) \\\\\n", " &= \\beta_i - C\n", "\\end{align}$$\n", "\n", "The lines in the parameters for $log(0.75)$ and $log(0.25)$ are\n", "\n", " lagged_choice_1_fishing,one_exp_fishing,-0.2877\n", " lagged_choice_1_hammock,one_exp_fishing,-1.3863\n", " \n", "For two experience in fishing, we have to make sure that fishing receives a probability of one which can be achieved by using a fairly large value such as 6. By default, all missing choices receive a parameter with value -1e300.\n", "\n", "The final parameters and options are the following:" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [], "source": [ "params = \"\"\"\n", "category,name,value\n", "delta,delta,0.95\n", "wage_fishing,exp_fishing,0.01\n", "nonpec_hammock,constant,1\n", "nonpec_hammock,not_fishing_last_period,-0.5\n", "shocks_sdcorr,sd_fishing,1\n", "shocks_sdcorr,sd_hammock,1\n", "shocks_sdcorr,corr_hammock_fishing,0\n", "initial_exp_fishing_0,probability,0.33\n", "initial_exp_fishing_1,probability,0.33\n", "initial_exp_fishing_2,probability,0.34\n", "lagged_choice_1_fishing,zero_exp_fishing,1\n", "lagged_choice_1_hammock,zero_exp_fishing,1\n", "lagged_choice_1_fishing,one_exp_fishing,-0.2877\n", "lagged_choice_1_hammock,one_exp_fishing,-1.3863\n", "lagged_choice_1_fishing,two_exp_fishing,6\n", "\"\"\"\n", "\n", "options = {\n", " \"n_periods\": 10,\n", " \"simulation_agents\": 1_000,\n", " \"covariates\": {\n", " \"constant\": \"1\",\n", " \"not_fishing_last_period\": \"lagged_choice_1 != 'fishing'\",\n", " \"zero_exp_fishing\": \"exp_fishing == 0\",\n", " \"one_exp_fishing\": \"exp_fishing == 1\",\n", " \"two_exp_fishing\": \"exp_fishing == 2\",\n", " }\n", "}" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
value
categoryname
deltadelta0.9500
wage_fishingexp_fishing0.0100
nonpec_hammockconstant1.0000
not_fishing_last_period-0.5000
shocks_sdcorrsd_fishing1.0000
sd_hammock1.0000
corr_hammock_fishing0.0000
initial_exp_fishing_0probability0.3300
initial_exp_fishing_1probability0.3300
initial_exp_fishing_2probability0.3400
lagged_choice_1_fishingzero_exp_fishing1.0000
lagged_choice_1_hammockzero_exp_fishing1.0000
lagged_choice_1_fishingone_exp_fishing-0.2877
lagged_choice_1_hammockone_exp_fishing-1.3863
lagged_choice_1_fishingtwo_exp_fishing6.0000
\n", "
" ], "text/plain": [ " value\n", "category name \n", "delta delta 0.9500\n", "wage_fishing exp_fishing 0.0100\n", "nonpec_hammock constant 1.0000\n", " not_fishing_last_period -0.5000\n", "shocks_sdcorr sd_fishing 1.0000\n", " sd_hammock 1.0000\n", " corr_hammock_fishing 0.0000\n", "initial_exp_fishing_0 probability 0.3300\n", "initial_exp_fishing_1 probability 0.3300\n", "initial_exp_fishing_2 probability 0.3400\n", "lagged_choice_1_fishing zero_exp_fishing 1.0000\n", "lagged_choice_1_hammock zero_exp_fishing 1.0000\n", "lagged_choice_1_fishing one_exp_fishing -0.2877\n", "lagged_choice_1_hammock one_exp_fishing -1.3863\n", "lagged_choice_1_fishing two_exp_fishing 6.0000" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "params = pd.read_csv(io.StringIO(params), index_col=[\"category\", \"name\"])\n", "params" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "C:\\Users\\tobia\\git\\respy\\respy\\state_space.py:84: UserWarning: Some choices in the model are not admissible all the time. Thus, respy applies a penalty to the utility for these choices which is -400000 by default. For the full solution, the penalty only needs to be larger than all other value functions to be effective. Choose a milder penalty for the interpolation which does not dominate the linear interpolation model.\n", " \"Some choices in the model are not admissible all the time. Thus, respy\"\n" ] } ], "source": [ "simulate = rp.get_simulate_func(params, options)\n", "df = simulate(params)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The parameters produce previous choice probabilities on the left-hand-side and the correct conditional probabilities on the right-hand-side. The differences in the probabilities occur due to the low number of simulated individuals." ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAuIAAAFhCAYAAAA1L8/UAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdeZhcVZn48e/bHbaELZCEJSGELezIEvZNUFYRVFAQFGTEiAOCo6OiP0dlRAdUUBhAjIIoIyCIIpuAooDsiyCQhCWyJSyGQFhCAibd7++Puh0rnV6qO119u5rv53nuk6p7zz3nvdXV6bdOnXtOZCaSJEmS+ldT2QFIkiRJ70Qm4pIkSVIJTMQlSZKkEpiIS5IkSSUwEZckSZJKYCIuSZIklcBEvGQR8YmIuK3sODoSEatFxK0R8UZEnF52PL0RERdGxClLcH5GxPq9OG9sRMyJiObett3oIuK8iPivXp57c0Qc09cxSZI0kJiI94OI2CUi7oiI1yLilYi4PSK2LTuuGkwEZgErZuYX2h9c0iR3IIiINSLi/Ih4ofjA8WhEnBwRw5ak3sx8NjOXz8yWvoq13iJiXPHBY06xPR0RJ/W2vsw8NjO/1ZcxSpI0mJiI11lErAhcA/wvsAowGjgZeLsObQ3p4yrXBqbkIF31KSJWAe4ElgN2zMwVgL2AlYH1yoyt3rp5r6ycmcsDHwW+HhH79qL+d+w3AZIk1cpEvP7GA2TmJZnZkpnzMvPGzHyoulBEfD8iZkfEUxGxX9X+oyNiatFb+2REfLrq2LsjYkZEfDkiXgR+Vuw/ICIejIhXi574LToLLiJ2ioh7i976eyNip2L/hcBRwJeK3tH39uSiI+LMiJgeEa9HxP0RsWvVseUi4ufF9U6NiC9FxIyq41tHxAPFNV8eEb+q7nnv6voiYquI+Gtx7q+AZbsI8/PAG8DHMvNpgMycnpkntvv5vDcinijiPSciomirKSK+FhHPRMTMiPhFRKxUHGvrXR5SPF8lIn4WEc8X9VxZy/V08LpmRJxQvBdmRcT3IqKp6vi/Fa/p7Ii4ISLWbnfucRHxBPBEF68LxWtxJzAZ2Kw4f6OI+EPxrc5jEfGRqrovjIgfRcR1EfEmsEf7b0wi4lMRMa04/6qIWLPq2F7FtxGvRcTZQHQXnyRJjc5EvP4eB1qKxHO/iBjeQZntgceAEcB3gfPbkj1gJnAAsCJwNPCDiNi66tzVqfS0rw1MLI5dAHwaWBX4MXBVRCzTvtGiR/ha4Kyi7BnAtRGxamZ+Avgl8N1iiMUfe3jd9wJbFrFdDFweEW1J8TeAccC6VHqgP1YV09LAb4ELi3MvAT5YdbzT6yvOvRK4qDj3cuDgLmJ8L/CbzGzt5loOALYF3gV8BNin2P+JYtujuJblgbM7qeMiYCiwKTAK+EF319NFPB8EJgBbAwcB/1bU9QHgq8CHgJHAX6i8ftU+QOX9tklXFxwVOxfxPhCVoTp/oPKzHEWlt/zciNi06rTDgW8DKwC3tatvT+B/qLx+awDPAJcWx0YAVwBfo/I78Hdg567ikyRpMDARr7PMfB3YBUjgJ8BLRW/galXFnsnMnxTjiX9OJVFZrTj/2sz8e1bcAtwI7Fp1bivwjcx8OzPnAZ8CfpyZdxc98D+nMgxmhw7Cex/wRGZelJkLMvMS4FHg/X1w3f+XmS8X9Z4OLANsWBz+CPCdzJydmTOofBBoswMwBDgrM+dn5m+Ae6qOd3V9OwBLAT8szv01lQ8EnVkVeKGGyzk1M1/NzGeBP1P5gAFwBHBGZj6ZmXOArwCHRbthHxGxBrAfcGxxzfOLn2V319OZ0zLzlSKeH1JJiqGSzP9PZk7NzAXAd4Atq3vFi+OvFO+VzswCXgF+CpyUmTdR+TDydGb+rPiZ/pVK8nxI1Xm/y8zbM7M1M99qV+cRwAWZ+dfMfLt4rXaMiHHA/lSGQP06M+cX1/RiF/FJkjQomIj3gyIx+kRmjqHyNf+aVJKNNi9WlZ1bPFweoOhFv6v4Ov9VKknLiKpzX2qX9KwNfKEY5vBqcc5aRZvtrUmlZ7LaM1TGsS+RiPhCMUTitSKGlariXhOYXlW8+vGawHPtxqVXH+/q+jo6t/31VXuZyoee7lQnhXMpfjYs/vo9Q+VDRPWHLIr4XsnM2R3U3ZOfV5vq1+OZqrJrA2dW1fMKlSEeozs5tzMjMnN4Zm6cmW0fktYGtm8X5xFUvpGppe5FXqvig8vLRWyLvB+Kn18tcUqS1NBMxPtZZj5KZdjFZt2VLYYnXAF8H1gtM1cGrmPR8bPtb6ScDnw7M1eu2oYWvd3tPU8lwao2FniupovpPO5dgS9T6fkeXsT9WlXcLwBjqk5Zq+rxC8DoqqE57Y93dX0dnTu2i1D/CHyweox1D7V//cYCC4B/tCs3HVglIlbuoI6e/LzaVL8eY4s42ur6dLu6lsvMO6rK9/bG2+nALe3qXj4zP1Nj3Yu8VsVQl1WpvNdeqL6m4ue3VvsKJEkabEzE66y4we0LETGmeL4WlaEEd9Vw+tJUhnS8BCyIyk2ce3dzzk+AYyNi+2Kc77CIeF9ErNBB2euA8RFxeEQMiYhDqYwdvqbGywNojohlq7alqYwRXlDEPSQivk5ljHuby4CvRMTwiBgNHF917E6gBTi+iOkgYLsar+/Oot0TinM/1O7c9s4o4vp52/CNiBgdEWdEFzdMVrkE+I+IWCcilqcyFORXxbCQhTLzBeD3VMZUD4+IpSJitxqupzNfLOpZCzgR+FWx/zwqr+umxbWsFBEfruE6anENlffKx4v4l4qIbSNi4xrPvxg4OiK2LD5gfge4u7hJ9lpg04j4UDGs5wQW7WmXJGlQMhGvvzeo3Bx3d1Rmk7gLeARYbF7u9jLzDSpJyWXAbCo3w13VzTn3URl3fHZxzjQqNxR2VPZlKmN/v0BlmMCXgAMyc1YN19XmJGBe1fYn4AYqiefjVIYjvMWiQw3+G5gBPEWlV/rXFNM5ZuY/qdxs+EngVSo3cl5TdbzT66s69xPFsUOB33QWeGa+AuwEzKfy83kDuIlK7/20Gq79Aio3Yd5aXMtbwGc7Kfvxop1HqdyA+7nurqcLvwPuBx6kksSeX9T1W+A04NKIeJ3K+2y/zirpieK9uDdwGJXe7ReLtrq6qbT6/JuA/6LyDc8LVKaHPKw4Ngv4MHAqlffhBsDtfRG3JEkDWeTgnCJaDSQiPgMclpm7d3L8buC8zPxZ/0Y28EREAhtkZi0fFCRJ0gBmj7j6XVRWs9w5KvNwb0ilR/63Vcd3j4jVi+ElRwFbANeXFa8kSVI9mIirDEtTmS/7DSpDWX4HnFt1fEPgb1SGiHwBOKQYZy1JklSKiLggKgv4PdLJ8YiIs6KyeN1Dsei6Lx3X6dAUSZIkqWvFRAtzgF9k5mKz30XE/lTuFdufyv2BZ2bm9l3VaY+4JEmS1I3MvJXKGh2dOYhKkp6ZeRewclQW9euUibgkSZK05Eaz6CxxM+hmkcQhXR2spxEjRuS4cePKal6SBpT7779/VmaOLDuOMs2f+bRjJd/hZi71jv4VUJXRw4dF96U6d2yM6/H/Jz/mmU8DE6t2TcrMST2ooqOYu4yjtER83Lhx3HfffWU1L0kDSkQ8U3YMkjRYNPcijc/WnAT0JPFubwaLrgw9hn+tft0hh6ZIkiRpUGmO6PHWB64CjixmT9kBeK27Wd9K6xGXJEmS6qE3PeLdiYhLgHcDIyJiBvANYCmAzDwPuI7KjCnTgLnA0d3VaSIuSZKkQaWPergXkZkf7eZ4Asf1pE4TcUmSJA0q9egRrwcTcUmSJA0q9egRrwcTcUmSJA0q9ohLkiRJJWiUHvFupy+MiAsiYmZEPNLJ8YiIsyJiWkQ8FBFb932YkiRJUm2aerGVoZZ2LwT27eL4fsAGxTYR+NGShyVJkiT1TknziPdYt0NTMvPWiBjXRZGDgF8UU7bcFRErR8Qa3U1gLkmSJNXDO2mM+GhgetXzGcW+xRLxiJhIpdecsWPHdl3rN1da8si++doSnr+EMTR6+wMhhrLbHwgxlN3+QIih7PYHSgySpJoMmjHiNejoSrOjgpk5KTMnZOaEkSNH9kHTkiRJ0qKao+dbGfoiEZ8BrFX1fAzwfB/UK0mSJA1afZGIXwUcWcyesgPwmuPDJUmSVJZBc7NmRFwCvBsYEREzgG8ASwFk5nnAdcD+wDRgLnB0vYKVJEmSujNobtbMzI92czyB4/osIkmSJGkJNMrNmq6sKUmSpEFl0PSIS5IkSY3ERFySJEkqgUNTJEmSpBLYIy5JkiSVwB5xSZIkqQT2iEuSJEklsEdckiRJKoE94pIkSVIJ7BGXJEmSStBkIi5JkiT1v2iQsSkm4pIkSRpUmkzEJUmSpP4XzU1lh1ATE3FJkiQNKg5NkSRJkkrg0BRJkiSpBNHk0BRJkiSp39kjLkmSJJXAMeKSJElSCZw1RZIkSSqBQ1MkSZKkEkSTibgkSQJuu/teTj3zPFpaWzj4gP045mOHLnL8ngf+xglf+Saj11gdgPfutjOfOfpjAOz94SMZNnQ5mpqaaG5u5rKfnt3v8atv3HPn7Zz9g+/T2trC/gd+kMOPPLrDco9OmczxxxzFf51yKrvv+V4APvqB9zF02LCF74PzLvxlf4becJocmiJJklpaWjjljHP4yQ/+h9VHjuDQT32WPXbegfXWWXuRcltvsRnnfvdbHdZxwZnfZfjKK/VHuKqTlpYWzvz+aXzvrHMZOWo1PnP0x9hp190Zt866i5WbdM6ZTNh+x8XqOOOcH7PSysP7K+SG5s2akqR3jIjYCDgIGA0k8DxwVWZOLTWwAeDhqY8xdvSarLXmGgDs955386fb7lwsEdfg9uiURxg9Zgxrjh4DwJ577cMdt968WCL+28svZbc93sNjU6aUEeag0SiJeGP020uSBqyI+DJwKRDAPcC9xeNLIuKkMmMbCGa+9DKrjxq58PlqI0cwc9asxcr9bfJUPvSJYzn2P/8f0556euH+CJj4+a/ykU8ex+VXXdcfIasOZr30EqNGrb7w+YhRo3jppZmLlHlp5kxuu+XPvP+Dhyx2fkTwxROO49NHHc41V15R93gbXVNzU4+3MtgjLklaUp8ENs3M+dU7I+IMYDJwamcnRsREYCLAud/7NscceXg94yxFkovtCxbtrdtk/Pr84fKLGDp0OW698x5O+OrJXHfJzwC46NwfMGrEqrw8+1U+9R8nsc7YtZiw5eb9Erv6Tmb374Nzfvh9Jh53As3NzYuVPWvSzxgxciSzX3mFL57wGdZaexzv2mqbusXb6BqlR9xEXJK0pFqBNYFn2u1fozjWqcycBEwCmD/z6cUzlUFgtZEjeHHmSwuf/+OlWYwcseoiZZYfNmzh49123I5Tzjib2a++xvCVV2JUUXbV4Svznt125uGpj5qIN6CRo0Yxc+aLC5/PmjmTESNHLlLm8alT+NbXvgLAa6+9yt133kZzczO77L7HwrLDV1mFXXbfg0enTDYR70JTg8ya4tAUSdKS+hxwU0T8PiImFdv1wE3AiSXHVrrNNtqQZ2c8x4znX2T+/Pn8/qab2WOXHRYpM+vlVxb2mD485VFaW1tZeaUVmTvvLd6cOxeAufPe4o5772eDdcf19yWoD2y08aY8N306Lzz/HPPnz+dPf7iBHXfdfZEyF//2Gi658louufJadt/jvZz4xa+wy+57MG/ePOa++SYA8+bN47577mKdddcr4zIaRjQ39Xgrgz3ikqQlkpnXR8R4YDsqN2sGMAO4NzNbSg1uABgypJmv/sdxfPoLX6WltZUPvm9v1l9nHL+68hoADv3AAdx481/41ZXX0NzczLLLLMP3vvkVIoKXZ8/mxK+eDFRm09h/rz3YZfttS7wa9VbzkCF89j+/zJdPPI6W1lb2O+BA1ll3Pa76za8BOPBDi48LbzP7lZf5+pe/AFTeB+/Ze1+223Hnfom7UbmgjyTpHSMzW4G7yo5joNptx+3YbcftFtl36AcOWPj48IMP4vCDD1rsvLXWXIPfXHhe3eNT/9hhp13YYaddFtnXWQL+5a+fvPDxmqPH8NP/+1VdYxtsGmWMuENTJEmSNKjUY2hKROwbEY9FxLSOZoSKiJUi4uqI+FtETI6IjldsqmKPuCRJkgaVvh6aEhHNwDnAXhRD7yLiqsysnvD9OGBKZr4/IkYCj0XELzPzn53VayIuSZKkQSX6ftaU7YBpmfkkQERcSmURs+pEPIEVIiKA5YFXgAVdVWoiLkmSpEGlDgv0jAamVz2fAWzfrszZwFVUVhZeATi0uH+mU44RlyRJ0qASzdHzLWJiRNxXtU2srrKDZtqvfbAP8CCVdRW2BM6OiBW7itMecUmSJA0qvZkXvHqBsQ7MANaqej6GSs93taOBU7OyKMC0iHgK2Ai4p7M27RGXJEnSoBJNTT3eunEvsEFErBMRSwOHURmGUu1Z4D0AEbEasCHwZFeV2iMuSZKkQaWvx4hn5oKIOB64AWgGLsjMyRFxbHH8POBbwIUR8TCVoSxfzsxZXdVrIi5JkqRBpR5L1mfmdcB17fadV/X4eWDvntRpIi5JkqRBpR6JeD3UFGU9VhKSJEmS6qEOY8Trotse8XqtJCRJkiTVQzQ3lx1CTWoZmlKXlYQkSZKkemiUoSm1JOJ1WUlIkiRJqoemkoaa9FQtiXhPVhLaE1gP+ENE/CUzX1+kosoKRRMBxo4d2/NoJUmSpG40So94LVHWupLQb7JiGtC2ktAiMnNSZk7IzAkjR47sbcySJElSp6K5qcdbGWpptS4rCUmSJEn1MGhmTanXSkKSJElSPTTK0JSaFvSpx0pCkiRJUj00SiLeGFFKkiRJg4xL3EuSJGlQaWqQHnETcUmSJA0qZd182VMm4pIkSRpUGmWM+IBNxMe9dfES1/H0kochSZKkBmMiLkmSJJXAoSmSJElSCZqam8sOoSYm4pIkSRpUHJoiSZIklcBEXJIkSSqBY8QlSZKkEtgjLkmSJJXARFySJEkqgUNTJEmSpBJEk9MXSpIkSf3PRFySJEkqgUNTJEmSpP4XrqwpSZIklcChKZKkRhURZ3Ww+zXgvsz8XX/HI0k9YiIuSWpgywIbAZcXzw8GJgOfjIg9MvNzfd3g3Kt/2tdVqsFctO5RZYegAeKkPTZYovOdvlCS1MjWB/bMzAUAEfEj4EZgL+DhMgOTpG7ZIy5JamCjgWFUhqNQPF4zM1si4u3ywpKkGpiIS4PDuLcuXqLzn27w9vWO9V3gwYi4GQhgN+A7ETEM+GOZgUlSdxyaIklqWJl5fkRcB2xHJRH/amY+Xxz+YnmRSVIN7BGXJDW4JuAlKn8r1o+I9TPz1pJjkqTumYhL0uDxThsiFBGnAYdSmSmltdidgIm4pAHPBX0GgbL/8Jbd/kCJQVIpPgBsmJnemCmp8ThGXJLUwJ4ElgJMxCU1HoemSBos/GbkHWkulVlTbqIqGc/ME8oLSZJqEybikqQGdlWxSVLjcWiKJKlRZebPy45BknrLHnGpDzgkQupfEXFZZn4kIh6mMkvKIjJzixLCkqSeMRGXJDWgE4t/Dyg1CklaEg5NkSQ1msx8ofj3mbJjkaTeapR5xBvj44IkqV9FxIci4omIeC0iXo+INyLi9bLjkqSaNDX3fOtGROwbEY9FxLSIOKmTMu+OiAcjYnJE3NJdnfaIS5I68l3g/Zk5texAJKnH+niMeEQ0A+cAewEzgHsj4qrMnFJVZmXgXGDfzHw2IkZ1V6+JuCSpI/8wCZfUqKLvx4hvB0zLzCcBIuJS4CBgSlWZw4HfZOazAJk5s7tKTcQlSQtFxIeKh/dFxK+AK1l0QZ/flBKYJPVEL3rEI2IiMLFq16TMnFQ8Hg1Mrzo2A9i+XRXjgaUi4mZgBeDMzPxFV22aiEsa8JzGsl+9v+rxXGDvqucJmIhLGvii5z3iRdI9qZPD0dEp7Z4PAbYB3gMsB9wZEXdl5uOdtWkiLklaKDOPLjsGSVpivUjEuzEDWKvq+Rjg+Q7KzMrMN4E3I+JW4F1Ap4m4s6ZIkhYTEd+NiBUjYqmIuCkiZkXEx8qOS5JqkdHU460b9wIbRMQ6EbE0cBhwVbsyvwN2jYghETGUytCVLu+1qSkRr8d0LZKkAW3vzHydysI+M6iMffxiuSFJUo2iqedbFzJzAXA8cAOV5PqyzJwcEcdGxLFFmanA9cBDwD3ATzPzka7q7XZoSr2ma5EkDWhLFf/uD1ySma9EdDREUpIGoDr8f5WZ1wHXtdt3Xrvn3wO+V2udtYwRr8t0LZKkAe3qiHgUmAf8e0SMBN4qOSZJqk2DLHFfS5QdTdcyul2Z8cDwiLg5Iu6PiCP7KkBJUv/LzJOAHYEJmTkfeJNKJ4wkDXh1GCNeF7X0iPfZdC3V8zOOHTu259FKkuoqIvbMzD9VzSdOuyEpTl8oaeArKbHuqVoS8T6brqV6fsYJEya0T+YlSeXbHfgTi84n3sZ5xCU1hkGUiC+crgV4jsp0LYe3K/M74OyIGAIsTWW6lh/0ZaCSpPrLzG8U/zqfuKTG1SCJeLdR1mu6FknSwBMRF1Y9PqrEUCSp1wbTGPG6TNciSRqQ3lX1+ETg52UFIkm91iA94i5xL0mq5v07khpfg6x7YCIuSao2JiLOojJjVtvjhTLzhHLCkqQesEdcktSAqpexv6+0KCRpCZQ15runTMQlSQtlZk1jwiPifzPzs/WOR5IGMxNxSVJv7Fx2AJLUqQZZ4t5EXJIkSYOLQ1MkSZKkEpiIS5IGscaYG0zSO5OJuCSp0UXEsMx8s4NDZ/Z7MJJUI2dNkSQ1rIjYCfgpsDwwNiLeBXw6M/8dIDMvLDG8hnPnky9y+k0P0JrJQVusy1E7bLTI8Yvufozrpz4DQEtr8vTLr3PD8Qex0nJL863f38ttf3+B4UOX4dJ/26eM8NVHZky+n7svm0S2tjJ+573ZYt8PL3L873f/mYdvvAKAIcssy06H/zurjFkXgMk3/Y7Hb78BEsbvsg+bvuegfo+/oZiIS5Ia2A+AfYCrADLzbxGxW7khNaaW1uS7f/wrZ39kN0atMJSjfvFHdl1/TdYdseLCMh/ffkM+vv2GAPxl2vNcfN/jrLTc0gC8b7NxfHir9fnmdfeUEr/6RmtrC3dd8iP2OfEUhg5flav/5z8Yu8X2rLzm2IVllh+xOvt9/lSWGbY8Mx65j9v/72zef9IZzH7uaR6//Qbef9IZNDUvxY3/+3XGbDaBlVYbXeIVDXANsrJmY3xckCT1u8yc3m5XSymBNLjJL7zCmJWXZ/TKy7NUcxN7b7wWt057rtPyN0x9ln02/ldytvVaI1mxSMrVuGY9/TgrjFqDFUauTvOQpVh329149qG7Fimz2nobs8yw5QEYuc5GzJ09C4BXX5zByHU2YsjSy9LU3MzqG2zGsw/e2e/X0FCiqedbCUzEJUkdmV4MT8mIWDoi/hOY2puKIuLovg2tsbw0Zx6rrTB04fNRKwzlpTfmdVj2rfkLuOupF9lj/Jj+Ck/9ZO7slxk2fOTC50NXHsGbs1/utPzjt9/I6M0mADB8zbX5xxOP8Nac11nwz7eY8ch9vFkk6epYRlOPtzKYiEuSOnIscBwwGpgBbFk8742TOzsQERMj4r6IuO/CW/7ay+oHtsxcfGcnX5v/ZdoLbDF6xMJhKRo8OngXEJ28D1547CGeuONGJnzwEwCsvMZabL7PIdxw5n9x41nfYJUx6xBNzfULdjBokB5xx4hLkhaTmbOAI2otHxEPdXYIWK2LdiYBkwBeO/9rHeUqDW/UCkP5xxtzFz6f+cZcRi6/bIdlb3z0WfauGpaiwWPY8FV5c/ZLC5/PfXUWQ1deZbFyr8x4itsvOou9Pnsyyy7/r/sIxu+8N+N33huA+6/8OUNXHlH/oBtYOkZcktSoIuLnEbFy1fPhEXFBF6esBhwJvL+DrfPv398BNlljONNnz+G5V99kfksrN06dzq7rr7lYuTlvz+eB6S+xewfH1PhGrD2e12c+zxuzXqRlwXyevPdW1tpi+0XKzHllJn/68XfY9egvLHYj5rzXX11Y5pkH7mTdbXfvt9gbUWbPtzLYIy5J6sgWmflq25PMnB0RW3VR/hpg+cx8sP2BiLi5DvE1jCFNTXzxvVtxwuW30prJ+zdfh/VGrMQVD/wdgIO3Wg+Amx9/ju3Hrc5ySy/6p/lrV93F/dNf4tV5b3PAudfwqV025aAt1un369CSaWpuZodDj+XGs75OtraywU57MXzNtXn01usA2Gi3/Xnw2kt5+83XueuScwGIpmYO/OoPAfjzpO/w1pw3KvV89NiFN3WqY61lZdY9ZCIuSepIU0QMz8zZABGxCl38zcjMT3Zx7PA6xNdQdl5vDXZeb41F9rUl4G0O2HwcB2w+brFzTzlwh3qGpn601ubbstbm2y6yb6Pd9l/4eJePn8AuHz+hw3P3/8/v1jW2waYx0nATcUlSx04H7oiIXxfPPwx8u8R4JKlmrQ2SiZuIS5IWk5m/iIj7gT2o3HD5ocycUnJYklSTDmcrGoBMxCVJnXkUmE3xtyIixmbms+WGJEnds0dcktSwIuKzwDeAf1BZUTOoDLvcosy4JKkWDZKHm4hLkjp0IrBhZr6jpx6U1JjsEZckNbLpwGtlByFJveEYcUlSI3sSuDkirgXebtuZmWeUF5Ik1aa17ABqZCIuSerIs8W2dLFJUsNokA5xE3FJ0uIy82SAiBiWmW+WHY8k9USjjBFvKjsASdLAExE7RsQUYGrx/F0RcW7JYUlSTTKzx1sZTMQlSR35IbAP8DJAZv4N2K3UiCSpRq292Mrg0BRJUocyc3pEVO9qKSsWSeoJx4hLkhrZ9IjYCciIWBo4gWKYiiQNdK0Nkok7NEWS1JFjgeOA0cAMYMviuSQNeNmLrQz2iEuSFpOZs4Ajyo5DknqjUWZNMRGXJC0UEV/KzO9GxP/SQSdRZp5QQliS1CMNMjLFRFyStIi2ceD3lRqFJC2B1tIGm/SMibgkaaHMvDoimoHNMvOLZccjSb3RKD3i3qwpSVpEZrYA25QdhyT1Vmv2fOtOROwbEY9FxLSIOPhZNJkAAB7ESURBVKmLcttGREtEHNJdnfaIS5I68kBEXAVcDixc4j4zf1NeSJJUm77uES++KTwH2IvKTFL3RsRVmTmlg3KnATfUUq+JuCSpI6tQWVVzz6p9CZiISxrw6jBGfDtgWmY+CRARlwIHAVPalfsscAWwbS2VmohLkhaTmUeXHYMk9VZvesQjYiIwsWrXpMycVDweDUyvOjYD2L7d+aOBD1LpwDARlyT1TkSMB34ErJaZm0XEFsCBmXlKyaFJUrd6s7JmkXRP6uRwdHRKu+c/BL6cmS0RHRVfXE03a9ZjcLokaUD7CfAVYD5AZj4EHFZqRJJUo5bWnm/dmAGsVfV8DPB8uzITgEsj4mngEODciPhAV5V22yNer8HpkqQBbWhm3tOuV2dBWcFIUk/0pke8G/cCG0TEOsBzVDomDq8ukJnrtD2OiAuBazLzyq4qrWVoSl0Gp0uSBrRZEbEexVevxTedL5QbkiTVpqWPE/HMXBARx1PpcG4GLsjMyRFxbHH8vN7UW0si3meD06sHwY8dO7ansUqS+s9xVMZKbhQRzwFPAUeUG5Ik1aYOPeJk5nXAde32dZiAZ+YnaqmzlkS8zwanVw+CnzBhQoOseSRJ7zzFt6DvjYhhQFNmvlF2TJJUqxrGfA8ItSTiPRmcDjAC2D8iFnQ3LkaSNDBFxKrAN4BdgIyI24D/zsyXy41MkrpXjx7xeqglEa/L4HRJ0oB2KXArcHDx/AjgV8B7S4tIkmrU12PE66XbRLxeg9MlSQPaKpn5rarnp3Q3DZckDRStjZGH17agTz0Gp0uSBrQ/R8RhwGXF80OAa0uMR5Jq1tIgmbgra0qSOvJp4PPARcXzZuDNiPg8kJm5YmmRSVI3BtMYcUnSO0xmrlB2DJLUWy2NkYfXtsS9JOmdJSI+2e55c0R8o6x4JKknWjN7vJXBRFyS1JH3RMR1EbFGRGwO3AXYSy6pIbS0Zo+3Mjg0RZK0mMw8PCIOBR4G5gIfzczbSw5LkmrSKGPE7RGXJC0mIjYATgSuAJ4GPh4RQ0sNSpJq1JI938pgj7gkqSNXA8dl5k1RWTb581QWeNu0Xg1++Zj/q1fVahAXbNl+4W69U530wAVlh9AvTMQlSR3ZLjNfh8pchcDpEXFVyTFJUk0cmiJJajgR8SWAzHw9Ij7c7vDRJYQkST3W2po93spgIi5JqnZY1eOvtDu2b38GIkm95RhxSVIjik4ed/RckgakRhmaYiIuSaqWnTzu6LkkDUgtJuKSpAb0roh4nUrv93LFY4rny5YXliTVrqwx3z1lIi5JWigzm8uOQZKWVFljvnvKRFySJEmDimPEJUmSpBI4RlySJEkqQYtjxCVJkqT+ZyIuSZIklcBEXJIkSSqBibgkSZJUAhNxSZIkqQQm4pIkSVIJTMQlSZKkEpiIS5IkSSUwEZckSZJKYCIuSZIklWCBibgkSZLU/+wRlyRJkkpgIi5JkiSVoCVNxCVJkqR+Z4+4JEmSVAITcUmSJKkEJuKSJElSCVpaW8sOoSYm4pIkSRpU7BGXJEmSSmAiLkmSJJXAlTUlSZKkEjRKj3hT2QFIkiRJfamlNXu8dSci9o2IxyJiWkSc1MHxIyLioWK7IyLe1V2d9ohLkiRpUOnrHvGIaAbOAfYCZgD3RsRVmTmlqthTwO6ZOTsi9gMmAdt3VW9NPeL1+AQgSZIk1UMdesS3A6Zl5pOZ+U/gUuCg6gKZeUdmzi6e3gWM6a7SbnvE6/UJQJIkSaqH3vSIR8REYGLVrkmZOal4PBqYXnVsBl3nup8Eft9dm7UMTVn4CaAIsu0TwMJEPDPvqCpf0ycASZIkqR6yF4l4kXRP6uRwdHRKhwUj9qCSiO/SXZu1JOJ1+QQgSZIk1UNr38+aMgNYq+r5GOD59oUiYgvgp8B+mflyd5XWkoj32SeA6i7/sWPH1tC0JEmS1DOZfZ6I3wtsEBHrAM8BhwGHVxeIiLHAb4CPZ+bjtVRaSyLeZ58Aqrv8J0yY0BgTPEqSJKmh9GZoSpf1ZS6IiOOBG4Bm4ILMnBwRxxbHzwO+DqwKnBsRAAsyc0JX9daSiNflE4AkSar4+PnfZfMD9uSNmS/zrc33KTsclWTvnTbjjC8eTlNT8LMr/8L3fnZd2SE1rDoMTSEzrwOua7fvvKrHxwDH9KTObqcvzMwFQNsngKnAZW2fANo+BbDoJ4AHI+K+ngQhSdI72Z0X/pr/3feossNQiZqagjNP+hjvP/4HvOvgr3Hovtuz8bprlh1Ww8rWnm9lqGlBn3p8ApAkDS4RsRGVG/zvzsw5Vfv3zczry4ts4Jv2l3tYdW0nHHsn23azdfn79Jk89dxLAFx2w928/91bMvXJxUYDqwZ1GCNeFy5xL0laYhFxAvA74LPAIxFRvdDFd8qJSmoco0etzIx/vLLw+XP/mM2aI4eXGFFja23NHm9lcIl7SVJf+BSwTWbOiYhxwK8jYlxmnknHs28Bi86mtSursAkr9Ees0oATHfyaZMeT1KkGfX2zZr2YiEuS+kJz23CUzHw6It5NJRlfmy4S8erZtI6NcY3xl1OqgxkzZzNmtVUWPh+92nBeeOnVEiNqbI2SiDs0RZLUF16MiC3bnhRJ+QHACGDz0qKSGsR9k59i/bGrMW7NESw1pJmP7LM919z8YNlhNazWzB5vZbBHXJLUF44EFlTvKGbdOjIiflxOSI3jkxefxfh378DyI4bzP9Pv5Opv/IA7Lris7LDUj1paWvncaf/Hted+nqamJn7+u9uY4o2avdYoPeIm4pKkJZaZM7o4dnt/xtKIzj/8hLJD0ABw/W0Pc/1tD5cdxqBgIi5JkiSVoKxZUHrKRFySJEmDSqPMI24iLkmSpEGlrJUye8pEXJIkSYNKowxNcfpCSZIkqQT2iEuSJGlQcdYUSZIkqQQm4pIkSVIJylops6dMxCVJkjSo2CMuSZIklcBEXJIkSSpBo0xfaCIuSZKkQcWVNSVJkqQSODRFkiRJKoFDUyRJkqQSZGtL2SHUxERckiRJg4qJuCRJklQCE3FJkiSpBNliIi5JkiT1O3vEJUmSpBKYiEuSJEklMBGXJEmSSmAiLkmSJJXARFySJEkqQauJuCRJktT/7BGXJEmSSmAiLkmSJJXABX0kSZKkEtgjLkmSJJXARFySJEkqgYm4JEmSVIJsbS07hJqYiEuSJGlQsUdckiRJKoGJuCRJklQCV9aUJEmSStAo84g31VIoIvaNiMciYlpEnNTB8YiIs4rjD0XE1n0fqiRJktS9bG3p8dadeuTD3faIR0QzcA6wFzADuDcirsrMKVXF9gM2KLbtgR8V/0qSJEn9qq/HiNcrH66lR3w7YFpmPpmZ/wQuBQ5qV+Yg4BdZcRewckSsUUPdkiRJUp+qQ494XfLhWhLx0cD0quczin09LSNJkiTVXR0S8brkw7XcrBkd7MtelCEiJgITi6dzIuKxGtrvzAhgVlcF4rQlqL0PYngHtD8QYii7/YEQQ9ntD4QYym6/L2JYuy+DaUTn5dMd/S15R4mIiZk5qew4ynJe2QEMEO/090Ff+OcDF/T4/5N2eSrApKqfQ5/lw9VqScRnAGtVPR8DPN+LMhQX0ydvrIi4LzMn9EVdjRpD2e0PhBjKbn8gxFB2+wMhhrLbHygxaFCYSB/9nVRD831Qgm7y1D7Lh6vVMjTlXmCDiFgnIpYGDgOualfmKuDI4m7RHYDXMvOFGuqWJEmSBrq65MPd9ohn5oKIOB64AWgGLsjMyRFxbHH8POA6YH9gGjAXOLpn1yZJkiQNTPXKh2ta0Cczrysqr953XtXjBI6r7VL6zED4yqbsGMpuH8qPoez2ofwYym4fyo+h7PZhYMSgxuf7SOD7YECqRz4clXMkSZIk9aeaVtaUJEmS1LcaMhHvbonRfmj/goiYGRGP9HfbRftrRcSfI2JqREyOiBP7uf1lI+KeiPhb0f7J/dl+u1iaI+KBiLimhLafjoiHI+LBiLivv9svYlg5In4dEY8W74cd+7HtDYtrb9tej4jP9Vf7VXH8R/E+fCQiLomIZfu5/ROLtieXcf0aWCLihOJ3cXZXf58i4hMRcXYnx66LiJXrF6WWVESMKysH6Knib9WIsuNQx2oaIz6Q1LjEaL1dCJwN/KIf26y2APhCZv41IlYA7o+IP/Tja/A2sGdmzomIpYDbIuL3xSpS/e1EYCqwYgltA+yRmV3OZ19nZwLXZ+YhxV3cQ/ur4cx8DNgSFv5ePgf8tr/aL9odDZwAbJKZ8yLiMip3sl/YT+1vBnyKyopr/wSuj4hrM/OJ/mhfA9K/A/tl5lO9rSAz9+/DeCQNYI3YI17LEqN1lZm3Aq/0Z5vt2n8hM/9aPH6DSiLabyuZFku3zimeLlVs/X6zQUSMAd4H/LS/2x4IImJFYDfgfIDM/GdmvlpSOO8B/p6Zz5TQ9hBguYgYQuWDSJdztvaxjYG7MnNuZi4AbgE+2I/tawCJiPOAdYGrim9qzi72f7j41uRvEXFr1SlrRsT1EfFERHy3qp6nI2JE0es6NSJ+UnzjcmNELFeU2TYiHoqIOyPie43SOzvINLf/2UTEpyLi3uJnfUVEDAWIiAsj4kfFt9lPRsTuxbfrUyPiwrYKI2JORJwWEfdHxB8jYruIuLk458CizLIR8bPiG9kHImKPYn9zRHy/2P9QRHy2Otgivusj4lP9+BqpG42YiPd4+dDBLCLGAVsBd/dzu80R8SAwE/hDZvZr+4UfAl8CWktoGyofPm4s/sOc2G3pvrcu8BLws+I/459GxLAS4oBKL/Ql/d1oZj4HfB94FniBypytN/ZjCI8Au0XEqsUf3P1ZdDEHvYNk5rFUPgjuAcyuOvR1YJ/MfBdwYNX+LYFDgc2BQyOio/fOBsA5mbkp8CpwcLH/Z8Cxmbkj0O3a3KqLjn42v8nMbYuf9VTgk1XlhwN7Av8BXA38ANgU2DwitizKDANuzsxtgDeAU6iMAPgg8N9FmeMAMnNz4KPAz4sheROBdYCtMnML4JdVbS9ftHlxZv6k714CLalGTMR7vHzoYBURywNXAJ/LzNf7s+3MbMnMLamsGrVd8RV9v4mIA4CZmXl/f7bbzs6ZuTWwH3BcROzWz+0PAbYGfpSZWwFvAmXcM7E0leTi8hLaHk7lG7F1gDWBYRHxsf5qPzOnAqcBfwCuB/5GZeiYVO124MKiJ7K5av9NmflaZr4FTAHW7uDcpzLzweLx/cC4Yvz4Cpl5R7H/4noFri4t9rMBNouIv0TEw8ARVBLtNlcX09s9DPwjMx/OzFZgcnEuFEPciscPA7dk5vzicVuZXYCLADLzUeAZYDzwXuC84ts5MrP6m/vfAT/LzLKG1KoTjZiI93j50MGoGJt9BfDLzPxNWXEUQyFuBvbt56Z3Bg6MiKepDE/aMyL+rz8DyMzni39nUhkbvV1/tk/ld2FG1bcRv6aSmPe3/YC/ZuY/Smj7vVT+GL5U/LH6DbBTfwaQmedn5taZuRuVIWuOD9ciip7yr1H52/VgRKxaHHq7qlgLHd+31VGZjjqk1P86+tlcCBxf9FafDCzbQfnWdue28q+f/fz817zSC8sVCXtbmc5+/kHnHZO3A/tFhO+dAaYRE/Falhgd1IpfpPOBqZl5Rgntj2y7o78Yr/he4NH+jCEzv5KZYzJzHJX3wJ8ys996QiNiWHGjLMVwkL2pDFPoN5n5IjA9IjYsdr2HSq9af/soJQxLKTwL7BARQ4vfi/dQ+Tq430TEqOLfscCHKO+10AAVEetl5t2Z+XVgFks4fCkzZwNvRGUJbaj8H6iBYQXghaKz7Ig6tXFrW90RMR4YCzwG3AgcW9wvQ0SsUnXO14GXgXPrFJN6qeES8eIrl7YlRqcCl2Xm5P6MISIuAe4ENoyIGRHxye7O6WM7Ax+n0gvcNnVcf95lvwbw54h4iMoHoz9kZr9PH1iy1ajMFvM34B7g2sy8vptz6uGzwC+Ln8WWwHf6s/FiXPReVHqi+13xbcCvgb9S+eq2if5fke6KiJhCZfzlcUWSJFX7XnED3SNUkqi/9UGdnwQmRcSdVHpCX+uDOrXk/ovKPVt/oH4dVOdSuVH0YeBXwCcy820qExc8CzxU/G06vN15nwOWrb4xWOVzZU1JkhpMRCzfNntVVOYrXyMz+3VNCUlLruHmEZckSbwvIr5C5e/4M8Anyg1HUm/YIy5JkiSVoOHGiEuSJEmDgYm4JEmSVAITcUmSJKkEJuIaFCKipZjG8ZGIuLyYVq8n5/80IjbpQflPRMTZPY9UkiSpwkRcg8W8zNwyMzejskTwsbWeGBHNmXlMZpaxGI4k1V1EzCk7Bug+jogYHxHXRcS0iJgaEZdFxGq96fwo6ll5ySJeWNdGEXFnRLwdEf/ZF3VKYCKuwekvwPoAEfGxiLin6C3/cUQ0F/vnRMR/R8TdwI4RcXNETCiOfbRt8Y2IOK2t0og4OiIej4hbqCyqJEnqIxGxLHAt8KPMXD8zNwZ+BIzsTX2ZuX9mvtpH4b0CnAB8v4/qkwATcQ0yxdK++wEPR8TGwKHAzpm5JdDCv5YcHgY8kpnbZ+ZtVeevCZwG7EllpcptI+IDEbEGcDKVBHwvoOZhLJI0EEXE+yPi7oh4ICL+GBGrFftHRsQfIuKvRQfGMxExojj2XxHxaHH8krbe4YhYLyKuj4j7I+IvEbFRsX+doif53oj4VjchHQ7cmZlXt+3IzD9n5iPF0zWLNp6oXh2yi86Tp6viPjIiHoqIv0XERVXXeUUR270R0WkHS2bOzMx7gfm1v8JS90zENVgsFxEPAvdRWeL3fOA9wDbAvcWx9wDrFuVbgCs6qGdb4ObMfCkzFwC/BHYDtq/a/08qywpLUiO7DdghM7cCLgW+VOz/BvCnzNwa+C0wFqD41vBgYCvgQ8CEqromAZ/NzG2A/6SyDDvAmVR6uLcFXuwmns2A+7s4viWVzpXNgUMjYq3OOk+qT4qITYH/B+yZme8C2lYgPRP4QRHbwVSWiJf6lStrarCYV/R6LxQRAfw8M7/SQfm3MrOlg/3RRRuufiVpMBkD/Kr4xm9p4Kli/y7ABwEy8/qImF21/3eZOQ8gIq4u/l0e2Am4vPLfLgDLFP/uTCXJBbiIStLcWzdl5mtFm1OAtYFVKTpJiv1tnSdXVp23J/DrzJxVXNMrxf73AptUxbxiRKyQmW8sQYxSj9gjrsHsJuCQiBgFEBGrRMTa3ZxzN7B7RIwoxpN/FLil2P/uiFg1IpYCPlzPwCWpH/wvcHZmbg58Gli22N9Zh0Rn+5uAV4sb5tu2jauO19qJMZnKt5idebvqcQuVzsSuOk/aRCcxNAE7VsU82iRc/c1EXINWMQvK14AbI+Ih4A/AGt2c8wLwFeDPwN+Av2bm74r93wTuBP4I/LWOoUtSf1gJeK54fFTV/tuAjwBExN7A8Kr974+IZYte8PcBZObrwFMR8eHinIiIdxXn3A4cVjxuu0enMxcDO0XE+9p2RMS+EbF5F+d01nlS7SbgIxGxalHnKsX+G4Hjq9raEqmfRabftkuSNJhFRCvwfNWuM4C/Az+gkozfBWybme8uvkW8hEoCfguVcdnrZObbEfFNKsnuM8BLVIaF/CQi1qEyw8kawFLApZn538X+i6n0Xl8BfC0zl+8izo2AHwLrUbkx8iEqY7r3AyZk5vFFuWuA72fmzRFxOJUOlACuy8wvFWWeLs6ZFRFHAV+k0pP+QGZ+oriR8xxg4yK+WzOzw6lvI2J1KvcgrQi0AnOATYoPIVKvmYhLkqSFImIZoCUzF0TEjlRuttyyOLZ8Zs6JyqJptwITM9NvCKVe8mZNSZJUbSxwWUQ0UVkg7VNVxyZFZRXiZancDG8SLi0Be8QlSVK/KcZ8X9Ru99uZuX0Z8VSLiKP51/SGbW7PzOPKiEeDn4m4JEmSVAJnTZEkSZJKYCIuSZIklcBEXJIkSSqBibgkSZJUAhNxSZIkqQQm4pIkSVIJTMQlSZKkEpiIS5IkSSUwEZckSZJKYCIuSZIklWBI2QGo5+6///5RQ4YM+SmwGX6Y0sDQCjyyYMGCY7bZZpuZZQcjSVIjMBFvQEOGDPnp6quvvvHIkSNnNzU1ZdnxSK2trfHSSy9t8uKLL/4UOLDseCRJagT2pjamzUaOHPm6SbgGiqamphw5cuRrVL6lkSRJNTARb0xNJuEaaIr3pP+nSJJUI/9oqldOOeWUUeuuu+6mK6644pZf/epXV++s3FlnnbXqkUceObajY7vvvvv6s2bNaq5flAPfY489tvQGG2ywadlx1GL06NGbv/DCCw5nkySpj/hHdRAYd9K12/RlfU+f+r77uytz/vnnj/z973//xEYbbfTP3rZzyy23TOvtuXXxzZX69HXkm691+zpKkqR3LnvE1WOHH3742BkzZixz4IEHrn/yySePauvxvuCCC4ZvsMEGm2644YabTJgwYcO28i+++OJSu+666wZrr732Zscee+yYtv1tPayPPfbY0uuuu+6mhx122Nrrr7/+pjvvvPMGc+bMCYBbbrll6Pjx4zfZcsstN/r0pz89plF6j3uipaWF9td++umnj9hss8023nDDDTfZZ5991nvjjTeaAA4++OBxRxxxxNjtt99+/JgxYza/9tprl//whz88bt1119304IMPHtdW59ChQ7f6zGc+M3rTTTfdeKeddhr/5z//eeh222234ZgxYzb/5S9/uRLA3Llz45BDDhk3fvz4TTbeeONNrr766hUAFixYwMSJE8eMHz9+k/Hjx2/y7W9/e1R1vHPmzIldd911g9NPP31EP75MkiQNOibi6rGLL7742VGjRs2/5ZZbHh8+fHhL2/5TTz11jRtvvPHxxx57bMr111+/sLd7ypQpQ6+88sonp06dOvmqq64aPm3atKXa1/nss88ue8IJJ8ycNm3a5JVWWqnlF7/4xXCAY445Zp1zzjnnmQcffPDR5ubmQTkuvqNrP+KII2Y/8sgjUx977LEpG2644byzzjprYdL72muvDbnzzjsfP/XUU6cfeuihG3zxi1/8xxNPPDH50UcfXe6OO+5YDmDevHlNe+yxxxuTJ0+eOmzYsJavfe1ro//yl788fvnll0/71re+NRrgtNNOGwXw+OOPT7n44oufnDhx4ri5c+fG6aefPvKZZ55ZZvLkyVMef/zxKcccc8zLbW2//vrrTXvvvfcGhx566Ctf+MIXZvX3ayVJ0mBiIq4+M2HChDlHHHHEuNNPP33EggULFu7fZZddXl911VVbhg4dmuuvv/5bf//735dpf+7o0aPf3mmnneYBbLXVVnOffvrpZWbNmtX85ptvNu21115vAhx11FGv9NvF9KOOrv3+++9fbpttttlw/Pjxm1xxxRWrTp48edm28u973/tebWpqYuutt5676qqrzt9uu+3mNTc3M378+Hltr+1SSy2VhxxyyOsAm2666bxddtnljWWWWSa32267ec8999zSAHfcccfyRx555MtFu2+tueaa/3z44YeX/dOf/rTiscce+9JSS1U+L6222moLP2wdeOCB63/84x+fdfzxxy9MziVJUu+YiKvPXHzxxc+ecsopz0+fPn3pLbfcctMXX3yxGWDppZde2JPd3Nyc8+fPj/bnti+zYMGCyByUHeCL6ejaJ06cuM7ZZ5/97OOPPz7ly1/+8vNvv/32wt/VZZddNouyi5zb1NTEggULAmDIkCHZ1NS0cP8yyyyz8JyWlpYA6Oz1zUwiosOD22677Zzrr79+pdbW1iW9bEmS3vFMxNVnJk+evMyee+755g9/+MPnhw8fvuDJJ59ceknqGzlyZMuwYcNab7rppmEAF1100Sr/v507ZlEbjOM4/ii0xUEC0kLLiS2UROsV4yQZdPM9OLldBscsIr4FcTfQ4CIivpuIL0A4OONlOGyLvdQkdimlpQe9Uo/H9L6fPckvkOH3JE/+x0l6+na7XbpQKOyDIEjNZrMHue96vf55MpnkhBBisVg8W6/XTyuVym2z2fw4Go1e7Pd7IYQQm83mx2SbwWBwlcvlwna7feckHAAAcH8UcRyNZVl5TdPKqqqeG4bxyTCML/96Ttu2V51O53W1Wi0dDgeRzWajPx+VfL1e76pWq71rNBqaqqq3D3GNbrd7HUVRStO0cqvVemvb9iqTyRwsy/Lz+fzXUql0XiwWy47j/LIQcBznMgiC9M8/3gIAgL/3aD7//09c113puv4ofpTbbrdpRVFiIYTo9/sv1+v1k/F4fCk7F+7muu5zXdffyM4BAEASMEccJ20+nyvD4fBVFEWps7OzYDqdrmRnAgAAOAaKOE6aaZo3pmneyM4BAABwbOwRBwAAACSgiCdTHMfxbyMAAZm+P5PMNQQA4J4o4sm09H1foYzjVMRxnPJ9XxFCLGVnAQAgKdgjnkBhGF54nvfB87z3gsUUTkMshFiGYXghOwgAAEnB+EIAAABAAt6mAgAAABJQxAEAAAAJKOIAAACABBRxAAAAQAKKOAAAACDBN4kDkOsieGPTAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, axs = plt.subplots(1, 2, figsize=(12.8, 4.8))\n", "\n", "(df.groupby(\"Period\").Lagged_Choice_1.value_counts(normalize=True).unstack().plot\n", " .bar(ax=axs[0], stacked=True, rot=0, title=\"Share of Lagged Choice per Period\"))\n", "\n", "sns.heatmap((df.query(\"Period == 0\").groupby([\"Experience_Fishing\"]).Lagged_Choice_1\n", " .value_counts(normalize=\"rows\").unstack().fillna(0)), ax=axs[1], cmap=\"RdBu_r\", annot=True)\n", "\n", "axs[0].legend(loc=\"upper center\", bbox_to_anchor=(0.5, -0.15), ncol=6)\n", "\n", "plt.show()\n", "plt.close()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Observables\n", "\n", "Now, we have talked enough about experiences and lagged choices. The last section is about observables which are other characteristics such as whether the island has rich fishing grounds or not. We can also express the distribution of observables as probability mass functions or softmax functions.\n", "\n", "Note that observables are sampled first before all other characteristics of individuals. After that, experiences and lagged choices in chronological order (highest lag first) follow. Every group itself is sorted alphabetically. Keep this in mind if you want to condition a distribution on other characteristics.\n", "\n", "For observables, it is also possible to use labels instead of numbers to identify levels which might be more convenient and intuitive. The evenly distributed levels are ``\"rich\"`` and ``\"poor\"``. The keyword is ``observable_*_*`` and everything after the last ``_`` is considered the name of the level." ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [], "source": [ "params = \"\"\"\n", "category,name,value\n", "delta,delta,0.95\n", "wage_fishing,exp_fishing,0.01\n", "nonpec_fishing,rich_fishing_grounds,0.5\n", "nonpec_hammock,constant,1\n", "nonpec_hammock,not_fishing_last_period,-0.5\n", "shocks_sdcorr,sd_fishing,1\n", "shocks_sdcorr,sd_hammock,1\n", "shocks_sdcorr,corr_hammock_fishing,0\n", "initial_exp_fishing_0,probability,0.33\n", "initial_exp_fishing_1,probability,0.33\n", "initial_exp_fishing_2,probability,0.34\n", "lagged_choice_1_fishing,zero_exp_fishing,1\n", "lagged_choice_1_hammock,zero_exp_fishing,1\n", "lagged_choice_1_fishing,one_exp_fishing,-0.2877\n", "lagged_choice_1_hammock,one_exp_fishing,-1.3863\n", "lagged_choice_1_fishing,two_exp_fishing,6\n", "observable_fishing_grounds_rich,probability,0.5\n", "observable_fishing_grounds_poor,probability,0.5\n", "\"\"\"\n", "\n", "options = {\n", " \"n_periods\": 10,\n", " \"simulation_agents\": 10_000,\n", " \"covariates\": {\n", " \"constant\": \"1\",\n", " \"not_fishing_last_period\": \"lagged_choice_1 != 'fishing'\",\n", " \"zero_exp_fishing\": \"exp_fishing == 0\",\n", " \"one_exp_fishing\": \"exp_fishing == 1\",\n", " \"two_exp_fishing\": \"exp_fishing == 2\",\n", " \"rich_fishing_grounds\": \"fishing_grounds == 'rich'\",\n", " }\n", "}" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
value
categoryname
deltadelta0.9500
wage_fishingexp_fishing0.0100
nonpec_fishingrich_fishing_grounds0.5000
nonpec_hammockconstant1.0000
not_fishing_last_period-0.5000
shocks_sdcorrsd_fishing1.0000
sd_hammock1.0000
corr_hammock_fishing0.0000
initial_exp_fishing_0probability0.3300
initial_exp_fishing_1probability0.3300
initial_exp_fishing_2probability0.3400
lagged_choice_1_fishingzero_exp_fishing1.0000
lagged_choice_1_hammockzero_exp_fishing1.0000
lagged_choice_1_fishingone_exp_fishing-0.2877
lagged_choice_1_hammockone_exp_fishing-1.3863
lagged_choice_1_fishingtwo_exp_fishing6.0000
observable_fishing_grounds_richprobability0.5000
observable_fishing_grounds_poorprobability0.5000
\n", "
" ], "text/plain": [ " value\n", "category name \n", "delta delta 0.9500\n", "wage_fishing exp_fishing 0.0100\n", "nonpec_fishing rich_fishing_grounds 0.5000\n", "nonpec_hammock constant 1.0000\n", " not_fishing_last_period -0.5000\n", "shocks_sdcorr sd_fishing 1.0000\n", " sd_hammock 1.0000\n", " corr_hammock_fishing 0.0000\n", "initial_exp_fishing_0 probability 0.3300\n", "initial_exp_fishing_1 probability 0.3300\n", "initial_exp_fishing_2 probability 0.3400\n", "lagged_choice_1_fishing zero_exp_fishing 1.0000\n", "lagged_choice_1_hammock zero_exp_fishing 1.0000\n", "lagged_choice_1_fishing one_exp_fishing -0.2877\n", "lagged_choice_1_hammock one_exp_fishing -1.3863\n", "lagged_choice_1_fishing two_exp_fishing 6.0000\n", "observable_fishing_grounds_rich probability 0.5000\n", "observable_fishing_grounds_poor probability 0.5000" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "params = pd.read_csv(io.StringIO(params), index_col=[\"category\", \"name\"])\n", "params" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [], "source": [ "simulate = rp.get_simulate_func(params, options)\n", "df = simulate(params)" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAu8AAAFiCAYAAABGRTbKAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3de5wU5Z3v8e93ZgAFRFHwAoigzgwOKiiIN0w0RgPJRjfRrCjRNVkkmBg1atTj5mRN1BzdLNmEqEFFjtGImqjxfj0mq9moq+AFGW6OCIKIgncuIgO/80fXYNv00A0001Mzn/frxes1XfVU9a+b6ae+/dRTNY4IAQAAAGj9KspdAAAAAIDiEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjs2i+1Lbf9hM7edaPt/l7qmUrHdz3bYrtrM7cP23s2sG237sXxtC70vti+xPWlzagKAzbElfX1LKVSj7fm2v7wVn3+57T1LtK/1/fyWHovy7LtvUmtlKfaH8iG8o1m2T7E9Nfmwv2X7YdvDt3S/ETEuIi4rRY1Nks57TVLrB7aftn1oKZ+jFCLi1og4tpl1698X20faXpSz/hcRMaYl6gTQPtg+3fYrtlfaXmL7d7Z3KHddrUHSD69LjivLbS+y/UfbB2W3i4iuETGviH0t2libZF8l6+dzv7RExBtJrWtLsX+UD+Ededk+T9KvJf1C0i6S+kq6VtLx5ayrgDsioquknpL+W9Ldtp3bqFSjGACQZrbPl3SVpB9L2l7SIZL2kPS47Y4tWEdr7pMXJ8eV7ZR5f2ZL+pvto0v9RK38fUArQnjHBmxvL+nnkn4QEXdHxIqIWBMR90fEj7OadrR9s+2PbdfbHpq1j31s/1cyCl5v+7isdTfZvjzr8fG2X7L9ke3XbI9oqsP2jcmo/5u2Ly/mdF9ErJH0e0m7StopGVn6u+3/tP2epEttV9j+ie0Ftt9JXsf2Obv6ru3FyfOfn1XvMNvPJK/tLdtX5znQfdX2PNvLbP/SdkWy7em2/7uZ9/2m5DV2kfSwpF5ZIz69ck8N2z4kOcPwge2XbR+Zte705Pk/tv267dGF3jcA7YftbpJ+JumHEfFI0sfPl/RPygT4b2c138b2HUl/8oLtQVn7uSjpnz+2Pacp1CZ97MVJn/5uMmK9Y7KuaTrIv9h+Q9JfbD9i+6ycGl+2/c3k59/YXpgcJ6bZPiLnJTVbY84+m61rYyJjUUT8VNIkZb70NO0ze/rjV23PTOp40/YFBfr0O23/wfZHkk7P7ecTzR2Lco+l60f3bd+izKDb/cnzXeicaThJDffZfs92g+0zsvZ1afLe5D3Go7wI78jnUEnbSPpzgXbHSbpd0g6S7pN0tSTZ7iDpfkmPSdpZ0g8l3Wq7NncHtodJulmZkZ8dJH1B0vxk9e8lNUraW9IBko6VVPB0ou1Okk6XtCgiliWLD5Y0L6nnimT96ZKOkrSnpK5N9Wc5SlJ18rwX+7PTj2sl/UhSD2Xeq6MlfT9n229IGirpQGXOVny3UN1NImKFpJFKRnySf4tzXmNvSQ9KulzSjpIukHSX7Z7JgWKCpJERsZ2kwyS9VOzzA2gXDlOmn787e2FELFcmaB6Ttfh4SX9Spq+ZIuke2x2SPv0sSQclfc1X9Fn/fbakf5T0RUm9JL0v6ZqcGr4oaZ9kuymSTm5aYbtOmS8RDyaLnpc0OKuGP9neplCNeV53MXUVcrekA5O+NteNkr6XvB/7SvpLgT79eEl3KnP8u7WZ52vuWNSsiDhV0huSvp4837/naXabpEXKvA8nSvqFP39GIe8xHuVHeEc+O0laFhGNBdr9d0Q8lMyfu0VS00jHIcqE4Ssj4tOI+IukB5TVMWf5F0mTI+LxiFgXEW9GxGzbuyjT2Z2bjPy/I+k/JY3aSD3/ZPsDSQslDVGmg26yOCJ+GxGNEbFK0mhJv4qIecnB6n9JGuXPn7b8WfLcr0j6v031R8S0iHg22dd8SdcpcyDIdlVEvBcRbygz/Sjfa98S35b0UPL+r4uIxyVNlfTVZP06Sfva3jYi3oqI+hI/P4B066Hm+/m3kvVNpkXEnclZzV8pE/oPUWYgo5OkOtsdImJ+RLyWbPM9Sf+ajFavlnSppBNz+thLkz52lTKDRYNt75GsGy3p7mRbRcQfIuLdpN8dnzxv9oBQczXmKqauQhZLsjKhNtea5P3oFhHvR8QLBfb1TETck/Tjq5ppk/dYtCVs7y5puKSLIuKTiHhJmTMKp2Y1a+4YjzIjvCOfdyX1KKIzW5L180plTltWKfMtfmFErMtav0BS7zz72F3Sa3mW7yGpg6S3kmkhHygTknfeSD1/jIgdImLniPhSREzLWrcwp22vpKbs+qqUmd+fb5sFyTayXWP7AWcu7vpImesCsg90zW5bQntI+lbTe5O8P8Ml7ZaM8pwkaZwy79+DtgeU+PkBpNsyNd/P75asb7K+P0v69UWSekVEg6RzlQnA79i+3XZTX7eHpD9n9U+zlAn7efvYiPhYmVH2pgGaUcoaibZ9vu1Ztj9M9re9Pt/v5q0xz2srpq5CeksKSR/kWXeCMoMoC2w/6cI3Tsg9NhVqU6rjSS9J7yXve/a+s4/TzR3jUWaEd+TzjKRP9PmR602xWNLuTuZ5J/pKejNP24WS9mpm+WpJPZJAvkNEdIuIgZtZU+SpcY+sx32VmaLzdtay3XPWN53m/J0yFy1VR0Q3SZcoMwqjIrbd3HpzLZR0S9Z7s0NEdImIKyUpIh6NiGOUOQjPlnTDJj4/gLbtGWX62G9mL0ymgoyU9ETW4t2z1ldI6qOkT4uIKRExXJn+NPTZXPCFykzdy+6jtomI7ONAbj93m6STk8C7raS/Js95hKSLlJmP3z0idpD0oT7f7zZbY45i6irkG5JeSAZKPicino+I45UZaLpH0h+bea0qsDxbc8eTFZI6Z63bdRP2vVjSjra3y9n3prwPKBPCOzYQER9K+qmka2z/o+3OyfzGkbbzzZvL9T/KdCoXJtsdKenrysydy3WjpO/YPjq5kKi37QER8ZYyc+bH2+6WrNvLdu70lM11m6Qf2e5vu6syo+d35JxC/t/Jax8o6TuS7kiWbyfpI0nLkxHtM/Ps/8e2uyenJs/J2rZYbytzsW3uRbRN/iDp67a/YrvS9jbJxUp9bO9i+7jkILxa0nJlRpYAQNL6fv5nkn5re0TSV/dTZt74ImWmSTQZYvubyajrucr0K8/arrX9peQ6o08krdJnfc1ESVc0TYNJrscpdLeyh5T5EvBzZfrjprO32ykzuLJUUpXtn0rqlrNt3hrzPMfm1CVn9Lb9b8pce3VJnjYdnflbHtsn03c+0mfvR6E+fWOaOxa9pMzNEXa0vasyrzvb28pc07WBiFgo6WlJ/yc5fuyvzDTW5ubdoxUhvCOviPiVpPMk/USZDnOhMhcm3VPEtp8qc6HLSGVOvV4r6bSImJ2n7XPKdEb/qcxIypP6bET8NEkdJc1U5qKiO5UZSS6FycocnJ6S9LoyB54f5rR5UlKDMiNQ/xERTX9c6QJJp0j6WJkR7XzB/F5J05TpXB9U5ktK0ZL36jZJ85LTu71y1i9U5kKnS/TZ/8+PlflMV0g6X5mRlfeUmY+fe0EtgHYuuYjxEkn/oUzQ/B9l+pKjm+aaJ+5VZire+8rMif5mEk47SbpSmX5+iTKjzU2h9jfKXOT4mO2PlQnSBxeoZ7UyF4N+WZmLTps8qsxFtHOVmdrxiTacbtJcjbk2ta5etpcrMwjyvKT9JB2ZdTzIdaqk+cmUynFK7tpTqE8voLlj0S2SXlbmIuHHtOGx6P9I+knyfBfk2e/Jkvopc6z4s6R/S66fQivniGLO2AAAAAAoN0beAQAAgJQgvAMAAAApQXgHAAAAUoLwDgAAAKQE4R0AAABIibL9pawePXpEv379yvX0ANCqTJs2bVlE9Cx3HeXEcQEAMjZ2TChbeO/Xr5+mTp1arqcHgFbF9oJy11BuHBcAIGNjxwSmzQAAAAApQXgHAAAAUoLwDgAAAKQE4R0AAABICcI7AAAAkBKEdwAAACAlCO8AAABAShQM77Yn237H9oxm1tv2BNsNtqfbPrD0ZQIAAAAoZuT9JkkjNrJ+pKTq5N9YSb/b8rIAAAAA5CoY3iPiKUnvbaTJ8ZJujoxnJe1ge7dSFQgAAAAgo6oE++gtaWHW40XJsrdyG9oeq8zovPr27VuCp97KLt2+3BWk36UflruC9OP3cMvxe1hStkdI+o2kSkmTIuLKnPVHSrpX0uvJorsj4uctWuRW0O/iB8tdQpsw/8qvlbsEINVKEd6dZ1nkaxgR10u6XpKGDh2atw0AoPWyXSnpGknHKDNY87zt+yJiZk7Tv0XEP7R4gQDQxpXibjOLJO2e9biPpMUl2C8AoPUZJqkhIuZFxKeSbldm+iQAoAWUIrzfJ+m05K4zh0j6MCI2mDIDAGgTmpsqmetQ2y/bftj2wJYpDQDavoLTZmzfJulIST1sL5L0b5I6SFJETJT0kKSvSmqQtFLSd7ZWsQCAsitmquQLkvaIiOW2vyrpHmXuSLbhztJ2LRQAlFnB8B4RJxdYH5J+ULKKAACtWcGpkhHxUdbPD9m+1naPiFiWuzOuhQKATcNfWAUAbIrnJVXb7m+7o6RRykyfXM/2rrad/DxMmWPNuy1eKQC0QaW42wwAoJ2IiEbbZ0l6VJlbRU6OiHrb45L1EyWdKOlM242SVkkalZylBVAC3LZ0y6X5lqWEdwDAJomIh5S53il72cSsn6+WdHVL1wUA7QHTZgAAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASImiwrvtEbbn2G6wfXGe9dvbvt/2y7brbX+n9KUCAAAA7VvB8G67UtI1kkZKqpN0su26nGY/kDQzIgZJOlLSeNsdS1wrAAAA0K4VM/I+TFJDRMyLiE8l3S7p+Jw2IWk725bUVdJ7khpLWikAAADQzhUT3ntLWpj1eFGyLNvVkvaRtFjSK5LOiYh1JakQANCqFJpKmdXuINtrbZ/YkvUBQFtWTHh3nmWR8/grkl6S1EvSYElX2+62wY7ssban2p66dOnSTS4WAFBeRU6lbGp3laRHW7ZCAGjbignviyTtnvW4jzIj7Nm+I+nuyGiQ9LqkAbk7iojrI2JoRAzt2bPn5tYMACifYqZSStIPJd0l6Z2WLA4A2rpiwvvzkqpt908uQh0l6b6cNm9IOlqSbO8iqVbSvFIWCgBoFQpOpbTdW9I3JE1swboAoF2oKtQgIhptn6XMqc9KSZMjot72uGT9REmXSbrJ9ivKTLO5KCKWbcW6AQDlUcxUyl8rcxxYm7mPwUZ2Zo+VNFaS+vbtW5ICAaAtKxjeJSkiHpL0UM6yiVk/L5Z0bGlLAwC0QsVMpRwq6fYkuPeQ9FXbjRFxT+7OIuJ6SddL0tChQ3O/BAAAchQV3gEASKyfSinpTWWmUp6S3SAi+jf9bPsmSQ/kC+4AgE1HeAcAFK3IqZQAgK2E8A60cv0+mVLuElJvfrkLaGMKTaXMWX56S9QEAO1FMXebAQAAANAKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKVJW7gNas3ydTyl1C6s0vdwEAAABtCCPvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKiqPBue4TtObYbbF/cTJsjbb9ku972k6UtEwAAAEBVoQa2KyVdI+kYSYskPW/7voiYmdVmB0nXShoREW/Y3nlrFQwAAAC0V8WMvA+T1BAR8yLiU0m3Szo+p80pku6OiDckKSLeKW2ZAIDWotDZWNvH256enI2dant4OeoEgLaomPDeW9LCrMeLkmXZaiR1t/1ftqfZPq1UBQIAWo+ss7EjJdVJOtl2XU6zJyQNiojBkr4raVLLVgkAbVfBaTOSnGdZ5NnPEElHS9pW0jO2n42IuZ/bkT1W0lhJ6tu376ZXCwAot/VnYyXJdtPZ2PVTKSNieVb7LtrwmAEA2EzFjLwvkrR71uM+khbnafNIRKyIiGWSnpI0KHdHEXF9RAyNiKE9e/bc3JoBAOVTzNlY2f6G7dmSHlRm9B0AUALFhPfnJVXb7m+7o6RRku7LaXOvpCNsV9nuLOlgSbNKWyoAoBUo5mysIuLPETFA0j9KuqzZndljk3nxU5cuXVrCMgGgbSoY3iOiUdJZkh5VJpD/MSLqbY+zPS5pM0vSI5KmS3pO0qSImLH1ygYAlEkxZ2PXi4inJO1lu0cz6zkjCwCboJg574qIhyQ9lLNsYs7jX0r6ZelKAwC0QuvPxkp6U5mzsadkN7C9t6TXIiJsHyipo6R3W7xSAGiDigrvAABImbOxtpvOxlZKmtx0NjZZP1HSCZJOs71G0ipJJ0UEF60CQAkQ3gEAm6TQ2diIuErSVS1dFwC0B8VcsAoAAACgFSC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICWKCu+2R9ieY7vB9sUbaXeQ7bW2TyxdiQAAAACkIsK77UpJ10gaKalO0sm265ppd5WkR0tdJACg9Sg0oGN7tO3pyb+nbQ8qR50A0BYVM/I+TFJDRMyLiE8l3S7p+DztfijpLknvlLA+AEArUuSAzuuSvhgR+0u6TNL1LVslALRdxYT33pIWZj1elCxbz3ZvSd+QNHFjO7I91vZU21OXLl26qbUCAMqv4IBORDwdEe8nD5+V1KeFawSANquY8O48yyLn8a8lXRQRaze2o4i4PiKGRsTQnj17FlsjAKD1KDigk+NfJD28VSsCgHakqog2iyTtnvW4j6TFOW2GSrrdtiT1kPRV240RcU9JqgQAtBbFDOhkGtpHKRPehze7M3uspLGS1Ldv31LUBwBtWjEj789Lqrbd33ZHSaMk3ZfdICL6R0S/iOgn6U5J3ye4A0CbVMyAjmzvL2mSpOMj4t3mdsYZWQDYNAXDe0Q0SjpLmbvIzJL0x4iotz3O9ritXSAAoFUpOKBju6+kuyWdGhFzy1AjALRZxUybUUQ8JOmhnGV5L06NiNO3vCwAQGsUEY22mwZ0KiVNbhrQSdZPlPRTSTtJujaZTtkYEUPLVTMAtCVFhXcAAJoUGtCJiDGSxrR0XQDQHhT1F1YBAAAAlB/hHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClBeAcAAABSgvAOAAAApAThHQAAAEgJwjsAAACQEoR3AAAAICUI7wAAAEBKEN4BAACAlCC8AwAAAClRVHi3PcL2HNsNti/Os3607enJv6dtDyp9qQAAAED7VjC8266UdI2kkZLqJJ1suy6n2euSvhgR+0u6TNL1pS4UANA6FDGgM8D2M7ZX276gHDUCQFtVzMj7MEkNETEvIj6VdLuk47MbRMTTEfF+8vBZSX1KWyYAoDUockDnPUlnS/qPFi4PANq8YsJ7b0kLsx4vSpY1518kPbwlRQEAWq1iBnTeiYjnJa0pR4EA0JZVFdHGeZZF3ob2UcqE9+HNrB8raawk9e3bt8gSAQCtSL4BnYPLVAsAtDvFjLwvkrR71uM+khbnNrK9v6RJko6PiHfz7Sgiro+IoRExtGfPnptTLwCgvIoe0ClqZ/ZY21NtT126dOkWlAUA7UMx4f15SdW2+9vuKGmUpPuyG9juK+luSadGxNzSlwkAaCWKGtApFoM6ALBpCk6biYhG22dJelRSpaTJEVFve1yyfqKkn0raSdK1tiWpMSKGbr2yAQBlsn5AR9KbygzonFLekgCg/Shmzrsi4iFJD+Usm5j18xhJY0pbGgCgtSlmQMf2rpKmSuomaZ3tcyXVRcRHZSscANqIosI7AABNihjQWSJuGQwAW0VRf2EVAAAAQPkR3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJYoK77ZH2J5ju8H2xXnW2/aEZP102weWvlQAQGvAMQEAyqdgeLddKekaSSMl1Uk62XZdTrORkqqTf2Ml/a7EdQIAWgGOCQBQXsWMvA+T1BAR8yLiU0m3Szo+p83xkm6OjGcl7WB7txLXCgAoP44JAFBGxYT33pIWZj1elCzb1DYAgPTjmAAAZVRVRBvnWRab0Ua2xypzClWSltueU8Tzo3k9JC0rdxEb46vKXQFaSKv+XUzJ7+Ee5S6gSCU7JkgcF7aCVv1ZlFLzecSW4fdwyzV7TCgmvC+StHvW4z6SFm9GG0XE9ZKuL+I5UQTbUyNiaLnrAPhdbFdKdkyQOC6UGp9FtAb8Hm5dxUybeV5Ste3+tjtKGiXpvpw290k6LbnDwCGSPoyIt0pcKwCg/DgmAEAZFRx5j4hG22dJelRSpaTJEVFve1yyfqKkhyR9VVKDpJWSvrP1SgYAlAvHBAAoL0fknYaIFLA9NjnlDJQVv4tA68BnEa0Bv4dbF+EdAAAASImi/sIqAAAAgPIjvAPYJLYfsr3DRtbfZPvElqwJaI/4LALtUzG3ikQbYLsyItaWuw6km21L+oeIWFfuWoD2jM8i2gryyaZj5L0Vst3P9mzbv7c93fadtjvbPtr2i7ZfsT3ZdqekfXPL59v+qe3/lvStsr4opFby+zjL9rWSXpC01naPZN1pye/oy7ZvydrsC7aftj2PkT+gNPgsotzIJ60D4b31qpV0fUTsL+kjSedJuknSSRGxnzJnTc60vU2+5Vn7+SQihkfE7S1ZPNqcWkk3R8QBkhZIku2Bkv5V0pciYpCkc7La7yZpuKR/kHRlC9cKtGV8FlFu5JMyI7y3Xgsj4u/Jz3+QdLSk1yNibrLs95K+oMyHKN/yJne0RLFo8xZExLM5y74k6c6IWCZJEfFe1rp7ImJdRMyUtEtLFQm0A3wWUW7kkzIjvLdexd7D0wXWr9jSQgDl/z2ymv89XZ3TDkBp8FlEuZFPyozw3nr1tX1o8vPJkv6fpH62906WnSrpSUmzm1kObG1PSPon2ztJku0dy1wP0F7xWURLIp+UGeG99Zol6Z9tT5e0o6T/VOZPjP/J9iuS1kmaGBGf5FtepprRjkREvaQrJD1p+2VJvypzSUC7xGcRLYx8Umb8hdVWyHY/SQ9ExL5lLgUAAEAS+aS1YOQdAAAASAlG3gEAAICUYOQdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3rFV2V5r+6Wsf/1sD7U9YSPbHGn7gWbWTbJdV+Iaq20/YPs129Ns/9X2FwpvWVq2l7f0cwIAgHSpKncBaPNWRcTgnGXzJU3dnJ1FxJgtriiL7W0kPSjpgoi4L1m2r6Shkp7KaVsVEY2lfH4AAIBNwcg7Wlz2yLrtL2aNyr9oe7ukWVfbd9qebftW207a/5ftocnPy21fYftl28/a3iVZvlfy+HnbPy8woj1a0jNNwV2SImJGRNyU7OtS29fbfkzSzba3sf1/bb+S1HtU0u5021dnvcYHbB9ZoM7+tp9J6rwsa9vdbD+VvCczbB+xZe84AABoKwjv2Nq2zQrnf86z/gJJP0hG54+QtCpZfoCkcyXVSdpT0uF5tu0i6dmIGKTMKPkZyfLfSPpNRBwkaXGB+gZKeqFAmyGSjo+IUyT9QJIiYj9JJ0v6fTJ6vzEbq/N3SZ1LstqfIunR5D0ZJOmlAvsHAPwiplEAAA9jSURBVADtBOEdW9uqiBic/PtGnvV/l/Qr22dL2iFrWspzEbEoItYpE1775dn2U0lNc+OnZbU5VNKfkp+nbEqxtv+cjHbfnbX4voho+lIxXNItkhQRsyUtkFRTYLfN1Xm4pNuSn2/Jav+8pO/YvlTSfhHx8aa8BgAA0HYR3lFWEXGlpDGStpX0rO0ByarVWc3WKv/1GWvisz8R3FybQuolHZhVzzcknS5px6w2K7J+djP7adTnP0/Zo/Ebq3ODP3EcEU9J+oKkNyXdYvu0jb8EAADQXhDeUVa294qIVyLiKmUuYh1QaJsiPCvphOTnUQXaTpF0uO3jspZ13kj7p5SZJy/bNZL6SpqjzEW4g21X2N5d0rAi6vx7Vn2jmxba3kPSOxFxg6QblfXlAgAAtG+Ed5Tbuck0lZeVme/+cCn2Kek8289J2k3Sh801TKbD/IOkcbbn2X5G0k8kXd7MJtdKqrT9iqQ7JJ0eEauVCeKvS3pF0n+o8Dx6STpH0g9sPy9p+6zlR0p6yfaLynwJ+U0R+wIAAO2APzubD7QNtjsrM9c+bI+SdHJEHF/uugAAALYU93lHWzRE0tXJ7SU/kPTdMtcDAABQEoy8o12wvZ8+f0cXSVodEQeXox4AAIDNQXgHAAAAUoILVgEAAICUILwDAAAAKUF4BwAAAFKC8A4AAACkBOEdAAAASAnCOwAAAJAShHcAAAAgJQjvAAAAQEoQ3gEAAICUILwDAAAAKVFV7gKw6aZNm7ZzVVXVJEn7ii9gaB3WSZrR2Ng4ZsiQIe+UuxgApcVxp02i304pwnsKVVVVTdp111336dmz5/sVFRVR7nqAdevWeenSpXVLliyZJOm4ctcDoLQ47rQ99NvpxbfndNq3Z8+eH9GBorWoqKiInj17fqjMqByAtofjThtDv51ehPd0qqADRWuT/E7SpwBtE8edNoh+O534D8Nmufzyy3fec889B3br1m3wJZdcsmtz7SZMmLDTaaed1jffui9+8Yt7L1u2rHLrVdn6zZkzp2N1dfXActdRjN69e+/31ltvMdUOAIAy4kDcBvS7+MEhpdzf/Cu/Nq1QmxtvvLHnww8//OqAAQM+3dznefLJJxs2d9ut4tLtS/o+6tIPC76PAJBG5TjuVFZWDqmurl7V9Pjee+9tePvtt6smT56800033bQw3zYPPPDAduPHj9/lr3/96wbHm5NOOmmPCy+88O0hQ4Z8smXVf+aVV17pdPbZZ+/e0NCwTbdu3dZ27dp17aWXXrp45MiRy0v1HMXo3LnzAStXrnyxJZ8TLYeRd2yyU045pe+iRYs6HXfccXv/7Gc/27lpZH3y5Mndq6urB9bW1tYNHTq0tqn9kiVLOhxxxBHVe+yxx77jxo3r07S8aSR3zpw5Hffcc8+Bo0aN2mPvvfceePjhh1cvX77ckvTkk092rqmpqRs8ePCA733ve33SMkq9KdauXavc1z5+/Pge++677z61tbV1X/nKV/b6+OOPKyTphBNO6Dd69Oi+Bx98cE2fPn32e/DBB7t+61vf6rfnnnsOPOGEE/o17bNz584HnHnmmb0HDhy4z2GHHVbz17/+tfOwYcNq+/Tps9+tt966vSStXLnSJ554Yr+ampq6ffbZp+7+++/fTpIaGxs1duzYPjU1NXU1NTV1V1xxxc7Z9S5fvtxHHHFE9fjx43u04NsEoJ3r1KnTutmzZ89s+ldbW/vpF77whZXNBfdC7rjjjgWlDO4rV67017/+9eoxY8YsXbhw4Yz6+vpZV1999Ruvvvpqp9y2a9asKdXToh0ivGOTTZky5Y2dd955zZNPPjm3e/fua5uWX3nllbs99thjc+fMmTPzkUceWT/KMXPmzM733HPPvFmzZtXfd9993RsaGjrk7vONN97Y5uyzz36noaGhfvvtt1978803d5ekMWPG9L/mmmsWvPTSS7MrKyvb5HzLfK999OjR78+YMWPWnDlzZtbW1q6aMGHC+qD84YcfVj3zzDNzr7zyyoUnnXRS9Y9//OO3X3311frZs2dv+/TTT28rSatWrao46qijPq6vr5/VpUuXtT/5yU96/+1vf5v7pz/9qeGyyy7rLUlXXXXVzpI0d+7cmVOmTJk3duzYfitXrvT48eN7LliwoFN9ff3MuXPnzhwzZsy7Tc/90UcfVRx77LHVJ5100nvnn3/+spZ+rwAg2wMPPLDdUUcdtbckPfjgg10HDBhQN2DAgLp99tmn7v3336+QpBUrVlSOGDFiz/79+w887rjj+q9bt06SNGzYsNqnnnqqs5QZ8PjhD3/Yu7a2tm7QoEEDFi5cWCVJ9fX1nQYNGjRg33333efcc8/t1blz5wOaq+W6667b6cADD1w+evToD5uWHXTQQZ+cffbZ70rSeeed1+vkk0/e4/DDD6/+5je/2b+5AZTc6aZHHXXU3g888MB2G6tz9uzZHQcPHjxg33333eecc87p1bTtggULOgwdOrR2wIABddXV1QMfeeSRriV661FGhHeUzNChQ5ePHj263/jx43s0NjauXz58+PCPdtppp7WdO3eOvffe+5PXXnttg1GI3r17rz7ssMNWSdIBBxywcv78+Z2WLVtWuWLFiopjjjlmhST98z//83st9mJaUL7XPm3atG2HDBlSW1NTU3fXXXftVF9fv01T+6997WsfVFRU6MADD1y50047rRk2bNiqyspK1dTUrGp6bzt06BAnnnjiR5I0cODAVcOHD/+4U6dOMWzYsFVvvvlmR0l6+umnu5522mnvJs/7Sa9evT595ZVXtvnLX/7Sbdy4cUs7dMh8x9pll13Wf0E77rjj9j711FOXnXXWWesDPQC0hNWrV1c0hfNjjjlmr9z148eP33XChAkLZs+ePfPZZ5+d3bVr13WSNGvWrG2vueaahQ0NDfVvvPFGp8cff3yDALtq1aqKQw89dPmcOXNmHnrooct/+9vf9pSks846a/fvf//778yYMWNWr169NjpcXl9fv80BBxywcmNtpk+f3vnRRx9tuP/++19vbgBlY9s3V+f3v//9vmPGjFk6Y8aMWbvuuuv6OidPnrzj0Ucf/eHs2bNnzpo1q/7ggw/eaH1IB8I7SmbKlClvXH755YsXLlzYcfDgwQOXLFlSKUkdO3ZcP2JeWVkZa9as2aBzym3T2NjoiDY50L6BfK997Nix/a+++uo35s6dO/Oiiy5avHr16vWf1W222SaStp/btqKiQo2NjZakqqqqqKioWL+8U6dO67dZu3atJam59zciZDvvyoMOOmj5I488sn3TyBUAtJTsaTOPP/74a7nrDznkkOUXXHDB7pdffvnOy5Ytq2wagNhvv/1W7LXXXmsqKys1cODAla+99lrH3G07dOgQo0aN+lCShgwZsmLBggUdJenFF1/s+t3vfvc9Sco+C1mMY445Zq/q6uqBxx577PovGiNGjPiga9euITU/gLKxfTZX5wsvvND1jDPOeE+Svve9762v85BDDllx22239TjvvPN6Pffcc9t2796dzrsNILyjZOrr6zt96UtfWvHrX/96cffu3RvnzZu3QQe5KXr27Lm2S5cu65544okuknTLLbfsWJpKW7+VK1dW9O3bd83q1at9++23b5XXPXz48OV/+MMfdpSk6dOnd3rrrbc67r///p98+ctf/mjixIk9m+Zkvv322+vvCPTLX/5y8Y477th46qmn5r2DEACUyy9+8YslkyZNWrBq1aqKww47bJ8XX3xxG0nrBy+kzABG0yBHtuwBj6qqqrxtChk4cOAnL774Yuemx48//vhrN9544+sffPDB+puDdOnSZX14bm4ApaqqKrIHSLIHbzZWZ75beY4cOXL5U089Nad3796fnn766f2vvvrqnTb1daH1IbyjZH70ox/1qampqauurh54yCGHfHzIIYesKrzVxl133XXzzzzzzD0GDx48ICK03XbbrS28VfpdfPHFi4cNG7bPEUccUVNdXV2yC6qyXXjhhe+sXbvWNTU1dSeddNJe11133fxtt902fvSjHy3t06fPpwMGDBhYW1tbd+ONN37uy8ONN964cPXq1RXZFx8DQLnV19d3GjZs2KorrrhiyX777bdixowZGx3FLsbgwYOX33TTTd2lzBSUjbU944wz3p06dWrXppsCSNKKFSuazVnNDaDstdden9bX13deu3atGhoaOkyfPr1LoToPPPDA5TfccMOOknTDDTesD+hz587t2Lt37zXnn3/+sm9/+9vLXnjhhc7N7wVpwa0i24BibrFVam+++eYrkpRciPOuJD322GMbnMbMXi9J2bfratrHbrvtpldffbW+afnPf/7zt5t+HjJkyKq5c+fOlKRLLrlk10GDBq0o/atJlOHWjrW1tZ8299ovuuiipbnt77rrrvnNbZu9LvsWYb/61a8WZ++jaV3nzp0je5smHTp00KRJkxZJWpS9vOn/S5LuvPPODbYD0H6U47hTyL//+7/v/PTTT3erqKiImpqaVSeeeOKHTzzxxBZdoPnb3/524ejRo/tPmDBh12OPPfaDrl27NjuA1LVr17j33nsbzj333D4XXXRR3x49eqzp0qXL2ksuuWRxvvYXXnjhO6eeeuoeNTU1dZWVlWoaQDnmmGOWX3PNNatra2sH1tbWrqqrqys4T/3aa699Y9SoUXtee+21uxx33HHvNy1/9NFHt5swYcKuVVVV0blz57W33nrr65v3TqA1aTfzituSl19+ef6gQYPaxZ0+brjhhu7jx4/fbe3ate7du/fqKVOmzO/Vq1dj4S1RDi+//HKPQYMG9St3HQBKqz0dd7J9/PHHFV26dFlXUVGh66+/vvsdd9yx4xNPPLHBQFWa0W+nDyPvaNXOOOOM988444z3C7cEAKC0/v73v3c+55xz+kaEunXrtvamm26aX+6aAMI7AABAHiNGjFg+Z86cmdnLnnvuuW1PO+20/tnLOnbsuG769OmzW7Y6tFeE93Rat27dOue7shwol3Xr1lkStyED0KYNGzZs1ezZs2cWbglsHdxtJp1mLF26dPskLAFlt27dOi9dunR7STPKXQuArWIdx5y2h0GXdGLkPYUaGxvHLFmyZNKSJUv2FV/A0DqskzSjsbFxTLkLAbBVzFi6dGldz549P+Ssb9vAoEt6cbcZAACwUdOmTdu5qqpqkiQGjdqO9YMuQ4YMeafcxaB4hHcAAAAgJfj2DAAAAKQE4R0AAABICcI7AAAAkBKEdwAAACAlCO8AAABASvx/pDcbFXG/VNwAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, axs = plt.subplots(1, 2, figsize=(12.8, 4.8))\n", "\n", "(df.groupby(\"Fishing_Grounds\").Choice.value_counts(normalize=True).unstack().plot\n", " .bar(ax=axs[0], stacked=True, rot=0, title=\"Choice Probabilities\"))\n", "\n", "(df.Fishing_Grounds.value_counts(normalize=True)\n", " .fillna(0).plot.bar(stacked=True, ax=axs[1], rot=0, title=\"Observable Distribution\"))\n", "\n", "axs[0].legend(loc=\"upper center\", bbox_to_anchor=(0.5, -0.15), ncol=2)\n", "axs[1].legend(loc=\"upper center\", bbox_to_anchor=(0.5, -0.15), ncol=2)\n", "\n", "plt.show()\n", "plt.close()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can see in the figure on the right-hand-side that poor and rich fishing grounds are evenly distributed. The figure on the left-hand-side shows that rich fishing grounds lead to higher engagement in fishing.\n", "\n", "## Conclusion\n", "\n", "We showed ...\n", "\n", "- how to express different distributions, probability mass functions or softmax functions, in the parameters of a respy model.\n", "- that you can use numbers and labels for discrete levels of observables.\n", "\n", "Happy modeling!" ] } ], "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.6" } }, "nbformat": 4, "nbformat_minor": 4 }