{ "cells": [ { "cell_type": "markdown", "id": "140849ff", "metadata": {}, "source": [ "# Simulation of Most Permissive Boolean networks\n", "\n", "The module `mpbn.simulation` implements random walks into the Most Permissive (MP) dynamics of a given Boolean networks, for estimating the stationnary distribution.\n", "\n", "The method is described in the paper [Variable-Depth Simulation of Most Permissive Boolean Networks](https://hal.archives-ouvertes.fr/hal-03704761/file/cmsb22.pdf) by Roncalli et al. (2022).\n", "The random walk is performed from a given initial configurations and transition probabilies are computed according to two parameters:\n", "- the *permissivity depth* of the transition, where depth 1 corresponds to asynchronous transitions, and further depth characterize the transitions that require permissive interpretation of transitions.\n", "- the numer of components that change of value in one transition.\n", "The parametrization allow scaling the probability of the transition as a function of these two quantities.\n", "\n", "This notebook demonstrates how to use to employ the `mpbn.simulation` module. Alternatively, the simulations can be performed using the command line utility `mpbn-sim`. See `mpbn-sim --help` for usage, and https://github.com/bnediction/mpbn/tree/master/examples for examples of configuration files (e.g., [simulation_bladder.json](https://github.com/bnediction/mpbn/blob/master/examples/simulation_bladder.json) and [simulation_i3ffl2.json](https://github.com/bnediction/mpbn/blob/master/examples/simulation_i3ffl2.json) files).\n", "\n", "Simulations are performed with `estimate_reachable_attractors_probabilities` and `parallel_estimate_reachable_attractors_probabilities` functions. These functions require\n", "\n", "- a Boolean network `f`, being either a `mpbn.MPBooleanNetwork` object or any object supported by its constructor (e.g., `colomini.minibn.BooleanNetwork` object, Python dictionnary, filename in BoolNet format, GINsim or bioLQM objects)\n", "- an initial configuration `x`, being a dictionnary associating components to a Boolean value\n", "- a list `A` of attractors reachable from `x`, as obtained with `f.attractors(reachable_from=init)` provided that `f` is a `mpbn.MPBooleanNetwork` object\n", "- a number `nb_sims` of trajectories to simulate\n", "- a function `depth()` which is called at each simulation step to determine the maximum permissivity depth of transitions to consider\n", "- a vector `W` of `n` reals (where `n` is the number of components of `f`) for determining the rate of a transition in function of the number of components that change of value.\n", "\n", "The simulation returns a dictionnary associating the index of the attrarctor in `A` with its estimated stationnary probability.\n", "\n", "The function for the `depth` parameter can be generated using:\n", "- `constant_maximum_depth(f)`: consider all MP transitions for the BN `f`\n", "- `constant_unitary_depth(f)`: consider only (general) asynchronous transitions of `f`\n", "- `poly_depth(f, power=1.2)`: random depth with polynomially decreasing probability: if $n$ is the number of components of `f`: it will draw a maximum depth $d$ with probability proportional to $(n-d-1)^{power}$\n", "- `reciprocal_depth(f)`: draw a maximum depth $d$ with probability proportional to $1/d$\n", "- `nexponential_depth(f, base=2)`: draw a maximum depth $d$ with probability proportional to `base`$^{-d+1}$\n", "\n", "The function for the `rate` parameter can be generated using:\n", "- `uniform_rates(f)`: constant uniform rate\n", "- `fully_asynchronous_rates(f)`: assigns a rate of 1 for transition changing exactly one components, and 0 otherwise\n", "- `reciprocal_rates(f)`: assigns a rate of $1/k$ for a transition changing the value of $k$ components\n", "- `nexponential_rates(f, base=2)`: assigns a rate of `base`$^{-k-1}$ for a transition changing the value of $k$ components" ] }, { "cell_type": "code", "execution_count": 1, "id": "1c317823", "metadata": {}, "outputs": [ ], "source": [ "import mpbn\n", "import mpbn.simulation as mpsim\n", "from colomoto_jupyter import tabulate # for pretty display" ] }, { "cell_type": "markdown", "id": "cad918ee", "metadata": {}, "source": [ "## Simple model\n", "\n", "Let us consider the following simple Boolean network we define with `mpbn.MPBooleanNetwork`:" ] }, { "cell_type": "code", "execution_count": 2, "id": "99b12217", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "a <- 1\n", "b <- a\n", "c <- c|(!a&b)" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f = mpbn.MPBooleanNetwork({\n", " \"a\": 1,\n", " \"b\": \"a\",\n", " \"c\": \"(!a & b)|c\"\n", "})\n", "f" ] }, { "cell_type": "markdown", "id": "55d1fd65", "metadata": {}, "source": [ "Its full most permissive dynamics from the configuration 000 is the following:" ] }, { "cell_type": "code", "execution_count": 3, "id": "4e8b93f8", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'a': 0, 'b': 0, 'c': 0}" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x = f.zero()\n", "x" ] }, { "cell_type": "code", "execution_count": 4, "id": "2ab7732a", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "# computing graph layout...\n" ] }, { "data": { "image/svg+xml": [ "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "000\n", "\n", "000\n", "\n", "\n", "\n", "100\n", "\n", "100\n", "\n", "\n", "\n", "000->100\n", "\n", "\n", "\n", "\n", "\n", "101\n", "\n", "101\n", "\n", "\n", "\n", "000->101\n", "\n", "\n", "\n", "\n", "\n", "110\n", "\n", "110\n", "\n", "\n", "\n", "000->110\n", "\n", "\n", "\n", "\n", "\n", "111\n", "\n", "111\n", "\n", "\n", "\n", "000->111\n", "\n", "\n", "\n", "\n", "\n", "100->110\n", "\n", "\n", "\n", "\n", "\n", "101->111\n", "\n", "\n", "\n", "\n", "\n" ], "text/plain": [ "" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f.dynamics(init=x)" ] }, { "cell_type": "markdown", "id": "7510debe", "metadata": {}, "source": [ "its fully asynchronous dynamics is the following:" ] }, { "cell_type": "code", "execution_count": 5, "id": "ff9b0ba7", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "# computing graph layout...\n" ] }, { "data": { "image/svg+xml": [ "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "000\n", "\n", "000\n", "\n", "\n", "\n", "100\n", "\n", "100\n", "\n", "\n", "\n", "000->100\n", "\n", "\n", "\n", "\n", "\n", "110\n", "\n", "110\n", "\n", "\n", "\n", "100->110\n", "\n", "\n", "\n", "\n", "\n" ], "text/plain": [ "" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f.dynamics(\"asynchronous\", init=x)" ] }, { "cell_type": "markdown", "id": "6ce58aea", "metadata": {}, "source": [ "Thus, we can deduce that the MP transition `000 -> 101` has requires a permissive depth strictly greater than 1 (see the method paper for details).\n", "\n", "The simulation functions requires as input the list of reachable attractors.\n", "These can be obtained as follows (two fixed point in this example):" ] }, { "cell_type": "code", "execution_count": 6, "id": "f2ea81cf", "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", "
abc
0110
1111
\n", "
" ], "text/plain": [ " a b c\n", "0 1 1 0\n", "1 1 1 1" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "A = list(f.attractors(reachable_from=x))\n", "A.sort(key=lambda a: [a[i] for i in f]) # stable ordering, for notebook reproducibility only\n", "tabulate(A)" ] }, { "cell_type": "markdown", "id": "45df790c", "metadata": {}, "source": [ "We first estimate their propability in the stationnary distribution from 1,000 simulations of a unform random walk in the full MP dynamics from `x`:" ] }, { "cell_type": "code", "execution_count": 7, "id": "90dd7ec0", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "100%|███████████████████████████████████████████████████████████████| 1000/1000 [00:00<00:00, 7634.96it/s]\n" ] }, { "data": { "text/plain": [ "{0: 52.1, 1: 47.9}" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nb_sims = int(1e3)\n", "urw = mpsim.estimate_reachable_attractors_probabilities(f, x, A, nb_sims,\n", " depth = mpsim.constant_maximum_depth(f), # full MP dynamics\n", " W = mpsim.uniform_rates(f)) # uniform random walk\n", "urw" ] }, { "cell_type": "markdown", "id": "be89cba6", "metadata": {}, "source": [ "Let us compare with exponentially decreasing depth and rate..." ] }, { "cell_type": "code", "execution_count": 8, "id": "54ad9380", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "100%|███████████████████████████████████████████████████████████████| 1000/1000 [00:00<00:00, 4209.22it/s]\n" ] }, { "data": { "text/plain": [ "{0: 94.9, 1: 5.1}" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "edw = mpsim.estimate_reachable_attractors_probabilities(f, x, A, nb_sims,\n", " depth = mpsim.nexponential_depth(f),\n", " W = mpsim.nexponential_rates(f))\n", "edw" ] }, { "cell_type": "markdown", "id": "36cd1ecb", "metadata": {}, "source": [ "... and considering only fully-asynchronous transitions:" ] }, { "cell_type": "code", "execution_count": 9, "id": "69d21078", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "100%|███████████████████████████████████████████████████████████████| 1000/1000 [00:00<00:00, 6010.25it/s]\n" ] }, { "data": { "text/plain": [ "{0: 100.0, 1: 0.0}" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "faw = mpsim.estimate_reachable_attractors_probabilities(f, x, A, nb_sims,\n", " depth = mpsim.constant_unitary_depth(f),\n", " W = mpsim.fully_asynchronous_rates(f))\n", "faw" ] }, { "cell_type": "code", "execution_count": 10, "id": "a096cbbb", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAs0AAAEXCAYAAABWAnlVAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABgP0lEQVR4nO3dd3wUdfrA8c/2ZDe9J7QkEHrvRQRR8EBRURQLNixnl7Oc6KlYT4HTw98himc9yx02REFRRFSkSe9NII30Xjfb5vdHJBJaAikzu/u8X6+8ILOzM0+SnZlnvvP9Pl+doigKQgghhBBCiFPSqx2AEEIIIYQQWidJsxBCCCGEEA2QpFkIIYQQQogGSNIshBBCCCFEAyRpFkIIIYQQogGSNAshhBBCCNEASZqFEEIIIYRogCTNQgghhBBCNECSZiGEEEIIIRogSXMLevfdd9HpdKSmptZb/vjjj9O+fXuMRiNhYWGqxNZYTz31FDqdDr1ez6FDh054vbKykpCQEHQ6HTfddFPd8tTUVHQ6Xd2XXq8nMjKSCRMmsHbt2lb8CYTwfcceb0899dRJ15k2bVrdOscaPXp0vWM1MDCQPn36MHfuXDweTytEL0TTLFy4kB49ehAYGIhOp2Pr1q2Nfu/xx8yPP/6ITqfjxx9/bPY4hfeTpLkFXXTRRaxdu5b4+Pi6ZYsXL+b555/nhhtu4KeffuL7779XMcLGCwoK4p133jlh+SeffILT6cRkMp30fffeey9r165l1apVvPDCC2zbto3zzjuPLVu2tHTIQvid4OBg3n333ROS3YqKCj755BNCQkJO+r7k5GTWrl3L2rVrWbhwIW3atOEvf/kLjz76aGuELcRZy8/P5/rrr6djx44sW7aMtWvX0rlzZ7XDEj5KkuYWFB0dzdChQ7FYLHXLdu7cCcB9993HiBEjGDhwYJP3U1VV1eRtNGTKlCm89957J1yM33rrLSZNmoTZbD7p+9q3b8/QoUMZMWIEt99+O++//z41NTXMnz+/xWMWwt9MmTKFtLQ0VqxYUW/5woULcbvdXHLJJSd9X2BgIEOHDmXo0KFccsklLF68mOTkZObNm4fT6WyN0IU4K/v378fpdDJ16lRGjRrF0KFDsVqtaoclfJQkzb+76aabSExMPGH50e4Jx9LpdNxzzz28//77dOvWDavVSp8+fViyZEm99Y7vnpGYmMjjjz8OQGxsbL3HQh6Ph9mzZ9O1a1csFgsxMTHccMMNZGZm1tvm6NGj6dmzJz///DPDhw/HarUybdq0usezc+bMYdasWSQmJhIYGMjo0aPrTiozZswgISGB0NBQJk2aRF5eXqN/P9OmTSMjI4Ply5fXLdu/fz+//PIL06ZNa/R2hg4dCkBaWlqj3yNEazl6vO/atYtrrrmG0NBQYmNjmTZtGqWlpfXWVRSF+fPn07dvXwIDAwkPD2fy5Mn1ujH973//Q6fTMW/evHrvnTlzJgaDoe54Onr8zp49m+eff5727dsTEBDAwIEDT0iAT6dLly4MHz6ct99+u97yt99+m8svv5zQ0NBGbcdkMjFgwACqqqrIz89v9P6FaE033XQT55xzDlB7w6jT6Rg9enTd18nWP9l1/nTef/99dDrdSbsVPvPMM5hMJrKysk67jYULFzJu3Dji4+MJDAykW7duzJgxg8rKynrrHTp0iKuvvpqEhAQsFguxsbGcf/75dd1NbrnlFiIiIk7aUDZmzBh69OhR931j8xSAvXv3cs011xAbG4vFYqF9+/bccMMN1NTUNOZX5FckaT5LS5cuZd68eTzzzDN89tlnREREMGnSpJP2+z1q0aJF3HLLLQB1j5FuvfVWAO68804eeeQRxo4dy5dffsmzzz7LsmXLGD58OAUFBfW2k52dzdSpU7n22mv5+uuvueuuu+pee/XVV1m9ejWvvvoqb775Jnv37mXixInccsst5Ofn8/bbbzN79my+//77un03RkpKCiNHjqx3MX777bdJTEzk/PPPb/R2fvvtN6C2FV4Irbriiivo3Lkzn332GTNmzOCjjz7iL3/5S711/vznPzN9+nQuuOACvvjiC+bPn8+uXbsYPnw4ubm5AFx99dXccccdPPjgg2zcuBGAH374geeee47HHnuMsWPH1tvmvHnzWLZsGXPnzuWDDz5Ar9czfvz4MxoHcMstt/DFF19QXFwMwL59+1izZk3duaexDh48iNFoJDw8/IzeJ0RreeKJJ3j11VcB+Pvf/87atWub/SnmlClTiIuLq9vPUS6XiwULFjBp0iQSEhJOu40DBw4wYcIE3nrrLZYtW8b06dP5+OOPmThxYr31JkyYwKZNm5g9ezbLly/ntddeo1+/fpSUlABw//33U1xczEcffVTvfbt372blypXcfffd9ZY3Jk/Ztm0bgwYNYt26dTzzzDN88803vPDCC9TU1OBwOM701+X7FKEoiqLceOONSocOHU5YPnPmTOX4XxOgxMbGKmVlZXXLcnJyFL1er7zwwgt1y9555x0FUA4fPnzC9vLz8+uW7dmzRwGUu+66q95+1q9frwDKY489Vrds1KhRCqCsWLGi3rqHDx9WAKVPnz6K2+2uWz537lwFUC655JJ660+fPl0BlNLS0tP8VurH+8477ygWi0UpLCxUXC6XEh8frzz11FOKoiiKzWZTbrzxxhPimTVrluJ0OhW73a5s2rRJGTRokAIoS5cuPe1+hVDD0c/77Nmz6y2/6667lICAAMXj8SiKoihr165VAOWll16qt15GRoYSGBio/PWvf61bZrfblX79+ilJSUnK7t27ldjYWGXUqFGKy+WqW+fo8ZKQkKBUV1fXLS8rK1MiIiKUCy644LRxH33/nDlzlPLyciUoKEiZN2+eoiiK8vDDDytJSUmKx+NR7r777hPOZ6NGjVJ69OihOJ1Oxel0KllZWcqMGTMUQLnyyivP4LcnROtbuXKlAiiffPJJ3bJRo0Ypo0aNOmHdk13nAWXmzJknbG/lypV1y2bOnKmYzWYlNze3btnChQsVQPnpp5/OKF6Px6M4nU7lp59+UgBl27ZtiqIoSkFBgQIoc+fOPe37R40apfTt27fesjvvvFMJCQlRysvL6/1cjclTxowZo4SFhSl5eXln9HP4K2lpPkvnnXcewcHBdd/HxsYSExNzVt0OVq5cCVCv+gTA4MGD6dat2wmPZ8PDwxkzZsxJtzVhwgT0+j/+rN26dQNqByUe6+jy9PT0Rsd55ZVXYjab+fDDD/n666/Jyck5IebjPfLII5hMJgICAhgwYADp6eksWLCACRMmNHq/QrS24/v+9u7dG7vdXtelacmSJeh0OqZOnYrL5ar7iouLo0+fPvVG3lssFj7++GMKCwvp378/iqLw3//+F4PBcMJ+L7/8cgICAuq+Dw4OZuLEifz888+43e5GxR4UFMSVV17J22+/jcvl4j//+Q8333zzCd3MjrVr1y5MJhMmk4mEhAReeuklrrvuOv797383ap9C+LI777wToN7xMG/ePHr16sW5554L1HaxPPZccOzxeujQIa699lri4uIwGAyYTCZGjRoFwJ49ewCIiIigY8eOzJkzh5dffpktW7actHrN/fffz9atW1m9ejUAZWVlvP/++9x4440EBQXVW7ehPKWqqoqffvqJq666Sp7+NpIkzWcpMjLyhGUWi4Xq6uoz3lZhYSFAvSobRyUkJNS9ftTJ1jsqIiKi3vdHB+idarndbm90nDabjSlTpvD222/z1ltvccEFF9ChQ4fTvuf+++9nw4YNbNq0iYMHD5Kdnc3tt9/e6H0KoYbjj++jg3mPHt+5ubkoikJsbGxdsnn0a926dSd0qerUqRMjR47Ebrdz3XXXnfIYjouLO+kyh8NBRUVFo+O/5ZZb2Lx5M88//zz5+fkN3tx27NiRDRs2sHHjRnbu3ElJSQkffPBBo/tAC+HLYmNjmTJlCgsWLMDtdrN9+3ZWrVrFPffcU7fOtGnT6p0HjnZbrKioYOTIkaxfv57nnnuOH3/8kQ0bNvD5558Df5xTdDodK1as4MILL2T27Nn079+f6Oho7rvvPsrLy+v2c+mll5KYmFjXXeTdd9+lsrLyhK4Z0HCeUlxcjNvtpm3bts30m/J9RrUD0IqAgICTdno//uLXEo5+sLOzs0/48GZlZREVFVVv2elajFratGnTePPNN9m+fTsffvhhg+u3bdu2WSqECKElUVFR6HQ6Vq1aVa86zlHHL3vzzTdZunQpgwcPZt68eUyZMoUhQ4ac8L6cnJyTLjObzSe0Ip3OiBEj6NKlC8888wxjx46lXbt2p13/6KBDIXxBQEDACQN3oWnX8/vvv5/333+fxYsXs2zZMsLCwrjuuuvqXn/qqafqJdFHW3h/+OEHsrKy+PHHH+tal4G6fsrH6tChA2+99RZQO9D+448/5qmnnsLhcPD6668DoNfrufvuu3nsscd46aWXmD9/Pueffz5dunQ5458pIiICg8FwQsEBcWrS0vy7xMRE8vLy6gbwADgcDr799tsW3/fRrhYffPBBveUbNmxgz549ZzTQrqUNGzaMadOmMWnSJCZNmqR2OEKo4uKLL0ZRFI4cOcLAgQNP+OrVq1fdujt27OC+++7jhhtuYNWqVfTu3ZspU6bUDdQ71ueff17v6U95eTlfffUVI0eOPGl3jtN5/PHHmThxIg8++ODZ/6BCeKHExET2799fryGssLCQNWvWnPU2BwwYwPDhw5k1axYffvghN910Ezabrd4+jz0HHE1ijzZyHX8jvWDBgtPur3Pnzjz++OP06tWLzZs313vt1ltvxWw2c91117Fv3756yfqZCAwMZNSoUXzyySet0kDoC6Sl+XdTpkzhySef5Oqrr+bhhx/Gbrfzf//3f43uR9gUXbp04fbbb+df//pX3Wj51NRUnnjiCdq1a3fCqH21Hb0TFsJfHa07fvPNN7Nx40bOPfdcbDYb2dnZ/PLLL/Tq1Ys777yTyspKrrrqKpKSkpg/fz5ms5mPP/6Y/v37c/PNN/PFF1/U267BYGDs2LE88MADeDweZs2aRVlZGU8//fQZxzh16lSmTp3aTD+xEN7j+uuvZ8GCBUydOpXbbruNwsJCZs+efcrJfRrr/vvvryttd2zVqtMZPnw44eHh3HHHHcycOROTycSHH37Itm3b6q23fft27rnnHq688kpSUlIwm8388MMPbN++nRkzZtRbNywsjBtuuIHXXnuNDh06nFCF40y8/PLLnHPOOQwZMoQZM2bQqVMncnNz+fLLL1mwYEG9PtFCkuY6SUlJLF68mMcee4zJkycTHx/PAw88QH5+/lldsM7Ua6+9RseOHXnrrbd49dVXCQ0N5U9/+hMvvPDCSfslCSHUtWDBAoYOHcqCBQuYP38+Ho+HhIQERowYweDBgwG44447SE9PZ8OGDXWtUsnJybz55ptceeWVzJ07l+nTp9dt85577sFut3PfffeRl5dHjx49WLp0KSNGjFDjRxTCK40YMYL33nuPF198kUsvvZTk5GRmzpzJ119/3aTpsS+77DIsFgvnnXceKSkpjXpPZGQkS5cu5cEHH2Tq1KnYbDYuvfRSFi5cSP/+/evWi4uLo2PHjsyfP5+MjAx0Oh3Jycm89NJL3HvvvSdsd8qUKbz22mvceeed9Qb/n6k+ffrw66+/MnPmTB599FHKy8uJi4tjzJgxp5y0zJ/pFEVR1A5CCCH8WWpqKklJScyZM4eHHnpI7XCEECfx1Vdfcckll7B06VLVK0A9+OCDvPbaa2RkZEjDWiuSlmYhhBBCiFPYvXs3aWlpPPjgg/Tt25fx48erFsu6devYv38/8+fP589//rMkzK1MkmYhhBBCiFO46667WL16Nf379+e9995TtYLVsGHDsFqtXHzxxTz33HOqxeGvpHuGEEIIIYQQDZCSc0IIIYQQQjRAkmYhhBBCCCEaIEmzEEIIIYQQDZCkWQghhBBCiAZI0iyEEEIIIUQDJGkWQgghhBCiAZI0CyGEEEII0QBJmoUQQgghhGiAJM1CCCGEEEI0QJJmIYQQQgghGiBJsxBCCCGEEA2QpFkIIYQQQogGSNIshBBCCCFEAyRpFkIIIYQQogGSNAshhBBCCNEASZqFEEIIIYRogCTNQgghhBBCNECSZiGEEEIIIRogSbMQQgghhBANkKRZCCGEEEKIBkjSLIQQQgghRAMkaRZCCCGEEKIBkjQLIYQQQgjRAEmahRBCCCGEaIBR7QDEyeWVO8kpc1FS7a77Kre7qXJ4qHJ6qHZ4cLgVzAYdFqMei1H3+5eeAJMOs1FPgFFHdJCRNmEm2oSasFkMav9YQvi1GpeHggoX+RUuCipcFFS6ya9wUe304PYoeBTq/vV4FNxK7fusZj1BFj3BFgMhAXoibEYirAYibEZig41YjNL+IYRanG6FwkoXhZUuiirdFBzz/6rfj223p/bYBjDodRj0YNTrsFn0RNqMv38Z6v4fbjVg1OtU/snE8SRpVllRpYvDhQ5SixwcLqyp/X+hg0qHp9n3FRpooE2oiTZhJhJ+/7dTlIWkSDM6nRycQjSXaqeHA3k17Muzsy+3hsOFNeRXuCizN/9xrddBm1ATSZFmkqIsJEeaSY6y0CbMhF6OayGajcutcLiwhv15R7/sHCl1Um73oDTzvvS62mt2+3AzKdEWOsda6BITQLtwOa7VpFMUpbn/1uI0DhXUsCG9io3pVezPq6G02q12SIQGGuidEECfNoH0aWulY5RZDkohGklRFPbn1bAr216XJKcXO/CofGa1GHUkRZrp0yaQQe2t9GoTKC3SQpyBzGIHW49Usy/Xzv68Gg4VOnC61T2wA006OkVb6BwTQJcYC/3aWYkOkvbP1iJJcwsrs7vZmF7FhrQqNqRXUVDhUjukBgVZ9PRKCKRPm0AGd7CSHGVROyQhNKXc7ubXtCrWp1ayIa2KYg3c/DbEYtTROyGQgR2sDGovx7UQx3N7FHZl21lzuJLVhyrIKHaqHVKDdEBKjIXhSTZGJNtIiQlQOySfJklzC8grd7JsTzlrD1WwL69G9RanpmofbuK8zsGM6RxMhwiz2uEIoYrCShe/HKzg598q2HqkGnfz97RoVZE2AyM7BvGn7iF0jZULrfBPVQ4PG9IqWXO4knWpVZp4+tsUMUFGhifbGJZko387KyaDPDVuTpI0NxOXW2H1oQqW7ipjY3qV1yfKp9Ip2sKF3YIZ2yWYMKs8EhK+zaMo/JpaxZc7SlmXWumzx3VSpJkLu4UwrmswETY5roVvUxSFrZnVLN1Vxs+/VeBQuctFS7GZ9VzQJZiLe4ZIC3QzkaS5idKLHSzdWcp3e8q94hFtczHqYWiSjct6hzGwvVXtcIRoVsVVLr7eVcaSnaVkl2m/S1VzMehhcAcbf+oWzIjkIIzSSiV8SLndzdJdZXy1o5QjpdrvetGcOsdYuKRXKGO7BsvYhiaQpPkseBSFnw5U8Pm2EnZk2dUOR3VdYy1cNzCCczrapAqH8GrbjlTzxbYSVh2swOXl3S+aKirIyOS+YVzSKxSrWS6ywnulFtbw2dYSlu8tx+7y75QnJEDPxT1Duax3KDHBJrXD8TqSNJ8Bj6Lww75yPthQTGqRQ+1wNCcx0sy1A8M5v3MwBqkvKbzIjqxq3lxTyLYj1WqHojlBFj2X9Q7lqv7hhARIrXfhPQ4X1vDmmkJWH6pUOxTNMejhwq4h3DQ0QpLnMyBJcyMoisLKAxW8s67QK0bTqi0h1MSU/mGM7x6CWR4DCQ3bn2fnrTWFrE+rUjsUzbOa9VzRN4yr+ocRLBMlCQ3LLXPy9rpClu8t99lxCM3FbNBxWZ9Qpg6KkJviRpCkuQEb0yp5Y00h+/Nq1A7F60QFGbnznCjO7xKsdihC1JNaWMPbawtZdbCy2Scl8HU2s57rB0cwuV+YzFgmNKW02s0HG4r4Ynup6vWUvY3NrOeaAeFM7hdGgEkau05FkuZTyCxx8MrKfDakSwtUUw1oF8j082JoFy7l6oS6yuxuXl9VwLI9ZdIC1URJkWamnxdDnzaBaoci/Fy108MnW0pYuKm4RWbT9SeRNgM3DI7gop6hclN8EpI0H8ftUfh0awlvry2kxs8HDDQnk0HH1QPCmTooXEbuClX8sL+cf/2Y71dVblrDuK7B3DkyinApQSlUsDmjitnf55LjR1VuWkPnGAszxsbKJEjHkaT5GIcLa5i9PJc9udIVo6UkhJq4b1Q0Q5Nsaoci/ER+hYt//pDHmsMyGKilBFn03DIskkt7h6KXCjqiFVQ5PCz4pYAvd5RKF6sWYjLouH5wBNcODJdW599J0kztxCQfbizigw3F0g+qlZzb0cYD58cSFigDD0TLUBSFL3eU8sbqQnlk20q6xVp4/E9xtAmTrlii5Wz5vXXZn2qoq0lanf/g90nz3lw7s5fncqhQSsi1tuggI0+Oj6NXgvSJFM0rt8zJ37/LlRJyKrCZ9Tx8QQyjU2QAsGhe1c7a1uXF26V1ubWZDDpuGBzBNX7e6uzXSfOibSW8+nO+309ioCaDHm4bHsXVA8LVDkX4iHWHK/n7dzmU2eXAVtPEniHcMypaxjCIZrEv187T3+SQ5Wcz+WlNlxgLT18UT1yIf9Z29suk2elWmLsyj6W7ytQORfxuRLKNGeNipf6rOGtuj8Lbawv5aGOxtEJpRMcoMzMnxNNeKueIJlixr5zZ3+fK4HyNCAs08PRF8X5ZOcfvkubCShczl2azM1umv9aahFATT02Io3NMgNqhCC9TbnfzzDc5UiJSgwJNOh4YE8PYriFqhyK8jEdReGtNIR9uLFY7FHEcox7uGx3DJb1C1Q6lVflV0rwnx84TS7IoqJSSU1plMui4f3Q0F/f0rwNRnL3DhTU8/lU2R+SxrabdOCSCm4dGqh2G8BJVDg/PLcuRqjcad2mvUO4dHe03/Zz9Jmn+ZncZL/+QJ9UxvMT1gyO4ZZhcYMXpbcmo4m9LsqmS6hheYWLPEKafF4PBTy6w4uwcKXHwt6+ySS2SAfreoF/bQJ6aEE+oH1TD8ouk+e21hfzn1yK1wxBn6OKeIfxFLrDiFNYermTm0mwcciPsVUYk23hyfJwMEBQntf1INY8vyZKBvF4mPsTI7Mva+PzMvz6fNL+5poAPNkh/KG91TkcbM8fHYzJI4iz+sHJ/Oc9/myOVb7xUz/gAXrgkgeAA32+ZEo23Kb2Kv32VhV0G/HmlCKuBly5vQ1Kk79Zz9umk+Y3VBXwkAwi83uAOVp69OF5apgQAX+8q5R8r8vD47JnLPyRGmJl9WQIxwf5ZukrUt+5wJU/KkyOvFxpo4B+T2pAS7ZuJs88mzQt+KeC/myRh9hX92wby/CUJBJokcfZnn20tYd5P+VJSzkckhJr4v8ltiQoyqh2KUNGaQxXM/DpHxhz5iGCLnpcvb0OKD1bC8skM5HVJmH3O5sxqZizOwiHP4/3W/zYV8y9JmH1KVqmTBxcdoaRaKhr5qw1plZIw+5jyGg8PfZHFoYIatUNpdj6XNM9flc//JGH2SduOVPPC8lx89OGIOI3v95ax4JcCtcMQLSCtyMHDi45QUSOJs7/ZklnF40uyJWH2QaXVbh5cdIT0Yt+qgOJTSfN76wv5eHOJ2mGIFrRyfwVvrC5UOwzRirZkVjHr+zxpYfZhB/JreHxJtjxJ8iO/5dfw2JdZMsufDyuucvPg50corHSpHUqz8Zk+zSv2lfPsshy1wzhB6jf/R9q38+otMwVHMfzZNXjcTlKXzqVoz09UF2ZgDAgmvPMwkiY+hCU09pTbrMw+QOo3r1CesYua4iN0vOwx2o6+qd46uRu/5PCSf+B2VBM3ZDIdL32k7jV7YSbbX59G/wc/xxgQ1Kw/b2t5YIz/zUTkj1ILa7jnk0wqaiSZ8gfnpQTx5Pg4dDqpluPLSqpc/Pl/GeSW+04yJU6te1wAc69og9kHBvP7xOiLXdnVvLg8V+0wTskal0Kfu979Y4G+tsySx2GnPHMX7cfdRVBCV1zVZfy26Hl2vnknAx78/JTbczurCYhsR3TfP3HwixdOeN1ZUcT+hX+jyzUvEhDVjp1v3E5Yp8FE9jgPgP2fzCTp4ge9NmEGmLsyj9hgI0MSbWqHIlpIYaWLRxZnScLsR1YeqCA6qIC7zo1WOxTRQlxuhSeXZkvC7Ed259h56Yc8Hh0Xp3YoTeb1aX9euZPHv9J2nyid3oA5JPqPr6AIAIyBwfS5611i+k3AGptMSGJfOl3xBBUZO7EXZ51yeyHte9Px0keI6X8xOsOJhcSrCzMwBAQT0/8iQtr3JqzTEKpyDwKQu+kr9EYz0X0ubJkftpV4FHjq62wO5NnVDkW0gCqHh0e/zJILqx/6eEsJK/eXqx2GaCFzf8xje5act/3Nt3vKWbjZ+8ebeXXSXOPy8MSSbIo1PvK6uiCNtU+ew/pnxrD7velUF6Sfcl13dTnodBgDQ856f4HRiXgc1ZRn7sZZWUJ5xg5s8V1wVpaQ+s0rdLriybPetpZUOxUe/TKL/ApJrHzNi8tz2Z/neyOvRePMWZFHho8NIBLw+bYSluwsUzsMoZIFvxSwPrVS7TCaxKuT5rkr89mn8QtrcIc+dL1uNr3ueIvOU57FUVbAlleuxll54h2Xx1nDoSUvEdN/YpO6TpisoXS9bhb7Pvwrm/85mdiBlxHRbSSHFs+izcjrsRdlsmnOpWx48SLyty5ryo+nuoJKNy98myMVNXzIlztK+fm3CrXDECqqcniYuTSbGhkY6DM2Z1Tx6s/5aochVORR4Jlvcry6oobX9mn+akcp3+zW/h1rZPdRx3zXhZDEfqx/7gJyfl1Eu/Om1b3icTvZ/d50UDykXPlUk/cb1XscUb3H1X1fcmA9ldn76DT5SX59bizdbngZc3AUm/85mdCOgzAHRzZ5n2rZnFnNws0lXD0gXO1QRBOlFtbIhVUAcKjQwT9/yGOGD/SD9Hd55U6e+jobt9wD+b1Kh4e/fZXFv69pT4AXTlbmfRFTWxD/1VXeeWE1WKzY4jtTnZ9Wt8zjdrL73fuxF2XS+853mn2Ansfl4MCnT5Fy1bNUF6SheFyEdRqMNTYZa3QiZWnbmnV/anhrbSG/5Wv7qYM4PYfLw7PLcqQElaizbE85S3eWqh2GaKJ/rMijzC4Zs6iVUezk32u8s3Ss1yXNiqIw5/tc7E7vvLB6XA6qcg9iDqkdHX40Ya7OT6P3Xe9hsjV/a2nat68S0e1cgtv1QPF4UDx/9AH3uF3g0Xaf8MZwuhWe+zZHHud6sdd/KeBggfc+thMt45Uf831yZjF/sWRnKb+mVakdhtCYz7eWsO1ItdphnDGvS5q/3FHKlkzv+UUfXPwiJb/9SnVhBmWp29j1zr247RXEDZ6E4nax+537qMjYSbfr/wEeN46yfBxl+XhcfyQPez94mENf/aPue4/LQUXmbioyd6O4ndSU5lKRubte6/VRldkHyN/yNYnj7wfAGpMMOh3Z6z6hcNdKqvIOEdy+V8v/IlpBaqFDJj7xUmsOVfD5NmlRFCdyuBXmrMjDI+MWvE5euZPXVslMnuJECjBreS7VTu9q6PKqPs25ZU4WeFlSVFOSw57/PICzshhTUDghHfrS7y+fEBDRBnthJoU7VwCwac6l9d7X5+73CUsZAoC9OBt0f9zfOErz2PSPy+q+z1z5Fpkr3yK042D63vtB3XJFUdj/8RN0nPQYBosVAIM5gK7XzuLAp0/jcTlIueJJLGG+02fw860lDE20MqiD1G/2FiXVbmZ9n6d2GELD9uTYWby9lEl9wtQORZyBOd/nUenwrqRItJ6sUidvrC7g/tExaofSaF41I+DDi46wIV0e84jTi7QZeO/6DgRZDGqHIhph1vJcrxjUK9RlM+v5zw0diLR5VVuP31qys5R/rJCbYXF6OuDlK9rQr61V7VAaxWu6ZyzdWSoJs2iUwko37/9apHYYohF2ZFWzTBJm0QiVDg//96N3DgD3N7llTuZLtwzRCAow24u6aXhF0lxY6ZIDUJyRz7eVklXqVDsMcRpuj8I/f8jDax51CdX99FsFaw5JDW+tm/dzPlXSLUM0UnaZy2saurwiaX53XaH0ixJnxOlWWPCL3Ghp2Vc7SjlUKNUyxJl55cd8r2mV8ke7sqtZddC7Z30Tre+zrSUUeMHsvppPmrNKndLfUZyVn36rYEeW91Ra8SflNW7eWe8dLQtCW3LLXSzcdOKMqkIbpIKROBs1LoX31mv/s6P5pPmddYVI6V1xtuavKpAptjXo/V+LKK32/vrgQh2fbCmhvEY+P1qzPrXSK2vvCm34encZGRqfYlvTSfPhwhpW7CtXOwzhxfbk2OUzpDGFlS4WSU1m0QSVDg+fbilROwxxDEVReGO1dIkTZ8/tgTc1PlOgppPmt9YW4pFGQtFE/15TKDMFasgnW0pwuuXAFk3zqbQ2a8r3+8plRk/RZD/9VsHeHLvaYZySZpPmvTl2fpHBBKIZ5Ja7pKyZRpTXuPlyh7Qyi6ardHj4ZHOJ2mEIwOVWeHuttlsIhffQ8hMLzSbN/16j3V+a8D6fbCmRaXg14IttpVKKSjSbz7aWUG6X1ma1LdtTRnaZ9isfCO+wObNas33jNZk078+zsylDm78w4Z0yS5yslicXqqpxefhsa4naYQgfUunw8LH0bVadHNeiuX22RZsVcjSZNH8lj29FC1i4WZsHob9YuquMEqmYIZrZ4u0lMmZBRZvSqzgs9dZFM/vlUCU5ZdqboExzSXOVw8OK/TLjk2h+O7PtHMjT7gADX+byKHwsNy2iBZTZPayUa4ZqpJVZtASPAou2lagdxgk0lzT/sL9c+jyKFrNYnmKo4tfUSnKkz6NoITK4VB05ZU7WpUq3N9EyvtldhkNjT5E0lzRL1wzRkr7fV06FlKlqdculVrZoQbtz7PyWX6N2GH7nm11lUhZWtJgyu4eff9PWUyRNJc0H8uzsy5MTn2g5dqfC8r2SwLWmyho3aw5Ja5RoWVJWsnW5PQpfy+9ctLAlu7T1GdNU0vzVTm39coRv+vGAtu5cfd3Pv1VQ45LmKNGyvt9XjksmzWk1G9KqyK+QLleiZW3NrCazRDsDTTWTNFc7PTLdsWgVO7KqKamSk31r+U5a9kUrKKl2sz5Nnmi0lp809thc+C4tddHQTNL8a2ollTIAULQCjwKrpbtAq8grd2q2SL3wPVJFo3V4FIV1h+UcKlqHlrr3aSZpXiMHoGhFqw7KxbU1fL+vXAYKiVazIb1KZv5sBXty7BRLzXXRSnbn2DVT418TSbPbo7AutUrtMIQf2ZRRLaUNW8EqmYVRtKLSajd7c2UweUuTJ3WiNXkUWKuRhlVNJM27su2UauQuQvgHp1thvdQXbVGVNW725cpkMqJ1yXHd8uTJsGhtaw9p4+mwJpLmX2XwhlCBdNFoWduOVEvXDNHqfpWkuUVllTpJlWmzRSvbkF6liYlOjGoHALVz14szYzPrCQ00YDbqUJTamZmkrNeZWZdahcujYNTr1A7FJ23OkAGAZ0qvg7BAA+FWAw6XQkWNh5JqN3JkN96+vBpKqlyEWTVxefM5azTS4udNLEYd4VYDZoMOBXC4FAorXWggB/Qa1U6FLZnVDEm0qRqH6meV8hq3TGjSAB3QMdrCoPZWBnew0jUugEBT/YcEiqJQUOnmSImDtCIHaw5VsjmzGqfULT2lKoeHw4UOUqItaofik7Zkys1wQzpFWxiaaKV/Oyttw0xE2Iwn3MS53AoZJQ42pVexMb2KrUeqsTvluD4VjwK/plcxrmuI2qH4pPVpclw3pF2YicGJNjrHWEiJttAhwozhuOPa6VY4VFDD/rwa9ufZWZtaRYHUvT6tX9OqJGnemimPcE8lPNDAFX3DGN8jhEjb6f9UOp2O6CAj0UFG+ra1cmnvMMpr3CzfU84X20tIL3a2UtTeZV+uXZLmFlBS5eJQgTzCPRmTQcfFPUO4ok8YbcPNDa5vNOhIirSQFGlhcr9wnG6FVQcrWLipWBocTmGTJM0tRsYpnJxBD+elBDOpTyg94gMbXN9k0NElNoAusQFAKG6Pwsb0Kj7fWiI3JqewVwOfPdWT5t3Z6v8StCbQpOPagRFc1T8Mi/Hsu50HWwxc3jeMS3qH8vnWEt5ZV0i1tFDVszfXzsU9Q9UOw+dsyayWLgXHMejhwm4h3Dg4gtgQ01lvx2TQMaZzMGM6B7PmUAVvrC4ktUhuUI4lNxMtI7vMSZld+hQcb0gHK9PHxBDfhOPaoNcxJNHGkEQbB/Jr+Mf3ufI5Ps5v+TW4PcoJrfatSfWk+XChfCiO1S0ugGcuiic6qPn+NEa9jqv6h3Ne52Dm/ZQvMzkdQ8pTtQyZ0KS+dmEmnpwQ3+xPNYYnBzGgvZX5PxeweEdps27bm2UUO6hxeZrU6CBOtF8DLX1aYjPruefcKMb3aN6Gl5RoC69Oacd/Nxbzn1+LpJvl72pcCmlFDpKj1Hs6rPoZ5aA8wq0zvnsIr1zRplkT5mNFBxl5+qJ4/nJeNDL2rdbhwhpqZDRGs/stX25GjhrZ0caCa9q3WDcgi1HPX8bE8NSEOGxm1U/pmuD2wGGp8NDs9kvLZ53ucQG8M7V9syfMRxn1Oq4fHMGCq9uREHr2Ldi+Ru3PoKpn2PIaN/nS8R29Du4dFc0jY2Mxt0LLyKW9w3j24nhMBsmc3R44IBeCZqUoiiQsv5vUJ5RnLorH2grJ7OiUYOZd1ZawQEOL78sbyI1b81M7YdGK/u0CeenyNsQEt3wymxxl4V9XtiUxsuHxD/5gf566TztUTZoPSyszAHeeE8UVfcNadZ8jkoN4ekIc8vRSG4MLfElOuYtKmW2RUZ2CuG9UNDpd692cJkVa+MekNtLijCTNLWGfygmLFvROCOCFSxJOqGDVkiJtRl6e1Ia2YdLirPaNm6pn1kMFclK7tFcoV/YPV2Xfw5ODuHlopCr71hLp19y8zmbiA5e9gt8+f551T49m1cO92DJ3CmXp20+67v6FT/DT9M5k/vjuabfpcTtJXTaP9c+ez88P9WTj7IkU7fm53jq5G79k3VPnsvqxQRxcPKvea/bCTH59fhwu+5mPAegZH8BjF8a2asJ8VKdoC09fFI/Bz/NmSZqbV44MAiQh1MQLlySo0lc+wmZkjtwQ1w0GVIuqv/2Dfp40d421cM+oaFVjuGZgOH3aNFwex5dll0k5vuaUWXzmSfP+//2N4v2r6Tp1DgP/uoTwLiPYPv8makpy6q1XsH05ZWnbMIfGNLjN1KVzyV77Pzpd8QSDZnxNwvBr2PX23ZRn7gbAWVHE/oV/I/mSR+h1x1vkblhE4a6Vf8T0yUySLn4QY0DQGf0s4VYDz01U58J61MD2Vq7qp87NuFb8VlCDR5EBVM3F37tc6YC/XhCDzaJe96f4EBN3joxSbf9aYHcpql6z1e2e4ccHoQ54+PxY1fsV63U6HhsX69d3r1JQvnlllpzZCc3tsJO//TuSJz5MWMdBBEZ3IHH8fQREtCVr9X/r1qspyeHAZ8/Q7fqX0OkbfkyZu3Ex7S+4g8juowmMak/COdcS3uUcMle+DUB1YQaGgGBi+l9ESPvehHUaQlXuwdr3bvoKvdFMdJ8Lz+hngdrxCVroV3zT0AjiQ1QvkKQau1OhuMqtdhg+w9/Pk5P6hNK3rVXtMLi4ZyiDOqgfh5rU/CyqmikdOcOLqy8Z0zmIjhqZVCM2xMTtI/y3m0ZRlRtFWqSaTWbJmd0MKx4XeNzoTfWPB70pgNJDm35fx8PeD/9KuzG3YotPadR2PS7HabcZGJ2Ix1FNeeZunJUllGfswBbfBWdlCanfvEKnK548o58DoE+bQMZ0Dj7j97WEo1U1/FmRJM3NprDSf5Pm2GAjt4/QTgvvw+fHEGDy34H8hZXqHdeqJc2KolBq988TmkEPNw/TVpI6vnsIkTb1W8fU4HQrlFb752exJZxpRRxjQBAhif1I+3Y+NaW5KB43uRsXU56+DUdZPgAZK95ApzfQ5twbGr3diK7nkPnjO1Tlp6J4PBTtW03hzhU4yvIAMFlD6XrdLPZ9+Fc2/3MysQMvI6LbSA4tnkWbkddjL8pk05xL2fDiReRvXdaofd42XFvH9eAONga2999WqSI/TvSamz8nzVcPCCegFQf+NSQm2MSE7v4746Wax7Vqz+4qHB6/nT77wm4htA3TVvkYs1HPZb3DeGttodqhqKKg0k2Y1X8fZTen8pozHyzUdeoc9v33UdbNHAl6A8FtuxPTfyIVmbsoz9hJ5s//YcBDi85oYF3Hyx9n///+xoa//wl0OgIj2xM35HJy1n9et05U73FE9R5X933JgfVUZu+j0+Qn+fW5sXS74WXMwVFs/udkQjsOwhx86qQ4MdJMzwTtjQ+4vE8oG9P9c1peaWluPv76uwyy6BmvwQR1cr9wFm0r9cuZVwtV/CyqliWU+/Eo3D91094BCLWtze+sK/TLm5nCShedNNJdxttVnEXSHBjVnr73foi7pgqXvQJLaAy7372fgMi2lB7ciLOikHVPj/7jDR43Bxe/SOZP7zF05sqTbtMcFEHPW1/D46zBWVmMOTSWw1/9g4DItidd3+NycODTp+g69R9UF6SheFyEdRoMgDU6kbK0bUT1HHPKn2FsF210yzjekEQb4VaDX/bvLa7y39bR5uavLc1jOgdrqpX5qIRQE33bBrIl0/9mX1Xzs6ha0lzmp4/Dw60GeiYEqB3GSUUFGekRH8COLP+rxSmT7DQPu9PTpClfDRYrBosVZ1UpRXt/IfmSh4nucyHhXYbXW2/769OIHXgpcYOvaHCbepMFS1gcHreT/O3fEt13/EnXS/v2VSK6nUtwux6UZ+5G8fxxjvK4a/tdn4oOGNtVm0mzQa9jdEoQi7b53zTbRSr2ffQ1BX76uxzT+cyq57SmMZ2D/TNpVvF6rV7SXOOfB2D/toHoVajd2lhdY/0zafbXVpTmVnaW4xSK9qwCFAJjkqguSOfQ4llYY5KIG3IFeoMJk61++TSd3oQ5OBprbHLdsr0fPIw5NJbkiQ/VxpK6jZrSHILadKOmNJe0Zf8CxUP7MbedsP/K7APkb/maAQ8vBsAakww6HdnrPsEcHEVV3iGC2/c6ZfwpMZZWmR3sbPWIC/DPpFlampuFR1EabLUvObiBjB/epCJjF46yPHpMe5Wo3mPrXlcUhbRl/yJ77ce4qksJbt+HlMkz6w3s9bgcHFz8Inmbl+Bx1hCeMoyUK5/CEhZ32n0f+eVDMn94i5qyPGxxKXSc9BhhHQfVvZ7xw1tk/PAmAO0vuJ22o2+ue60sdRsHPn2K/g98ik5ff1yPXgddYrXZyAXQNU67sbUkv+ye4a8tzX00ULLmdLpq+ATRkvy9aH9zOZuuGQAuezmHl7xETUkOJlsYUb3HkXTRA+gNjU9E7cXZoPvjMarHVUPq13NrS8tZrER2G0XXqXMwWut3j1IUhf0fP0HHSY9hsNQenwZzAF2vncWBT5/G43KQcsWTp71wp2i8a09KjH8e1yV+ep1pbpUOD+4GDm13TRVBCV2JG3w5u9+594TXM1b8m8wf36HLtS9ijUki7bv5bH/tZgY9tqyuFvpvnz9P4a4f6H7DPzHawjm0+EV2vHF77XgG/ckHqudtXsrBRX8nZfJMQpL6k71mITsW3MagR78mIDyBiqx9pH7zCj1vWwAo7Pz3nwnvMgJbfGc8bicHPnmSlCnPnnT77cLNrTrz35lKijBjMuia9HTPG6k5cF/Flmb/TFJig7U92KxrrLYv/i1FzRmGfMnZDAIEiOk3gZh+Exq9/sn6Mfe994N634d1GsygR79pcFs6nY5+9//vhOWRPc4jssd5jYonMVJbA3uP1y7cRIBRh93lX59zf0smWoqrEb/HyO6jiOw+6qSvKYrCkZ/fo/3YO+tqn3e9bjZrHh9G3qYlJIy4Gld1OTnrP6XrdbMJ7zKidp2pc1j31CiK960hotvIk24788d3iBsymfhhVwHQ6fK/Ubx3FVm/fETyxIeoyj2ILaEL4Z2HAWCL71K7LL4zGT+8SWjHQYS0733SbXeO0fb10GjQ0THK7Hez2rr8cUbAs22R8nahAdou69YmzIyKE5mpRpLm5lHpp92u4kO02zUDaicxahuu7cS+JTTUOioax9XE36O9MANHWT7hXc+pW6Y3mgnrNJiy1M0AlGfsRHE7661jCY3FFp9St87xPC4H5Zm7iOg6ot7y8K7nUJa6BQBbfGeq81OxF2dhLzpCdX4q1rgUqvPTyP11EYkTpp8y7g5ecMx0iNB+jM1NzaRZtWZPlSfCU02IBmYKa4jFqMfl8K+rjVxcRVNo+RHuUYF+OBmC3Aw3j6b+Hh3lBQAnlGw0B0diL8qqW0dnMGGyhh63ThSOsoKTbtdZWQweN6bg+hOPmIIj695ji+tE0kUPsH1+bT/mpIsfxBbXiW3zbyR54sMU7/2F1GX/Qm8w0vHyx+v1hfaGCUQC/LCVS83rtWpJs9rTR6slOED7H/CHdzyHp7JM7TBaVYhxFIy9T+0wvJ6/Htfe0A3gmn3/ZtKBXWqH0arMbZOAf6gdhtdrvglT658fFEWBBgbGN2ad47eLUn9RwohrSBhxTd33Oes/x2CxEZLUj1+fv5D+D35GTUkOe977C0Oe/AG9sbb11qjX/vnM6IfnXDVn8FUtaTb74R8awO0FF1fdT59BZYXaYbQqW2yY2iH4BL9Nmr2gRdOydSX2javVDqNVWfoMangl0aCm5o7m31uCHeUFWEL/mNrdWVFU1/psDo5CcTtxVpXWa212VhQSmtTvpNs12cJBb8BZnl9vubOisG6fx3NWFJH23Tz63vsRZWnbsMYkYo2u/VLcTqryDhOU0KU2Xi+4Xjua2nfGC6lZgUy1Zk9/vbhqvdalp6YGd7X/zR6mM2h7gKa38Nfj2uEFA+ycRfkNr+RjTlVxQZyZprZmBkS2wxwSTfG+P27aPC4HJb/9SkhifwCC2/VEZzDVW6emNI/K7AN16xxPbzQT3LYHxfvW1FtevG81IYknT7R/W/R32o66qbYajseN4v6jlJ7icderx362A5tbkzfE2NzUbF1Xr6XZ6LsX1yhPGR1d2bSzHyGm/AihpRkEFKSjz00lpPtcGDtR7RBPqWLvdvD430GoM2l7IJe38NekOavUqXYIp6W43dgz087oPe+VKPznuNLO4Xr4tN3J/8aFLoXXi2G/A464YFIw3B1Rf92N1Qr/VwQlbhhuhQcjwfR7q1GFR+GubJgTC7HNdH3QGeVmuDk05rB211RSnf/HZ8xelElF5m6MtjACwhNoc+6NpC9/HWt0BwKjE0lf/joGcyAxAy4GwBgYTNyQyRxa/CImWxhGaxiHFr+ILb5zvcmNtr16A1G9x9Jm5PUAtB19M3s//CtB7XoSktiX7LUfYy/Ortcd46iifaupzk+l63WzAQhu35uqvEMU7v6JmpJs0OsJjPmj9vtv+dqvSuENMTY3NS8z0qf5LJgVJ8nuPBIdWcRVZRJZlomtKB1TXhqe7DTcZSUnvMf9+1fV/h2aTprLd21ROwRVGAK1XT/bW3jzcd0UBwu0feGqOrQPj+PMY0w01SaxR53u0aQTCDXAdaHw2UmGRHgUhRcK4OpQGBQAT+fD0gq47PeJFP9dDBODmy9hBjCFRjTbtvyZuRGDzcrTd7Lt1evrvj/4xQsAxA6aRNfrZtHu/NvwOO0c+PRpnFWlhHToQ+87366r0QzQadJjHDQY2P3udDxOO2Gdh9Hz2ln1nhhUF2TgrCiu+z6m/0U4q0pI+/ZVHGV52OI70+vP/yYgok29+NwOO799+gzdb/wnOn3tz2MJi6PT5U+w77+Pojea6XrtLAzmP2qa78/T9kRflTVujpRo+4a9JajZ6Cp9mk8hxlNCJ1c2bauziK7IJLQkg4CCNHQ5aThzM8F9YjeLxnx0i9f+SOLdjzV/wM2kwk+TZnNUbMMriQb5a9K8I0vbU9mWbl57Vu8zABGN/JvGGXXc83uOuqzixO4qpR4o8cClwWDW6RhmVUhz1L62066w3wH3NXOOa4qIbHgl0SCrWd9gne+wlCGMmrv/lK/rdDoSx99H4vhTD7jWmyykXPEkKVc8ecp1Tlajvc0519HmnOtO+R6onbBo8N++PWF5/LCr6mo8H6+w0k1hpYtImzafWOzPr0H7HcOaX4RVvb+Hit0z1K0iYcFJR1cOHRxZxFceIeL31mJjbhqe7FTcFSc2lTTHhKwla1fiKCrAHHHyQQpqUhSFwlXL1Q5DFeaYeLVD8AlWLyi91hIKK90cLqwhKVKbkyHkr/v5rN53xAVXZSqYgK4WuCUMEs6yDFeYHiINsLEaBgQo7KiBcTZwKgpzi+DhSDA08wAfU0R0s27Pn0XYjJrvhtQSth2pZkznYLXDOKntR7R9s95SIm3qjVVQLWkObYXSa/GeEpKdWbStPkJ0RSYhJRkE5P/eWpx35IS+uwqNay1uCsXtJv/bRbS55rYW3tOZK17zA/b0Q2qHoQpz9KmnSBaNF241+OW0rgBLdpZx7yjtJWnlNW4eTn6coXddxpCM5QT/ugRXbmaD7+tqgUcioa0Jit3wYSnclwNvJSiEnsUTBZ1OxxNRCq8Vw6tFMDgQxgfBR6XQPwAsOrgvR6HUXdsf+rKQpifQ5kjt/T28VZTN4JdJ87e7yzSbNH+7p1ztEFShZsu/anuODm76wKsAxUkndzYdarKIq8wkoiwTa2E6xrw03NmpeE5SNq05WoubKvfL/2kyac7675tqh6Aai7Q0NwudTkd0kH+2SC3bXcYtwyKxmrXV2r5ibzkOD/xs6cnPnXqi6zidc5y7GZq+nJANS3DmZJz0fUMC6yet3S0K1x+B7yrhypCzi6VXgI75xxxqGU6F7ythQTxMz4UrgmFQINyaBb0CFDqam5Y4m8K190TPW2m1i0JL25BeRXaZU3Ozfm5Kr/LL8yz4adIcFmjAYtRR00CpprbuQpJdWbSpziK6LJPgknQs+emQm4YrL+uEquut0VrcVCVrV1K05gciho9RO5Q6jvxc8r9bpHYYqpGW5uYTF+yfSXOlw8PyvWVc2jtM7VDqVDs9fLChqN4yRadjlbkHqzr1QNfxfkY49jAsYzkhG5fgzE4/5bYC9TqSzApHmulPqygK/yyEO8LBA/zmgHOtEKDX0TtAYbsdOjZxhmCTtDQ3mwg/TZo9Cny6pURzT5EWbi5ueCUf5ZfdMwBigowUFVfQyZVN+5raShThJRlYC9Mx5KXhzk7Dc5KawVpoLW6q3557iEFLNtaN4lXbgWcfQHE41A5DFXqzBVNouNph+IyYYP+8uAJ8tLGYcd1CNDOt9sJNxaetDa/odPxi6c4vnbpDp/sZ7tjDsIzvCdv4Fc6s+iXqHIpCuhN6NVO37a8rIMQAw606yn/vznP03O5SahPpprLEJDTDVgSom6io7asdpVzWO5R24U28i2smG9Mq+TXN/+ZTOErNGzhVr24P/vooRV9/fEJrsYfmOWFqWcXureR8/j7xk29UOxQKf1xG7pf/VTsM1Ugrc/OK1dhjzNaUW+7i3XWF3DlS/VapwkoX/zvD1qg15m6s6dgNOt5L9ed/5bIoPcmH11OQfYQPSqHKAxf+XiHszWKFAjfMiPqjC8VvjtpzebUCpe7a741A4nHdLIrdCh+Wwiu/H3rBBh3tTQqfl8GAQIUt9trSdU1lTUpp+kYE4L/dM6B2ZsBZy3P5vyvbqjobHdSWmZuzIk/VGNTml90zAKzR0RSpOIe42g48+wChA4aremJ3V1ex7/G7VNu/FgQmd1Y7BJ8S68ctzQCfbClheHIQfdoEqhaDR1F48btc7M6zP78eLnfw1NaNOCuLsdki6JUQzuvBdmKLswAodEPecY/9/pz9x//3O2BFFcQa4KO29dd7tQiuCoHoY+qtPhIJswrg8/La17pampacmGMTMFhtTdqG+IO/H9c7s+18sqWEKf3VfSr52i8F5Jb7wvP2s6fmZ1HdpDnJv5MVV2kx22+5hAGL1mIKDWv1/SseD7vun4o9M7XV960lwd37qh2CT4lthkG+3syjwLPLcvjXlW1VGzz01ppCNqQ37fFt9xvnnrhdYIhzPyMyl/P4xq9wZh6u9/qKDo1LdB+PPnG9rhYd77Q5ycpnyZrcpfk2JugYrc1yiq3prTWF9EkIpGtcQMMrt4CfDpSzZOdJZg7yI+FWg6otzap2vJOTWu1MXTvvvgqPq3XvHBVFYd/jd1Hw3Retul8tCpKkuVklRmqj35+aCipcPPBZJnnlrT8g8ssdpXy4seUGCa03deblpLt54splfHvPEiouux9Tu+SG39jKgrr0VDsEnxJsMZAQ6t83xA63wl8XH+FwYevPALoxrZLnvs1t9f1qTZcYdW/eVE2ag7r1UXP3mlH8y/fsuO0yXCcpkdcSFLeb/U/cTdZHb7TK/rQuuEc/tUPwKZE2I1F+PGjoqOwyF3/57AgFFa13Q/ze+kL++UPr9Xdcb0rh5aS7eGLyNyy7ZwkVk6Zjatex1fZ/OjZJmpud2gmLFpTZPUz/7Aj7cltviu01hyp47Ktsv6x/f7zOMeq08h+latJsjoohUAZqAFC48ms2Tx5JVdrBFt2PPSuDzdeM4cgHr7fofryFwWrDKn2am53aJzatOFLq5M6FGS0+zXaNy8PT32Tzzroi1abV/dWUwsuJd/LE5K9Zes9Syi5/AFMH9c7vQV17q7ZvX9VZkmYASqvdPPD5Eb7d07JdJdwehY82FvHE0mwckjAD6n8GVa+LFDpghNohaEbFnm1svHggRz56A8V96jJRZyt/2SJ+ndCP0l9XNfu2vZWta2/NlP3zJV1j5eJ6VH6Fi+mfZfLmmoIWaSnam2vn7o8zWbm/dZ5UNcYmUyfmdvgzT1y+hKX3fk3Z5Q9iSmy9m1N9QKA8QWoBcjP8h0qHhxe+y+XRL7MorGz+p0lpRQ7u+TiDN1YX4vb1cmJnQO2kWfXhsGEDh5Pz6btqh6EZrvJS9j12Bxlvv0LHv/6d6HGXNnmbxWtWkvraixSvWt4MEfqW4B591Q7BJ3WPV69yhBa5PfDBhmJ+PFDB1QPCGdc1GLOxaTdr2WVO/rO+iG/3lOHRcCPUJmNHNnXoCB1up7/zEOdmf0/kpiU4U/e12D5D+w9Db5a+9c1N7YRFi9YeruSm99O4eWgkE3qEENDEGu2l1W4+31bCRxuLpTvGccIDDcSoPNBcpyjq1nyr/G0P6y/ooWYImhbceyCxl1xN9IWXE9gusdHvcxYXUrRqORlvv0LZ1vUtF6CX6zr7LRKuulntMHxOlcPDxa8f1HQyp6Zwq4FJfcL4U7fgM7oIVDk8rD1cyU8HyllzuBKXF7dA9XMd5tys74navATn4b3Nuu2kvzxF0v1PNus2Ra1r3031yxk/GyMkQM+fuocwvnsISZFndoOxI6uaZbvLWLGvHHsDMyX7qyEdrMy6rBlL7JwF1ZNmRVFY1S8aV0lRwyv7uaAe/QgfOoqAtokEtE0ksG0i+kArnho7zuJCHHnZlO3YSPGalVTs3nrCpDHiRCPWZWCJU/cg9FW3/zed/XmtP8rc23SIMNMrIYCusQFEWA3YLAasZj2BRh0l1W7yK1ykFTnYn1fDxvQqn+zb2MeVyuis74nasgTnoT1N3l6///5A+LDRTQ9MnOD5b3NYvrdc7TA0LyrISJcYC51jLHSOCSDKZsBs1KMoCg63Qk6Zi/15dvbl1bA/r4bS6ubvkulrpg2L5IbBEarGoHrSDLDtlksoXLFE7TCEnwnq1ofB32xROwyfNf/nfD7eUqJ2GMLL9HGlMSr7e6K3LMF5cPcZv19vtjByezGGAOl/2xJ+2F/OM9/kqB2G8EP/vrY9KSrXC9fECKiwwSPVDkH4ocjzxqsdgk8b0TFI7RCEF9pm7MD/tbuFJy5ZxBf3Lqdk8iOYOja+C19w38GSMLegIR2sNLE7vhBnLDbYqHrCDBpJmqMvnKR2CMIPRY65SO0QfFqvhADCAqVeszh7243t+b9203jiks9ZdO9yiq+cgSml12nfE3nuuFaKzj/ZLAZ6qzhFvPBPw5JsaocAaKB6BoA1sRNBPfpRsUselYvWYQyLILTfULXD8Gl6nY7hSTa+3u3f076K5rHD2J4dbW+GtjfTw5XBebkriN2yFOeB7fXWk0aYljciOYjNGS1be1yIYw1PlqS5npgJkyVpFq0m8twL0RmkFbSlndNRkmbR/HYZ27GrzU3Q5ia6uzMZk7OC2K1LMbmqsaV0Vzs8nzc82ca/fspXOwzhJ6xmPf3aWtUOA9BI9wyAmIuuVDsE4UekP3PrGNjeSqBJp3YYwoftNrRlXpsbeeKij0mduUztcPxCfIiJpEipgy1ax8D2VkwGbVxHNJM0WxM7EdS9r9phCD9gsNqIGtv0SWNEw8xGPUMStfFYTfi+IT3j1Q7Bb4zQyONy4fu09FnTTNIMEHPRZLVDEH4gZuIUjEHBaofhN87tJFU0RMtLjjSf8YQS4uyd30XOoaLlBRh1nCNJ88nFXHSV2iEIP5Bw9W1qh+BXzkm2SRUN0eIkiWtdSZEWesRLaT/RskanBGGzaOf6oamk2ZrYibBh56kdhvBhtq69CO03RO0w/IrZqOeiHiFqhyF8mMmgY4J8xlrdxfI7Fy3s4p6haodQj6aSZoD2t0xXOwThwxKuvlXtEPzSpb1D0WtjHIfwQed3DiLcqpliUH7jvM7B2MyaSyOEj0iMNNMzQVs1wTX3aY88/2ICk1LUDkP4IH1AIHGTrlc7DL8UE2zSTJ1N4Xsm9wtXOwS/FGDSM767tDaLljGpt7ZamUGDSbNOp6PdTfeqHYbwQTETJmMKDVM7DL81qXeY2iEIH9SvbSCdNDC9rr+6vG+YPEUSzS7YoufCbtq7IdNc0gwQf+XNGEPC1A5D+BKdjnbTpqsdhV8b0N5Khwip7Sqa1+R+YWqH4NcSQk2ameJY+I6LeoYSYNJeiqq9iKito5twjfQ9Fc0netxlBPfsp3YYfu/yPtp73Ca8VxtJ2DThqv7SPUY0H5NBp8muGaDRpBmg7Y33ojPKwA7RDHQ6kv7ylNpRCGB89xBig+W4Fs3jir5h6HXSN0BtfdoEMqiDNqY5Ft7vkl6hxIaY1A7jpDSbNAcktCNhyi1qhyF8QMyEyQR17aV2GILa8nPThkWqHYbwAbHBRi7uqb0+j/7q9hFRyO2LaCqrWc/1g7T75EKzSTNA0l+exiAzt4mm0OtJmj5T7SjEMcZ2DSY5Uvo2i6a5bXgkZqOmL2F+JSXawpjOMvunaJor+4URpuHykZo+45ijYuhw5wy1wxBeLOaiq7CldFc7DHEMvU7HrcOltVmcva6xFpkBUIOmDYtE7mPE2QoLNDBF4/3jNf/xbnfLX7AktFM7DOGFdAYDSfc/qXYY4iSGJwfRK0Gm4BVn586R0eikL7PmtAkza24GN+E9rh8cgVXjk+VoOzrAEBBAx4efVzsM4YXa3ngPtk5d1Q5DnMKfR0SpHYLwQiOSbfRpo61ZwsQfbhgcQYBJbmjEmYkPMXJJL+3fcGk+aQaIvew6gnsNUDsM4UUsbdqT/NBzaochTqNnQiDndpI+kKLxjHq44xy52dKyCJuRqzX+iF1oz63DozAZtH+z5RVJs06nI+WJl0Eex4lG6vLMPAxWqd+qdfeNisam8cdxQjsm9QmjXbgMItW66wZFkCSDfUUjDUuyec0YBa+5WoUNHknbG+9WO4wW81GpwvlpCq8WKXXLzk9TTvq1sFQ55XZcisJ/ShSmHlH4U5rCbVkKv1bXX//7CoWrMxUuy1BYUFz/tRyXwg1HFCo9p96H1sVcdCVR51+sdhiiEaKCjPxZWg5FI7QLM8kAUi9hMuiYMTYWg9dkGEItwRY9D50fo3YYjeZVH+mOM2ZhTe6idhjNbm+NwtJySD6ulvcnbet/PRwJOmDkaWrIv10CSyrg3gh4OwEmBsPMfDjgqE2CS90KLxXBn8PhxRj4rgLWVf2RIM8thNvCwab3zlZ9Y0gYKTNfUTsMcQYm9gyRPqritPQ6ePTCOCxSmsFrdIkN4JoB0k1DnN69o6KJtGm3xNzxvOoMZAgIpPvL76EzGNQOpdlUexT+XgAPRELwcX+NCIOu3tfqKugbAAmnGWTxfSVcGwpDAnUkmHRcEqxjYAB8Ulb7erYLbDo4z6ajq0VH3wBIc9a+tqJSwaSDkVbvTJgBOs54EUtMnNphiDOg0+l4ZGwsgTJ4SJzCtQPD6R4n1Va8zY1DIqUmuzilEck2xnXzrgmKvCppBgjpO5gOd/lO7eZXimBoIAwIPH3CUORWWF8N4xsYN+VQ4PhTlEUHO+21/29jhBqltuW5zK2wzwHJZihzK7xbUttC7a3CR5xPwjW3qR2GOAsJoSbuGhmtdhhCgzpGmblxiHTL8EYmg44Z46SbhjhRSICeB8Z4T7eMo7zyo5x435ME9eindhhN9kOlwm8OuLURT7C+qwCr/vRdMwAGBcCn5ZDpVPAoChurFdZUQ5G79vVgg45HomBWAdydA2NtMChQx4JiuCwYclzw5yyFW7IUfqr0nn7N5tgEerzyodRu9WITe4UyNLGBD7jwKyaDjr9dGOcVo+rFyXWOCeC6gV7cGiNahLd1yzjKK5NmvclE93/+B73ZonYoZy3PpfBqETwaBeZGJHrLKuB8W8Pr3h1R25p8cxZcmA7/KoILbbV9Ao86x6rjzQQd77fRcWOYjq12hcNOuCgIniuAuyLgqWh4qRCK3dpPnHUmE73mf4w5yvvuWkV9j46LIy7E+06komXcPDSC5CjvPc+LWjcMjqBnvHSvEbUu6BLM2K7e1S3jKK9MmgGCOvegy99fVzuMs7bfASUeuCMbxqYpjE1T2FYDi8prv3crfySr2+0KGS6Y0IiStmEGHc/G6FjaHj5qA+8mQKAe4k6RhzgUhVeKYHoEHHGBW4E+ATramXS0NcGemmb6gVtQpxmzCB0wXO0wRDMIDTTw3MUJBBilZdHfjUi2yUAyH2E06HjmonhiguSG2N91ibHw8AXe28DltUkzQPzkG2n/54fVDuOs9A+AN+PhjWO+uphrW5PfiAfDMS3K31RAZzN0NDc+kTDrdEQbdbiBVVUw/BTFCT4ogcEB0NmiwwO4j3nNpYDnbH64VhRz0VW0u2W62mGIZtQp2sIjY2PVDkOoqH24ib9dGCfdrXxIhM3IcxPjscgNsd+KtBl4bmKCV1fB8d7If9fxkReIGnup2mGcMateR5K5/leADkL0kHRMclzpUfi56tStzC8WKLx5TK3lPTUKq6oUspwK2+0KM/JAAa4+yeyUqQ6FH6vgprDa79sba0vafV2usK5KId1Zm8hrlbVjV7rOflPtMEQLOK9zsLQy+imbWc9zExOwyqQ3PqdzTIDcEPspk0HHsxcnEO3lTxu8/qyk0+vpPvd9grr2VjuUFrGysjbpPe8Uk9vluaDwmOZhh1Jbq3laVm195igDzI2FoOPqLiuKwstFcGc4BP7+mkWv46+R8H4p/KOwtpJGtEZbBYzBofR6/VOMNpmG2VfdNiKSIR1kYKA/Mejh6YviaS+z/vmsMZ2DmTpIboj9zUPnx/hE2UidoijaH+nVCPYj6Wy8dAiOgly1QxEtTG+20Oc/ywgfOkrtUEQLK69xc9f/MsgocaodimgFD50fw8U9T/JYTPgURVF4fEk2qw9Vqh2KaAVX9Q/zmZKiXt/SfFRAm/b0WvA5eov338mI09Dr6f5/H0rC7CeCLQZeuryNVNTwA9cODJeE2U/odLWlBFOipTKKrxuRbOOOc6LUDqPZ+EzSDBA6YBi93lgkibMP6/LcfGL+dLnaYYhWFBNs4p+Xt5WR9z7syn5h3D7Cdy6somFWs545k9qQJDMG+qxB7a3MHB+H3ocG9PpU0gwQOepCer+5GH3AKcpFCK+V8uQ/aXPt7WqHIVQQH2ri5SvaEGkzqB2KaGZX9Qvj7nN949GtODNhgQZemtSGduEmtUMRzaxv20CemxiP2YsrZZyMb/00v4sYOVYSZx/T8ZEXaDftfrXDECpqG2bm5cvbEh4oibOvuKpfGHdJwuzXImxGXr68Le3CJHH2Fb0TAvi7l5eWOxWfGQh4MkVrfmD7LZfgqa5SOxRxtnQ6Oj06i/a3P6R2JEIjDhXUMP2zTMrsWq8iLk7HlwYHiaYrrHTx4OdHSC1yqB2KaIL+7QL5+8QEAky+lzCDjyfNAMVrVrLtlomSOHshndlMtznvEHfpNWqHIjTmYH4NM77MIr/CpXYo4ixIwixOpqTKxYOLjnCwQBJnbzSkg5VnLo73yRbmo3w+aQYo2fALO26fhLO4UO1QRCMZg0Pp9cYiwoeNVjsUoVEFFS4e+yqL/XleMNe7AECvg1uHR3LtwAi1QxEaVV7j5qml2WzKqFY7FHEG/tQtmAfPj8Vk8J1BfyfjF0kzQFXqb2yfNpGqQ/vUDkU0wJLQjj7vLCWoS0+1QxEaV+308NyyHKn36gWsZj2PXxjL8GSZkEicntujMH9VAZ9tLVE7FNEAvQ7uOCeKq/r7x4Q1fpM0AzhLi9l555UUr/lB7VDEKQR17U2fd5diiWujdijCS3gUhddXFfDxlhK1QxGnEB9i5PmJCSRHSV1e0XhLd5Yy98d8nG6/SVO8SpBFz5Pj4xjc4RRTFvsgv0qaATwuFwdf+CsZb81VOxRxnKixl9L95fcwBoeoHYrwQou3l/B/P+XjlvGBmtKnTSBPXxRPmFQ9EWdhR1Y1Ty7JprjarXYo4hjtwk38fWIC7fxsynu/S5qPyln8X/bOuE0GCGqA3hJAp8f/Qdvr71I7FOHltmRU8fy3ORRUygVWCy7uGcL00TEYfbyfo2hZeeVO/vZVNgfyZfyCFgzpYOWJ8XEEWfzvRthvk2aAin072f3AjVTs2qJ2KH7LltKdHv/6L0Fde6kdivARJdVuXvwuh3WpckOslpAAPfePjuH8LsFqhyJ8hN3pYf6qAr7aUYrfJi0qM+rh+sERXD84wqdm+TsTfp00A3icTlLnPU/a/BdQnE61w/ErCdfcSsqTczEEWtUORfigz7eWsGB1ATUuvz7FtbqhiVYeviCWSJtMey6a38b0KuZ8n0tuuZSbbE2doi3MGBtLp2j/Hpfg90nzUeU7t7DnoZup2Ltd7VB8njEkjK4vvkHMhMlqhyJ8XEaxgxe+y2V3jl3tUHye1azn7pFRXNQzVO1QhI+rcnh4bVU+X+0sUzsUn2fUw9RBEUwdFCHdrJCkuR6Pw8Hh/3uG9Ndno7jkLrbZ6XTEXXEjnWa8iDkqRu1ohJ9wexQWbi7mP+uLsEurc4vo1zaQR8bGEhciUyGL1rMxrZI5K/Kk1bmFdIwyM2NcHCl+3rp8LEmaT6Js+0b2/e1OyndsUjsUnxHcawCdn/4Xof2Hqh2K8FP5FS7eWF3A93vLpU9kMwkLNDBtWCQTe4ag89M+jkJdlTVuXv+lgKW7yvDIgd0sTAYd1w4I5/rB0rp8PEmaT0FRFPK+WsjBfzyOPf2Q2uF4LVN4JMkPP0/C1bei0/vu1JrCe+zOsTPvp3zpstEERj1c3ieMG4ZE+OUIeqE9hwpq+PeaQtYelomOzpZeB+O6BnPz0Ehi5anRSUnS3ACPw8GRD18n9V/P4SwqUDscr6EzGkm4+laSH3oOU5hMmSu0RVEUvt9XzhurC8mvkEe7jaUDRqcEccuwSNr6WX1W4R22H6nmjdUF7MyWm+IzMSLZxq3DI0mKlK4YpyNJcyO5ystIf2MO6W/+U2o7n4Y+IJCEKdNod9tDBLbtoHY4QpyW3enhs60lfLq1hOIqqe18OgPbW7lteCRdYgPUDkWIBq0+VMG/1xSSWuhQOxRN650QwO0jouiZEKh2KF5BkuYzVJOXTeZ788j677+l5fkYxuBQ2lx/J+2mTZdBfsLrOFwevttbzsebi0kvltKTRxn1MDolmMn9wugqybLwMh5F4bs95fx3UzFpRZI8H6tHfABTB0UwLMl/psBuDpI0nyW33U7eV/8j8715lO/crHY4qjFFxdBu2nTaTr0TY4iUmhLeTVEU1h6uZOHmErYdqVY7HNWEBRqY2CuUy3qHSr1l4fUURWFjehWfbS1hfWqV3w4ErrsJ7htG1zi5CT4bkjQ3g5KNq8l851/kf/u5X5Sq0xmNRIwcR9zl1xM19lIMAXLwCd+zJ8fO59tK+OVgBdVO/zhNJkeauaJvGBd0DcZilIG7wvdkljhYurOMZbvLKK72jy5Z8SFGLuoRyvgeIXIT3ESSNDejmpwj5H71P/KXLaJ081rwsV9tcK8BxE2aSuwl10gXDOE3alwe1h6uZMW+ctalVuF0+9ZxnRBqYlSnIEanBEl/ZeE3XG6F1Ycq+HZPOZsyqnxu5lCbWc/QRCvje4QyoF2glIRsJpI0t5CavGzyv/2C/G8XUbLuR69tgbZ27Er0hZcRd/n12Dp1UzscIVRVWeNm1cFKVuwvZ3NGFW6P2hGdnbZhJkanBDGqUxApMZIoC/9md3rYmF7FmsOVrD1c6bWDguNDjAxPDmJ4ko0+bQKlxnILkKS5FThLiij4/isKV35D6Za11GRlqB3SKVkS2hE+/HwiRowhfPgYLLEJaockhCaV2d1syahiU0Y1mzOqyCzR7gDCYIueHvEB9IwPZFiSjY4yw5cQJ6UoCrtz7Kw9XMnqQ5Uc1nD1Db0OusYGMDzJxvBkG8lRcly3NEmaVVCTm0XplnWUbV5L6Zb1lO/YhMfe+oOOdGYz1qTO2Dr3JHzoKMJHnI81sVOrxyGEL8ivcLEru5qdWXZ2ZldzqMCBQ6WuHG1CTfRMCKBXQiA94gNIjDDL41khzkJptZsDeXb25dWwP6+G/Xl2ssta/8mxXgftws10jrH8/hVASrQFq1nGHrQmSZo1wONyUbl3O1WH9lOdmUrNkTSqM9OwZ6ZiP5LWtIRap8MYEkZg+2SsnbphO/qV0p3ADh3RGWQ2LyFagkdRyCt3kVnsIKPESWaJk4xiB5klTnLKnE2e8jckQE9MsInYYCNtw0y0CzfTNsxEYoSZMKsM9hGipZTZ3XUJdFapk8JKN4WVLgorXRRXuc/62DbqIcJqJNJmIMJmJMpmpF24ic4xAXSSBFkTJGn2Ao6CPGpyjuCpqcZTU4PHUYOnxl777+//R6fDaAvGEBSMwRaMKTQcU0Q0pogo9Ea5gAqhJW6PQpXDQ6XDQ9XvX5XH/Ov2KFiMegJMutp/jTosptp/A0x6wq0GAk1yARVCa9wehZLqo0m0myqHB7ei4PbUvqYDDHodBn3tvzaznihbbaIcGmiQJ0IaJ0mzEEIIIYQQDZCmCiGEEEIIIRogSbMQQgghhBANkKRZCCGEEEKIBkjSLIQQQgghRAMkaRZCCCGEEKIBkjQLIYQQQgjRAEmahRBCCCGEaIAkzUIIIYQQQjRAkmYhhBBCCCEaIEmzEEIIIYQQDZCkWQghhBBCiAZI0iyEEEIIIUQDJGkWQgghhBCiAZI0CyGEEEII0QBJmoUQQgghhGiAJM1CCCGEEEI0QJJmIYQQQgghGiBJsxBCCCGEEA2QpFkIIYQQQogGSNIshBBCCCFEAyRpFkIIIYQQogGSNAshhBBCCNEASZqFEEIIIYRogCTNQgghhBBCNECSZiGEEEIIIRogSbMQQgghhBANkKRZCCGEEEKIBhjVDkAIIY7l8XhwOBxqhyE0wGQyYTAY1A5DCCEASZqFEBricDg4fPgwHo9H7VCERoSFhREXF4dOp1M7FCGEn5OkWQihCYqikJ2djcFgoF27duj10nvMnymKQlVVFXl5eQDEx8erHJEQwt9J0iyE0ASXy0VVVRUJCQlYrVa1wxEaEBgYCEBeXh4xMTHSVUMIoSppyhFCaILb7QbAbDarHInQkqM3UE6nU+VIhBD+TpJmIYSmSN9VcSz5PAghtEKSZiGEEEIIIRogSbMQQgghhBANkIGAQghNG/3KgVbd34/3p5zV++bPn8+cOXPIzs6mR48ezJ07l5EjRzZzdA37IbH12kLGpJ55acCff/6ZOXPmsGnTJrKzs1m0aBGXXXZZ8wcnhBDNTFqahRCiiRYuXMj06dP529/+xpYtWxg5ciTjx48nPT1d7dA0p7Kykj59+jBv3jy1QxFCiDMiLc1CCNFEL7/8Mrfccgu33norAHPnzuXbb7/ltdde44UXXlA5Om0ZP34848ePVzsMIYQ4Y9LSLIQQTeBwONi0aRPjxo2rt3zcuHGsWbNGpaiEEEI0N0mahRCiCQoKCnC73cTGxtZbHhsbS05OjkpRCSGEaG6SNAshRDM4vp6woihSY1gIIXyIJM1CCNEEUVFRGAyGE1qV8/LyTmh9FkII4b0kaRZCiCYwm80MGDCA5cuX11u+fPlyhg8frlJUQgghmptUzxBCiCZ64IEHuP766xk4cCDDhg3jjTfeID09nTvuuEPt0DSnoqKC3377re77w4cPs3XrViIiImjfvr2KkQkhxOlJ0iyE0LSznWykNU2ZMoXCwkKeeeYZsrOz6dmzJ19//TUdOnRo9VjOZsKR1rRx40bOO++8uu8feOABAG688UbeffddlaISQoiG6RRFUdQOQggh7HY7hw8fJikpiYCAALXDERohnwshhFZIn2YhhBBCCCEaIEmzEEIIIYQQDZCkWQghhBBCiAZI0iyEEEIIIUQDJGkWQmiKjE0Wx5LPgxBCKyRpFkJogsFgAMDhcKgcidCSqqoqAEwmk8qRCCH8ndRpFkJogtFoxGq1kp+fj8lkQq+Xe3p/pigKVVVV5OXlERYWVndTJYQQapE6zUIIzXA4HBw+fBiPR9sTdIjWExYWRlxcHDqdTu1QhBB+TpJmIYSmeDwe6aIhgNouGdLCLITQCkmahRBCCCGEaIB0GhRCCCGEEKIBkjQLIYQQQgjRAEmahRBCCCGEaIAkzUIIIYQQQjRAkmYhhBBCCCEaIEmzEEIIIYQQDZCkWQghhBBCiAb8PybRDcOiWNK5AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "\n", "def make_pie(probs, ax=plt):\n", " labels=[n for n in names if probs.get(n,0) > 0]\n", " patches, _, _ = ax.pie([probs[n] for n in names if probs.get(n,0) > 0],\n", " wedgeprops=dict(width=.75),\n", " colors=[colors[n] for n in names if probs.get(n,0) > 0],\n", " autopct=lambda pct: f\"{pct:.1f}%\")\n", " return dict(zip(labels, patches))\n", "def make_plot(results):\n", " nb_cols = len(list(results))\n", " nb_rows = 1\n", " fig, axes = plt.subplots(nb_rows, nb_cols, figsize=(3*nb_cols, 3*nb_rows))\n", " patches = {}\n", " for col, (label, p) in enumerate(results):\n", " ax = axes[col]\n", " patches.update(make_pie(p, ax))\n", " ax.set_title(label)\n", "\n", " axes[1].legend(patches.values(), patches.keys(),\n", " bbox_to_anchor=(0.5, -0.2), # Legend position\n", " loc='lower center', \n", " ncol=2, \n", " fancybox=True)\n", "names = [0,1]\n", "colors = [\"#3f90da\", \"#bd1f01\"]\n", "make_plot([\n", " (\"uniform MP\", urw),\n", " (\"nexp MP\", edw),\n", " (\"fully-async\", faw)])" ] }, { "cell_type": "markdown", "id": "11483996", "metadata": {}, "source": [ "## Larger scale model\n", "\n", "We demonstrate how to import a model in GINsim format and perform simulations in different mutation conditions." ] }, { "cell_type": "code", "execution_count": 11, "id": "510d912d", "metadata": {}, "outputs": [], "source": [ "import ginsim" ] }, { "cell_type": "code", "execution_count": 12, "id": "b81594a9", "metadata": {}, "outputs": [ { "data": { "text/html": [ "Using local file SuppMat_Model_Master_Model.zginml
" ], "text/plain": [ "/home/pauleve/orga/projects/mpbn/code/examples/SuppMat_Model_Master_Model.zginml" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "" ], "text/plain": [ "" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "wt_model = ginsim.load(\"http://ginsim.org/sites/default/files/SuppMat_Model_Master_Model.zginml\")\n", "ginsim.show(wt_model)" ] }, { "cell_type": "markdown", "id": "373b8fb2", "metadata": {}, "source": [ "We convert the model to an `mpbn.MPBooleanNetwork` object and define its initial state according to the original publications (all nodes but microRNAs inactives, and DNA damage and ECMicroenv active." ] }, { "cell_type": "code", "execution_count": 13, "id": "d74e12bf", "metadata": {}, "outputs": [], "source": [ "f = mpbn.MPBooleanNetwork(wt_model)\n", "init_active = [\"miR200\", \"miR203\", \"miR34\", \"DNAdamage\", \"ECMicroenv\"]\n", "x = {node: node in init_active for node in f}\n", "nb_sims = int(1e3)" ] }, { "cell_type": "markdown", "id": "7f24322e", "metadata": {}, "source": [ "We compute the reachable attractors from the initial state. The first two attractors correspond to apoptosis while the latter to metastasis." ] }, { "cell_type": "code", "execution_count": 14, "id": "59a20afb", "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", " \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", "
AKT1AKT2ApoptosisCDH1CDH2CTNNB1CellCycleArrestDKK1DNAdamageECMicroenvEMTERKGFInvasionMetastasisMigrationNICDSMADSNAI1SNAI2TGFbetaTWIST1VIMZEB1ZEB2miR200miR203miR34p21p53p63p73
000110010110000000000100001001011
100110010110000000000100001101100
201001011111111111111111110000000
\n", "
" ], "text/plain": [ " AKT1 AKT2 Apoptosis CDH1 CDH2 CTNNB1 CellCycleArrest DKK1 \\\n", "0 0 0 1 1 0 0 1 0 \n", "1 0 0 1 1 0 0 1 0 \n", "2 0 1 0 0 1 0 1 1 \n", "\n", " DNAdamage ECMicroenv EMT ERK GF Invasion Metastasis Migration NICD \\\n", "0 1 1 0 0 0 0 0 0 0 \n", "1 1 1 0 0 0 0 0 0 0 \n", "2 1 1 1 1 1 1 1 1 1 \n", "\n", " SMAD SNAI1 SNAI2 TGFbeta TWIST1 VIM ZEB1 ZEB2 miR200 miR203 \\\n", "0 0 0 0 1 0 0 0 0 1 0 \n", "1 0 0 0 1 0 0 0 0 1 1 \n", "2 1 1 1 1 1 1 1 1 0 0 \n", "\n", " miR34 p21 p53 p63 p73 \n", "0 0 1 0 1 1 \n", "1 0 1 1 0 0 \n", "2 0 0 0 0 0 " ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "A = list(f.attractors(reachable_from=x))\n", "A.sort(key=lambda a: [a[i] for i in f]) # stable ordering, for notebook reproducibility only\n", "tabulate(A)" ] }, { "cell_type": "markdown", "id": "5945471b", "metadata": {}, "source": [ "We perform parallel simulations with the specified number of parallel jobs." ] }, { "cell_type": "code", "execution_count": 15, "id": "63a0a58b", "metadata": { "tags": [ "skip_test" ]}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "100%|█████████████████████████████████████████████████████████████████| 1000/1000 [00:25<00:00, 39.98it/s]\n" ] }, { "data": { "text/plain": [ "{0: 28.8, 1: 51.0, 2: 20.2}" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "mp_wt = mpsim.parallel_estimate_reachable_attractors_probabilities(f, x, A, nb_sims,\n", " depth = mpsim.nexponential_depth(f),\n", " W = mpsim.nexponential_rates(f), nb_jobs=8)\n", "mp_wt" ] }, { "cell_type": "markdown", "id": "572e6cfc", "metadata": {}, "source": [ "Now, we apply a mutation (p53 loss of function):" ] }, { "cell_type": "code", "execution_count": 16, "id": "10835137", "metadata": {}, "outputs": [], "source": [ "f_ko = f.copy()\n", "f_ko[\"p53\"] = 0" ] }, { "cell_type": "markdown", "id": "24677074", "metadata": {}, "source": [ "We compute the reachable attractors..." ] }, { "cell_type": "code", "execution_count": 17, "id": "ab9adbbb", "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", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
AKT1AKT2ApoptosisCDH1CDH2CTNNB1CellCycleArrestDKK1DNAdamageECMicroenvEMTERKGFInvasionMetastasisMigrationNICDSMADSNAI1SNAI2TGFbetaTWIST1VIMZEB1ZEB2miR200miR203miR34p21p53p63p73
000110010110000000000100001001011
101001011111111111111111110000000
\n", "
" ], "text/plain": [ " AKT1 AKT2 Apoptosis CDH1 CDH2 CTNNB1 CellCycleArrest DKK1 \\\n", "0 0 0 1 1 0 0 1 0 \n", "1 0 1 0 0 1 0 1 1 \n", "\n", " DNAdamage ECMicroenv EMT ERK GF Invasion Metastasis Migration NICD \\\n", "0 1 1 0 0 0 0 0 0 0 \n", "1 1 1 1 1 1 1 1 1 1 \n", "\n", " SMAD SNAI1 SNAI2 TGFbeta TWIST1 VIM ZEB1 ZEB2 miR200 miR203 \\\n", "0 0 0 0 1 0 0 0 0 1 0 \n", "1 1 1 1 1 1 1 1 1 0 0 \n", "\n", " miR34 p21 p53 p63 p73 \n", "0 0 1 0 1 1 \n", "1 0 0 0 0 0 " ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "A_ko = list(f_ko.attractors(reachable_from=x))\n", "A_ko.sort(key=lambda a: [a[i] for i in f]) # stable ordering, for notebook reproducibility only\n", "tabulate(A_ko)" ] }, { "cell_type": "markdown", "id": "7e83fe95", "metadata": {}, "source": [ "... in this example, they are a subset of the attractors reachable in the wild type model:" ] }, { "cell_type": "code", "execution_count": 18, "id": "0e773bf7", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[True, True]" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "[a in A for a in A_ko]" ] }, { "cell_type": "markdown", "id": "0100f64b", "metadata": {}, "source": [ "We perform simulations on the mutated model giving the wild type attractors for make comparison easier." ] }, { "cell_type": "code", "execution_count": 19, "id": "7d31ae26", "metadata": { "tags": [ "skip_test" ]}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "100%|█████████████████████████████████████████████████████████████████| 1000/1000 [00:19<00:00, 51.08it/s]\n" ] }, { "data": { "text/plain": [ "{0: 60.0, 1: 0.0, 2: 40.0}" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "mp_KO = mpsim.parallel_estimate_reachable_attractors_probabilities(f_ko, x, A, nb_sims,\n", " depth = mpsim.nexponential_depth(f),\n", " W = mpsim.nexponential_rates(f), nb_jobs=8)\n", "mp_KO" ] }, { "cell_type": "markdown", "id": "baed5840", "metadata": {}, "source": [ "We observe that the mutation substantially increases the probability of the metastasis attractor." ] }, { "cell_type": "code", "execution_count": null, "id": "90168ee2", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.0" } }, "nbformat": 4, "nbformat_minor": 5 }