diff --git a/figure_plots/20220218_new_shape_analysis_test.ipynb b/figure_plots/20220218_new_shape_analysis_test.ipynb index be85662e0f6256f7290abbc03b25410ecb4c26f3..59619945c36d9ca610eed6005dbe66560eb39455 100644 --- a/figure_plots/20220218_new_shape_analysis_test.ipynb +++ b/figure_plots/20220218_new_shape_analysis_test.ipynb @@ -2,7 +2,7 @@ "cells": [ { "cell_type": "code", - "execution_count": null, + "execution_count": 1, "metadata": { "collapsed": true }, @@ -18,6 +18,156 @@ "import warnings\n", "warnings.filterwarnings('ignore')" ] + }, + { + "cell_type": "code", + "execution_count": 2, + "outputs": [], + "source": [ + "# folder to save all panels for figure 2\n", + "# savefolder = r\"plots\\fig2\"\n", + "\n", + "# file containing the data for the controls\n", + "results_ctrl_file = r\"data\\shape_analysis\\histograms_HealthyControl_deformed_undeformed.txt\"" + ], + "metadata": { + "collapsed": false, + "pycharm": { + "name": "#%%\n" + } + } + }, + { + "cell_type": "code", + "execution_count": 24, + "outputs": [], + "source": [ + "def logistic_growth(x, A, x_0, k):\n", + " '''The logistic growth fct'''\n", + " return A/(1+np.exp(-k*(x-x_0)))\n", + "\n", + "def logistic_growth_simple(x, x_0, k):\n", + " '''A simplified logistic growth fct'''\n", + " return 1/(1+np.exp(-k*(x-x_0)))\n", + "\n", + "def hyperbola(x, a):\n", + " '''Hyperbolic function with growth rate a and limit 1 with\n", + " root at f(0)=0: f(x) = a/(x-a) + 1\n", + " '''\n", + " return a/(x-a) + 1\n", + "\n", + "def exponential_inverted(x, a):\n", + " return 1 - np.exp(-a*x)" + ], + "metadata": { + "collapsed": false, + "pycharm": { + "name": "#%%\n" + } + } + }, + { + "cell_type": "code", + "execution_count": 25, + "outputs": [ + { + "data": { + "text/plain": "<Figure size 360x252 with 1 Axes>", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWEAAAD1CAYAAACBf0Q2AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Z1A+gAAAACXBIWXMAAAsTAAALEwEAmpwYAABQ/0lEQVR4nO3dd3iUxdrA4d+7PZu2CZAAKdTQW5DeW0CIgDTFwrGggA2ODUH9sIJdDx49YjvYCyJIVUHwAAIiKBg6BEivpJfN1vn+WFgICSGB9Mx9XbnI7sz77jNZeJjMTlGEEAJJkiSpRqhqOgBJkqSGTCZhSZKkGiSTsCRJUg2SSViSJKkGySTcQLVv357Y2NirunbEiBHs2rWr1LJ9+/YxZsyYUusuW7aMp5566qpeszRvvfUWffv2ZeDAgZV2z4vNmDGD7777rlx1//zzT0aPHk14eDi//PJLlcRztVatWsUtt9xSqffcs2cPQ4YMcT8u6++EVDZNTQcgld+IESM4e/YsarUaDw8PhgwZwv/93//h6elZ06G59erVi59//rnUsjlz5ri/T0hIYOTIkRw+fBiNpuJ/DZOSkli+fDm//vorjRo1uup4K8vbb7/Nbbfdxh133FHToUh1jOwJ1zHLli1j//79rF69mkOHDvHee++VqGO322sgsuqVlJSEyWS6qgRcFT+fpKQkwsLCrurahvB+SZcnk3AdFRgYyODBgzl58iTgGl748ssvGT16NKNHjwZgxYoVRERE0KdPH+bMmUNqamqxe2zbto2RI0fSt29fXnnlFZxOJwBxcXH84x//oG/fvvTt25dHH32U3NzcYtcePHiQcePG0bt3bxYuXIjFYgFK/pp6sX//+9889thjANx+++0A9O7dm/DwcP744w/69OnD8ePH3fUzMjLo3r07mZmZxe6za9cu7r77btLS0ggPD2fBggUAbNmyhcjISHr16sWMGTM4deqU+5oRI0bwwQcfMH78eHr06FFq4tu5cyfXX3891113Hc8//zyXTqFfuXIlY8eOpXfv3sycOZPExEQARo0aRXx8PHPmzCE8PByr1Upqaipz5syhT58+REREsGLFimI/h7lz5/LYY4/Rs2dPVq9ezYwZM3jrrbeYPn064eHhzJkzh6ysLB599FF69uzJlClTSEhIcN/j1KlT3HXXXfTp04cxY8awceNGd1lWVhZz5syhZ8+eTJ06lbi4uFLfj/P27dvH9OnT6dWrF0OHDmXVqlUAWK1WXnnlFYYNG8aAAQNYtGgRRUVFZd5LugpCqjOGDx8udu7cKYQQIikpSYwbN0689dZbQggh2rVrJ+68806RlZUlzGaz2LVrl+jTp484dOiQsFgs4vnnnxe33nqr+17t2rUTt99+u8jKyhKJiYli9OjRYsWKFUIIIWJiYsRvv/0mLBaLyMjIELfeeqt48cUXi8URGRkpkpKSRFZWlrj55pvFm2++KYQQ4vfffxeDBw8uNea3335bPProo0IIIeLj40W7du2EzWZz133mmWfEq6++6n78ySefiNmzZ5f6s7j0dU6fPi26d+8ufvvtN2G1WsUHH3wgRo0aJSwWizuOCRMmiKSkJGE2m0vcLyMjQ/To0UP8+OOPwmq1iuXLl4uOHTu6fyabN28Wo0aNEtHR0cJms4l3331X3HzzzaW2Uwghbr31VvHMM8+IoqIiceTIEdG3b1+xa9cu98+hU6dOYvPmzcLhcAiz2Sxuv/12MWrUKBEbGytyc3PF2LFjxejRo8XOnTuFzWYTjz/+uFiwYIEQQoiCggIxZMgQsXLlSmGz2cThw4dFnz59xMmTJ4UQQvzzn/8Uc+fOFQUFBeL48eNi0KBBYvr06aX+HBMSEkSPHj3EunXrhNVqFZmZmeLIkSNCCCEWL14sZs+eLbKyskReXp6YPXu2eP3116/4PksVI3vCdcwDDzxAr169uPXWW+ndu3excdZZs2ZhMpkwGAysW7eOKVOm0LlzZ3Q6HY888ggHDhwo1pu69957MZlMNG/enH/84x+sX78egBYtWjBw4EB0Oh3+/v7cdddd7N27t1gct912G82aNcNkMnHfffexYcOGa27bpEmT2LBhg7sHumbNGiZMmFCuazdu3MjQoUMZOHAgWq2WmTNnUlRUxP79+911ZsyYQbNmzTAYDCWu3759O2FhYVx//fVotVruuOMOGjdu7C7/5ptvmDVrFm3atEGj0TBnzhyOHj3q7g1fLDk5mb/++ovHHnsMvV5Px44dmTZtGmvWrHHX6dGjB6NGjUKlUrnjmTx5MqGhoXh7ezNkyBBCQkIYMGAAGo2G66+/niNHjgDwv//9j6CgIKZMmYJGo6FTp06MGTOGn376CYfDwaZNm5g7dy5Go5F27doxadKky/7c1q9fz4ABA7jhhhvQarX4+fnRsWNHhBCsWLGCJ598EpPJhJeXF7Nnz66U91kqTn4wV8e8++67DBgwoNSyZs2aub9PS0ujc+fO7seenp6YTCZSU1MJDg4uUT8oKIi0tDQAzp49y+LFi9m3bx8FBQUIIfDx8bnsazVv3tx97bXo3r07BoOBPXv20KRJE+Li4hg5cmS5rk1LS6N58+buxyqVimbNmhUbgrk45tKub9q0qfuxoijF6iclJbFkyRJeeeUV93NCCFJTUwkKCipxL19fX7y8vNzPNW/enEOHDrkfX/xa512c9PV6fbHHBoOBwsJCABITE4mKiqJXr17ucofDwYQJE8jMzMRut5d4fy4nOTmZ0NDQEs9nZmZiNpuZPHlysfaeH7KSKo9MwvWIoiju7wMCAor10goLC8nOziYwMND9XHJysvvDpKSkJAICAgB48803URSFdevWYTKZ+OWXX3j++eeLvVZycrL7+4uvvZpYLzZp0iTWrl1LkyZNGDNmDHq9vlz3CwgI4MSJE+7HQgiSk5OLtfdyrwnQpEkTUlJSSlx/XrNmzZgzZ065euYBAQHk5OSQn5/vTsQVieVKmjVrRu/evVm+fHmJMofDgUajITk5mTZt2rhfu6x7RUVFlXjez88Pg8HAhg0bisUtVT45HFFP3XDDDaxatYqjR49itVp588036datm7sXDPDxxx+Tk5NDcnIyn332GePGjQOgoKAAo9GIt7c3qampfPTRRyXu/9VXX5GSkkJ2djbLli1zX1te/v7+qFQq4uPjiz0/YcIEfvnlF9auXcuNN95Y7vuNHTuWbdu2sXv3bmw2G//973/R6XSEh4eX6/qhQ4dy8uRJNm3ahN1u57PPPuPs2bPu8unTp/PBBx+4PwjNy8vjxx9/LPVezZo1Izw8nDfffBOLxcKxY8dYuXJluYdWrmTYsGHExMTwww8/YLPZsNlsREVFcerUKdRqNREREbzzzjuYzWaio6NZvXr1Ze81fvx4du3axcaNG7Hb7WRlZXH06FFUKhXTpk1jyZIlZGRkAJCamsqOHTsqpQ3SBTIJ11MDBgxg3rx5PPTQQwwaNIj4+HjeeuutYnVGjhzJ5MmTufHGGxk2bBhTp04F4MEHH+TIkSP06tWLWbNmuWdbXOyGG27g7rvvZtSoUYSGhnLfffdVKD4PDw/mzJnDLbfcQq9evThw4ADgSmCdOnVCUZRiv25fSevWrXnttdd44YUX6NevH7/++ivLli1Dp9OV63p/f3+WLl3KG2+8Qd++fYmNjaVnz57u8oiICO655x4eeeQRevbsyQ033MD27dsve78333yTxMREBg8ezIMPPshDDz102WGkivLy8uLjjz9m48aNDB48mEGDBvH6669jtVoBWLRoEYWFhQwcOJAFCxYUG1K4VPPmzfnwww9Zvnw5ffr04cYbb+TYsWMAPP7447Ro0YKbbrqJnj17cuedd3LmzJlKaYN0gSKE3MpSql0WLlxIQEAADz/8cE2HIklVTo4JS7VKQkICmzdvLvNXaEmqT+RwhFRr/Otf/2L8+PHMnDmTkJCQmg5HkqqFHI6QJEmqQbInLEmSVIPqzJiwEAK7ve5PFFerFRyOuv3LR31oA9SPdtSHNkD9aYdWq67wNXUoCUN2dmFNh3HNTCZjnW9HfWgD1I921Ic2QP1pR5Mm3hW+Rg5HSJIk1SCZhCVJkmqQTMKSJEk1SCZhSZKkGlTpSXjhwoX079+fG264odRyIQQvvvgiERERjB8/nsOHD1d2CJIkSXVGpSfhyZMnl7rr1nnbt28nJiaGTZs28cILL/Dss89WdgiSJEl1RqVPUevdu3ex0xsutWXLFm688UYURaFHjx7k5uaSlpZ2xf1oz5wt4PEVfxd7blT7Jkzr0Zwim4N5qw6VuOaGzoGM79KU7EIbT6w7UqJ8SvdmjO4QQEpuEc/8eLxE+W29ghnSphExmYW8tPlkifK7+4XSt4Ufx9PyefPXUyXK7x/Uku5BvvydmMN/fosBQKNRuec7PzK8De0DvNgTm8V/fy95DtjCiDBa+hvZfiqDL/eV/Jk+N7Y9TX0MbDqWxvd/l9wz9pXxnTAZtaw7lML6w6klypdO7oJBq+a7A0n8cjy9RPn7N3cH4PO98fx2+sI5bxqNCjXw9pSuAHy0O5a9cdnFrvX10PLqhE4AvLPjDAeTip9RF+Ct54VxHQB449dTnEjLL1Ye6ufBU6PbAbB40wnisszFytsFePHocNd+uf+38RhpeZZi5V2b+/Dg4FYAzF97hByzrVh571ATj43tCMDc7w9iuWQO+qDW/szo7Vo6Pfvb4n/vAEa1a8y0zk2wFBTy1Kq/0dptaO02NHYbGoeNYS18GBLqQ35+ER/tiEbjcKBy2lE7HKidDvqH+NAt0JOcfDNrDiShcjpQCScqp+urd7APbfwMZBUU8b/j6SjCiSIEKuFEcQp6BvsQ5K0j0+pgx4m0c+WgIFCEk/DmPjTx1JGeV0RUYg4KgBDnygU9gn3xM2hIzS3iaGo+invhrEAREB7si7deQ3KOmeizBa5rL1pce12wLx46NYnZZmIySk4t6x1qQqdWiM8yE59d/L1ThKBvCz/UKoUzmYUk5xShqBSE48J7MKCVPwCnzhaQesl7i3AS5u3E0+hJYoHgbIG1WLFOreK6EF8AjqXlk1VY/L3Xa1V0CzLhUBQOpRaQY3HgUFQ4FQWhqDDo1NjPxmG121AC2lAkVAiUc+UKRr2GVo29cKJwLL2AIodAoIAC924qucfzlVT7POHU1NRipwo0bdqU1NTUK28Krrj+8V/M6KHDZDJitjpKlAEYja5yp9Zaarmnp951vaKUWe5jc5Za7uXlKvc220st9/Y2uMpzLe5y5aLXOl/ulVFY6vU+58o9PfNLLff19cDk64Gnp7706309MHnqMBp1l7neiIdOjdGj9HKTyQiAxyXliqKg1ajc5QaDtsT1Wq36QrleU6Jcd1G5vpRyvV5T7nKdVl2i3HBRuVarRmNzoDideBXm4pWfQyt7AuqfY/HLyOT636LQ5udhNBfgUVSIwVJIiNZBI4OAggLejE9HZ7Wgs1nQW4vQ2qzo7FZU506Z+LTET+4CP+C5Msp9gAfLKPcCZpRR7glMLaO8xbmvy2l57utyWp/7upQA7Co1TTU6/NRaLBqd6+vc9yfitVg0WqxqLTaNFotai02tcT1Wazl9phCrSoNNrcem88Sm1mJVa7CpNNjVGtbnabCr1di9NNh81NjVGuwqNTaVBodKjU2txqGosTdW41CpsatU5/4891hR41CpsDdzPXYoKpwqFXaVGqFc9HelrMaXR/FDVbj3Km5RJXtHJCQkMGfOHPeZZRebPXs29957r3uv2DvuuIPHHnuMrl27lnlPp1OQkZFfZp26oD5MSq91bRACJT0ddXwsqsQE1MlJqJKTUaUko0pPQ5WWiupsOkpWFkolHs8jNBrQ6RF6HUKnB70eodWCTofQ6kCrAY3W9ZxG46qvOfe9Wg0XfQm1mr8PRuFUFMJ79XaVKypXuUoFiuJ6TqUAius5lQqDh44iq8NVriiuaxTF9RgFs6ImHzUFipoC1Gz9fQ95TmjfoxcYvSkUKgpRYRZqzKgwC4UiVJiFiiJUFAmFIqGiCIUioWA5971FKDi5+tNBapoinKiEQMX53zDOfw84HThsFhRAq1Zh1OvRa7WoFNf4rYJA7f4eVOce61QK3z17+fP8Lqfae8KBgYHFjpFJSUmRx6dI5aJkZ6E+fhzNiWOoT55AHXMa9ZnTqGNjUMp5FLvT14TT3x/h54e6cSOsnt4IHxPCxwfh7Y3T2xvh6XXuyxNh9ASjB8LoiTAYEAYPhMEAHh6uBFmJvn91CQBh858EXB9iF9md5Jht5BTZySuyk1vk+j7fYie3yI4VyMgtIt9iJ8dsI6/IRoHVToHVSZFdUKKH1b3the/NlxZWjEoBvUaFTq1Cr3F9adUqss6mokLQKjQEnUaFRqVCjRMcdpx2Kw6bBbvVgt1ahLXIjMNmQa9RcNjtqDiXDBHu5OZl9MDbyxNfLy+cDjtffv4JKpy0btWam6bdRECTxijCVVfBCUKgUlyJ9vx9EIBwolIEwulECHHuy4nz3GOn0/V9WloaL7zgOs7ruv4DeOyxhQQFBZf9w7gG1Z6ER4wYwRdffEFkZCR///033t7eFT6fTKpbXj2XXOafSy7loaSlof1rH5oDf6I5dBDNoYOok0qebHye088PR3AozqBgnM2b42jWHGdgU9dXkwCcTQIQ/v6g1bqvMZmM5FVzj95sc5BRYCWjwEpWoY3MQiuZhTayzTZ2OFpjdqjY+/HvFNgh22zDdo37KaiFAw0OtDgQVjO5WemoHFY0OPH39cLX6IEGVx01zgt/ivOPnajPPadVgUGrwdOgxUOnxdPDgE6nQ683oNfr0On06PUGPv74e3JyMunXZDyKoiI7O4uiMv6T1Gg0NGsagF7vicnkh5+fHyaTHyaTCR8fX9SX/Gc3aUBX4uPjCAkJrbLkuHv3TnJzc6o8AUMVJOFHHnmEP/74g6ysLIYMGcJDDz2E3W4H4JZbbmHo0KFs27aNiIgIPDw8WLJkSWWHINU1QqCKjUG36ze0v21H+8fvqONiS1YzGrGHtcfRrj2OsHbY27TF0bI1zlatEF4VX7NfmexOwdl8C6l5FtLyraTnW0jLs3K2wEJ6vpWzBVbO5lsptDkufxOd61TknOwLHzSphAM9dnTnvrTY0YkL3xs1AsVuQSvseOpUeBu0+Hjo8TXq8fX0wNvTiIeHN0ajkby8HJ544iWcTicdOnTiHzfejZ+fiaKiIqxWKxZLERbL+T8tFz3nxGKxuf4dW6DQAoVARilNyM3NZdeu7QghSE8/y4ABg/Dx8UGn0xVLsBcnWm9vH/z8PMs9xBUUFFzlifHFF1+u0vtfrNKT8JtvvllmuaIoPPPMM5X9slJdYzaj+20bui2b0W3ZjDo2plixMHpi63kd9vDrsHfrjr1LVxwtW1f6EEB5vfzqEszouWH63SRmF5GcW0RiThGpeRaScy2czbdQnk6rGicGbOiFFT1WDMKGHhtF2ekcj9qLxmHBUwOd2rSgXctQTF4GPDyMGI3Gc3/64OHhgdHoiYeHB82aNcZmAw8PIxrNlf85DxkynNzcHJ544qkKJzKHw3EuYVuwWC79ciXyAwf2o9cb8Pb2oXnzIPr27cfw4aPw9PS8phOm67M6s4uaVA+Yzeh+2YR+/Q/oNv2MquDCB61OPz9sAwZjHTgIW7+BODp2qvaE6xSClFwLsVmFxGW6plbFZZlJyDaToBqAUFRs/O7gZa/3UjvxVNkxOM1o7QV4CCseWDAIKx5YMWBFiwMF0Ol0+PqaMJlM+PoGYrP58cK2L1EUhe7dw3nikXmEhIReMeaKfkh6LT08tVqNp6cnnp6el63TokVLVq78FqfTSZMmTejZsxdeXl5X/ZoNgUzCUtUSAv/ok/Q++Df+y95FnZ/nLrJ17Y519PVYR0ZgD7+u2pKuwylIyiliT0IOB2LOcvpsAbFZRcTnWLBerjurqBCFWeht+QT5GvBSrOgdBXgKC0YseGBBfdG1iqLg4+NzLtEGnUu2JvefHh4eJXqGO3Zsv+peam0RFBTMgAGDzo2nLqiz7ahOMglLVaOoCP2q77C/9RoLLxpqiAsI4GjXHhzr0pW8xo1Rq9Uofx9AfeggarUatVqNSqVCpVK5H6vVahTl/GPVJY/VLF/+IQUFhUyffis+Pr7Y7TasVht2u438IhvJRSpSitSkWTWcdejJcupxUHrCNwgr3hTiLcx4UYSXMFOQeobdm9eg16jR+vjStFt3goNDMBgM5xJrs2IJ9nIfKF2Jj48PPj4+dT5xVed4an0gk7BUqZScbDw+XIbhw2Wos1yr7NIVhdWenmwIaIqxT1+Cg88d4pmTc82vl5uby549uxFCkHr2LJ0HXo/FqylZeJOpeJGHh2ve7CU8hAV/lRkfivBXW2mss9NY58TboEGj0aDT6dBoTGi1TcjN9eXgbz8jBLRq1ZqZM2fRrl0HPDw8rjn+i1Vk9ohUf8gkLFUKJTsLj2Xv4vHhMlR5riXKCYGB/K9bOC+djga9nm7dwpk37xGaNm2G0+nA4XDgcDjPzc88/9iB0+m8zJ8XlzvJLrLz66EzePQ3oApoQ74piD2q4n+l1QqE+upo42+gbWMj7QI8aR/gg7+3B40aeZOTU77JslFRf9f5oQKpdpJJWLo2NhuGTz/G87WXUGVlAXAytAW/9BuAavhIRo8Zx6CXXzyXwJ68pgSWmmfhz/hs/orP4a+EbOKzHUAodAjFCSAEIT5auof406mpN52aetO2sSf6UpZkAxX6tL6+DBVItY9MwtJV027djNdTT6A5FQ3AmRYt2dh/IMlt2jJy5Gg6d+7i/oDqahJYvsXO3rhs/ojN4o+47BKb+HhoVXRt5kMrb2gkchnYIZR2ra48o+BqyKECqarIJCxVXHo63g/NxbDqOwByAgNZ1X8QR9q0pWWr1tx1/Th8fHwrfFshBNFnC/jtdCa7z2QSlZRbbO6tp05NeLAvPYN9uS7ERLsALzQqOfdUqttkEpYqRL96JZonH0ebkYFDr2fr4GFs6doNlcFAxLAR9OjRs0K/5tsdTv5MyGFbdAY7TmWQctG2hWoFwoN86NPCj74t/OjY1FsmXanekUlYKhclPw+vBY9hWPE1AEmdu/BZ/4FkmPwIDg5h7NhI/Pz8S7320l/lbQ4nv8dkseXkWXacyiC3yO4u8zdqGdy6EQNa+9Mn1ISXXv4Vleo3+TdcuiJN1AG8770TzZnTOPQGtkRGsrlVWzRaLcMGDaVXr96oVGUf0uJwCv6Mz+bnY2n8ejKDPMuFxNuqkZHhbRsxpE0jOjb1RiWXt0oNiEzCUpn036/A++EHUYqKyAptwUcjIsgNbkZTU2PGjRtP48aNy7z+dEYB6w+l8vOxNNLyL2xM07axJxHtmzAirDEtGxmruhmSVGvJJCyVzuHAc/FzGN/5FwD7e/VmxaChOHU6Rg8bRufOPS+7IqzAamfzsXTWHkrhYPKFZcpBvgau7xjAmA4BtJKJV5IAmYSl0pjN+My+C/1PG3GqVKwZPopdPcJp3CSAyMjxtG/futRNY06m5/P938n8eCTNvWWjp07N6A5NiOwUSLfmPnInLUm6hEzCUjFKTjY+M6aj+30XZg8PPh1/I6dbtKRvn34MHDi4xHaJDqdg+6kMvv4rkf0JF5Yhhwf5MLFrM0a2a4xBWzPbT0pSXSCTsOSmSk3B5+ZJaI8cJtvbmw+n3Iy9XTtuHXtDiYUWhVYHaw6l8M2fCSTluqaVeerUjOsUyOTuzWjb+PLbHUqSdIFMwg3c+aOHFtxxN143jEYbG0Oavz8fTr2ZNsNHMmTIcLQXHQmUVWjlkz8T+fz3WPfUsmCTgenhQdzQJRBPnfwrJUkVIf/FSHjm56MfOxJ9QjyJAQF8e/e9jJlyMy1atHTXySiw8sW+BFYeSKLI7jqxuFtzH/7RO5jBbRrJaWWSdJVkEm7gCuPimL1uNV5mM8mNm7DzmRe5acIk9Ho9ANmFNj7dG893B5KwnEu+Q9s14fbw5vQIrvjSZEmSipNJuAHbs2UT965aQVu7ndMGA8df/xfDx40HXNPMvvozkS/3JVBgdc10GNKmEff0D6V/+8AKHakjSdLlySTcQP31x+80f3QePex24rRaloyIYKJ/I+xOwdqDyby/K5bMQhsA/Vv6cd+glnQMrNkTjSWpPpJJuIFxOp1s3bKJlq8uoWdSIlkqFXNCW+Lv58dZTWNu/exPzmS4erldm3nz4JBW9Aw21WzQklSPySTcgFgsFtat+4HQr7+k/98HcGp1fDg2Er13E+z9pvPcthQAmvsaeGhwK0a2aywXV0hSFZNJuIHIzc3h+++/w3/Pbsbu2IZQFDLf+whnQDhH9sRjSSjEQ6tiZr8W3NIzCN1lTqOQJKlyySTcAKSkJLNq1Ur0CXHctnE9KgQ7Hn+JRZnBnDkZC8CYDk2YO6Q1Ad76Go5WkhoWmYTruZMnT7B+/RooKOCBjetxOgRPznyJr5SukFFIqJ8HC0a1pXeoX02HKkkNkkzC9ZQQgr17/2Dbtq0Ip5M5e/dwUu3P/FmLifdqjFqlcEefEO7uG3rZgzAlSap6VZKEt2/fzuLFi3E6nUybNo1Zs2YVK09KSuKJJ54gLy8Ph8PBY489xtChQ6silAbJ4XCwZcsmDhzYD8BUh+BLj47899YbAWjXxJNnx7YnrIlXDUYpSRJUQRJ2OBw8//zzLF++nMDAQKZOncqIESNo27atu857773H2LFjufXWW4mOjmbWrFls3bq1skNpkIqKili7djUxMWfQaDRcF9aDRbtSOd57CBoEd/Vvwd19Q9GoZe9XkmqDSk/CUVFRtGjRgpCQEAAiIyPZsmVLsSSsKAr5+fkA5OXlERAQUNlhNEg5OdmsXLmCjIyzeHh4Yuw0nKf+TMfSuAUtrDk8d+dQOjfzqekwJUm6SKUn4dTUVJo2bep+HBgYSFRUVLE6Dz74IDNnzuSLL77AbDazfPnyK95XUcBkqvunMajVqippR0JCAitXfkNBQQGBQUGc8O3FpgMZoNYy7fRunnprLp4BjSrltaqqDdWtPrSjPrQB6k87rkaNfDC3YcMGJk2axN13383+/fuZP38+69evL/OwSCGoF/sVmEzGSm/HsWNH2bhxHXa7HWPztqwtaklcSgZGq5klm95l6CtPYdN5VNrrVkUbakJ9aEd9aAPUn3Y0aVLxpf2VPjAYGBhISkqK+3FqaiqBgYHF6qxcuZKxY8cCEB4ejsViISsrq7JDqfeEEPz++y7Wrl2N3W5H07IXKzKbE5ddRFhuMms/fZjRw3tg79O3pkOVJOkyKj0Jd+3alZiYGOLj47FarWzYsIERI0YUq9OsWTN2794NwKlTp7BYLPj7+1d2KPWaw+Hgp582sn37/0BRyG85lK/jDZhtTm6wp7D2wwdp6e9BwfwnazpUSZLKUOnDERqNhkWLFnHPPffgcDiYMmUKYWFhLF26lC5dujBy5EgWLFjA008/zSeffIKiKLz88styj4JSnD/1Yv4lidRsNrNmzSri4mJRNHpOBw5iX7wFlQLz2uh4aM5sFKeT7H/9BwyGmghdkqRyqpIx4aFDh5aY9ztv3jz3923btuWbb76pipeu97KyMlm16jsyMjJQjCb+8ryO6BQL3noNS8a2Y+zMCagcDgrnPCiHISSpDpAr5uqQhIR4Vq/+HrO5EJVfML86wkjLtBDka+Bfk7rQ8bv/ojlyCEdoCwoWPF3T4UqSVA4yCdcRR44c5scf1+NwONA068C67EDyrTa6NvPhjRs70Sg7HeMrruGL/JdeA2PDnO4jSXWNTMK1nBCCXbt+47fftgNgaN2bFQkeWBwOhoc15oVxHdBrVHjOfRJVQT6WceOxRlxfw1FLklReMgnXYjk52Rw5cpi8vDx8fX0xdBzGlydsOISTyd2aMX9kW9QqBe3WXzCsXY0wGsl/8eWaDluSpAqQSbiWSkxM4KefNlJUZCY/P4+I+17ks+OuM99m9gtl9oAWrhklViteTz4OQMFjC3EGh9Rk2JIkVZBMwrXU4cMHKSoyYzB4YOwRyXdnXM/PG9qa23sFu+t5/PcDNKdPYW8bhnn2/TUUrSRJV0sm4VrKbDa7erodRpLUtD8Aj49ow03hQe46SkYGxjdeBaDgucWg1dZIrJIkXT2ZhGshIQTp6Wk06j+Vwg7jUIAFEWFM7tasWD3P15agysnGOmwE1lFjaiZYSZKuidxUthZKS0tlX5aWwg7jgNITsPr4MQyf/hehUpH/3BLXNnOSJNU5MgnXQp/vOMKfimv/5R7OUyUSMIDns0+hOBwU/eMuHB07VXeIkiRVEpmEa5ntp86yIlYNikK7wmOEieQSdbQ7d6DfshmnlzcF85+qgSglSaoscky4Fvk7MYeF644iUOiuTSNMlw5cMswgBJ4vPgOA+YG5iMaNqz9QSZIqjUzCtcSpswU88sNhrA5BK2cKt3drwrBhJXu5uo3r0f65D2fjJhTOfqAGIpUkqTLJ4Yha4Gy+hXmrDpFbZCdYyaaniKZTp84lK9rteL70PAAFj84HL3lasiTVdTIJ1zCzzcEjPxwmNc9CmL+WPvYjBDRuUurhp4YVX6M5cRxHi5YUzbirBqKVJKmyySRcg5xCsGjjMY6m5tPc18BEv7OocdKxtNkOFgvG1137QhQseBp0umqOVpKkqiCTcA16Z/sZ/hedgZdezevj25Mccxyg1CRs+Opz1Anx2Dt2wjJpanWHKklSFZFJuIb8dDSNz/cloFYpvDK+E/bMRKxWK82bB2Ey+RWvXFSE8V+vA65NeijjVGpJkuoW+a+5BhxPzefFTScAeGRYG/q08OPYsSPAZXrBX36KOjkJe+euWCPHV2uskiRVrTKT8NatW3E6ndUVS4OQUWDlsTWHsdidTOgSyLQezTCbzZw+fQpFUWjfvmPxC8xmjP96A8B1crLsBUtSvVLmv+iNGzcyevRoXn31VU6dOlVdMdVbDqfg4RUHSMmz0KWZN/NHhqEoCidPHsfhcBAa2gKvS6adeXz2X9SpKdi69cB6/bgailySpKpS5mKN119/nfz8fNavX8/ChQtRFIXJkycTGRlZIllIV/bR7lh2n87E36jllfGd0Gtc/wceOXIYoOTcYLMZj3//C4DCxxfKTXokqR664u+2Xl5ejBkzhnHjxpGens7mzZuZPHkyn3/+eXXEV2/8HpPJx7/HoSjwYmQHArz1AOTn5xEfH4darSYsrH2xawxffY46LdXVCx4tz42TpPqozJ7wL7/8wurVq4mLi2PixIl89913NGrUCLPZTGRkJDNmzKiuOOu0tDwLizYeRwBzh7eld+iF2Q/Hjh1FCEHr1m0wGAwXLrJaMb7zLwAKH35c9oIlqZ4qMwlv3ryZO++8k969exd73sPDg8WLF1dpYPWFwyl4euMxssw2+oSauH9oG/Jyze7yo0fPz4ooPhRh+O4b1IkJ2Dt0xDo2slpjliSp+pQ5HNG4ceMSCfi1114DoH///lUXVT3y2d549ifk0MhTxwuRHVCrLvRos7IySU5OQqfT0aZN2wsX2e0Yl7pmRBTOe1TOiJCkeqzMf927du0q8dz27durLJj65nBKHu/vigXgmevb4W8svtT42LGjALRt2w7tRefD6X/4HnXMGeytWmOZOLn6ApYkqdqVOhzx1Vdf8fXXXxMXF8f48RcWBxQUFNCzZ88r3nT79u0sXrwYp9PJtGnTmDVrVok6Gzdu5J133kFRFDp06MAbb7xxDc2ofQqtDhZtPIbDKbilZxD9W/oXKxdCXDQr4qIFGk6nuxdsnvcoaORuo5JUn5X6L3z8+PEMGTKEN998k0cffdT9vKenJyaTqcwbOhwOnn/+eZYvX05gYCBTp05lxIgRtG174dftmJgYPvjgA77++mt8fX3JyMionNbUIm/+7xRxWWbaNvbkgcGtSpSnpaWRkXEWDw8jLVpcKNdt+gnN8WM4goIpmnpzdYYsSVINKHU4QlEUgoODWbRoEZ6enu4vgOzs7DJvGBUVRYsWLQgJCUGn0xEZGcmWLVuK1VmxYgW33XYbvr6+ADRq1KgSmlJ77DqTyZqDKWjVCi9EdnDPB77Y+WXK7du3R61Wu54U4kIv+L4H5U5pktQAlNoTfvTRR3n//feZPHkyiqIghHCXKYpSIqleLDU1laZNm7ofBwYGEhUVVaxOTEwMANOnT8fpdPLggw8yZMiQMgNVFDCZjFdsUE3LK7Lx0i/RAMwbEUavtk2KlavVKnx9PYiNjcZg0NK3by93u5Qd29H8uRfh74/+gfvQe9bO9qrVqjrxXlxJfWhHfWgD1J92XI1Sk/D7778PuPaOqAoOh4PY2Fg+//xzUlJSuP3221m3bh0+Pj6XvUYIyM4urJJ4KtMLPx8nJbeILs28mdIlsETMJpORw4dPkJKSjre3Dz4+jd11fJa8hAYonDmbQpsCtbS9JpOxTrwXV1If2lEf2gD1px1NmnhX+JpSk/Dhw4fLvKhz51KO3jknMDCQlJQU9+PU1FQCAwNL1OnevTtarZaQkBBatmxJTEwM3bp1q0jstc7OM5msPZSKTq3wzJj2aFSlL7A4PxTRoUNHlHOLMNQHo9Bv2YwwemKeWfKDTEmS6qdSk/DLL7982QsUReGzzz67bHnXrl2JiYkhPj6ewMBANmzYUGLmw6hRo9iwYQNTpkwhMzOTmJgYQkJCrrIJtUOB1c5Lm08CMGdgS1o2Kv1XK6fTybFjx4Die0UY33kLAPOMOxH+9WuMXJKkyys1CV/LvhAajYZFixZxzz334HA4mDJlCmFhYSxdupQuXbowcuRIBg8ezM6dOxk3bhxqtZr58+fj5+d35ZvXYst2xpKaZ6FjoBe3Xhd82XqnT5+msLCARo0aERDg+g1BFRuDfs1qhEaDeY48QVmSGpJSk/Du3bvp378/mzZtKvWi0aNHl3nToUOHMnTo0GLPzZs3z/29oigsXLiQhQsXVjTeWulISh4r9ieiVuCpiHbFVsVd6vxQT8eOnd1DER7vv4vidFI09WacQZdP4JIk1T+lJuG9e/fSv39/fv3111IvulISbkjsTsGSzSdxCrjtumDaB15+i0+bzeZeJdehg2vzdiUzA4+vXL95FN4/t+oDliSpVik1Cc+d60oGL730UrUGUxd9+1cix9PyaeqtZ/bAFmXWPX36FBaLhaZNm+F/btzXY/lHKIWFWEZG4Lh0P2FJkuq9MtfEZmVl8e677/Lnn3+iKAo9e/bkgQceqPPjt5UlLc/CB+f2hnhiVFs8tOoy6x89en4o4twyZbMZj49d0wHND8y73GWSJNVjZW7g88gjj+Dn58fbb7/N0qVL8ff35+GHH66u2Gq9t7efptDmYFjbRgxqXfaMhqKiIvc5cueHIgzffoXq7FlsPcKxDRxcHSFLklTLlJmE09PTeeCBBwgJCSEkJIT777+/Xu7zcDX2J+Tw87F09BoV/xzW+or1T548gd1up0WLFnh7+4DDgcd7/wbO9YLlpu2S1CCVmYQHDhzIhg0bcDqdOJ1ONm7cyKBBg6ortlrL4RS8ttW1NHlGr2CCfD2ueM35oYguXboAoPtpI5ozp3GEtsQSOaHqgpUkqVYrdUw4PDzcvWfEp59+yuOPPw64FhoYjUaeeOKJag2ytlkVlczJ9AKaeuu5o8+VF5nk5+cTFxeLSqWiQ4eOWK1g/M/bABTOuV9uVylJDVip//r3799f3XHUGblFNt7fGQPAw8PbYLjCh3EAJ04cw+l00qZNW4xGI87ffkW7dw9Ok4mi6bdXccSSJNVmV+yC5eTkEBsbi8VicT936ZFHDcnyPfHkFNnpGezL8LblW1586TlyxvNjwXfeA16Xn1csSVL9V2YS/u677/jss89ISUmhQ4cO/P333/To0aPMvSPqs8QcM9/uTwTgn8Nau1e8lSUnJ5vExAS0Wi1t24ZBdDS6jesQOh1FcqMeSWrwyvxg7rPPPmPlypU0b96czz//nNWrV5e53WR9958dMdgcgrEdA+gYWL4t644ePX+OXBg6nQ7V20tRhHAtUQ5seoWrJUmq78pMwjqdDr1eD4DVaqVNmzacOXOmWgKrbQ4l57LpeDo6tcL9g1qW+7oLCzQ6o2RkoPr0EwDMcx6sgiglSapryhyOaNq0Kbm5uYwaNYq77roLHx8fmjdvXl2x1RpCCN7edhqAW64LpqmPoVzXpaenk56ehsHgQatWrTE93h3FbMYyajSOcws2JElq2MpMwu+++y4ADz30EH379iUvL4/Bgxveyq5dMVnsT8zF16DhznJMSTvv/Obt7dq1R22zoV6XDIBZbtQjSdI5V5wdcfjw4WJ7R+ga2OGTTiF477cYAO7oE4KXvnxzeoUQxfaKMHz3DUqODdHWSy5RliTJrcwx4XfeeYcFCxaQnZ1NVlYWCxcu5D//+U91xVYr/HryLMfT8mnipWNaj/IPxSQnJ5GdnY2XlzchQcHuJcpi5hi5RFmSJLcyu3Xr1q1j7dq17g/nZs2axcSJE7n//vurJbia5nAKlp1bmDGzX2i5Fmacd74X3KFDRwxbNqOJPokjOATn3C+gwFYV4UqSVAeV2RMOCAgotkjDarWWOLSzPvvxaCoxmWaa+xqY0KX808kuPkeuY8dOeJxbomyedR9otVUSqyRJdVOpPeEXXngBRVHw9vYmMjKSgQMHoigKO3furPMnIpeX3eHkw91xAMwe0AKtusz/r4qJi4uloCAfPz8/gpOS0O3eidPHF13nTah3bIOu31VV2JIk1TGlJuHzO3117tyZiIgI9/N9+vSpnqhqgR+PppGUU0SonwdjOgRU6NqLlyl7nusFF91xNxqPveBwVnqskiTVXaUm4UmTJrm/t1qtxMTEANCqVSu0DeDXabtTsHyPqxc8s19oqQd3vvrqEgDmz3+y+LV2OydPHgegq4cR3fo1CK0W871z8E7YW8WRS5JU15T5wdyePXtYsGABQUFBCCFITk7mlVdeqfcb+Gw6lkZ8dhHBJgOjK9gLPnPmNEVFRQQEBBL07ZcoQmCeNh1n02aQUEUBS5JUZ5WZhF955RU+/vhjWrd2nRxx5swZHn30UVatWlUtwdUEx0W94Lv6hKIp4/j60pyfFdG9eRCGRQsBuThDkqTLKzMJ22w2dwIG13CEzVa/p1dtOZFOTKaZZj56xnWqWC/YYrEQHX0SgF57dqMUFWEZfT2Odu0BsDa+HrVH/R/OkSSp/MpMwl26dOGpp55iwgTX8Tvr1q1zf2hXHwkhWL4nHoA7+4aiqcCMCIDo6JPY7XZaNgnA9MqLQPFTlM0t56I3GSG7sPKCliSpTiszCT/33HN8+eWXfP755wD06tWLW2+9tVoCqwm7zmQRfbaAJl46buhU8fnQ54ciRsacRpWZie263tj6DajsMCVJqkcum4QdDgcTJkzgp59+4q677qrQTbdv387ixYtxOp1MmzaNWbNK37z8559/Zu7cuaxcuZKuXbtWLPIq8NleVy/4lp5B6DQV6wUXFhYSE3MGjRC0XeMaMy+c+0ixJcq++8ah1qigx/rKC1qSpDrtsplGrVbTqlUrkpKSKnRDh8PB888/z0cffcSGDRtYv3490dHRJerl5+fz2Wef0b1794pHXQUOJefyV0IOnjo1k7o1q/D158+RG5WehiYpCXu79ljHjK2CSCVJqk/KHI7Izc0lMjKSbt264eFx4Vj3ZcuWXfaaqKgoWrRoQUiIa8vHyMhItmzZQtu2bYvVW7p0Kffeey8ff/zxtcRfaT7f65o/NqV783LvlHaxo0ePoAhBv22/AlD44D9BVbHetCRJDU+Z2WbevHllFZcqNTWVpk0v7LMQGBhIVFRUsTqHDx8mJSWFYcOGlTsJKwqYTMYKx1MeMRkF/Bp9Fq1aYdawNpjKsWm7weCa5WAyGcnJySE9PZnwpHg8Y2MQwcF43H0HHpds+6nWqFBQqqwd1UWtVtX5NkD9aEd9aAPUn3ZcjVKTsMVi4euvvyYuLo527doxdepUNJqK9w5L43Q6efnll3nppZcqdJ0QkF1Fswre2xqNEDCuYyB6p7Ncr1NU5Jqql51dyB9//EWR2cqIP/4AoGD2A5gL7VBoL3aNr92JRqOqsnZUF5PJWOfbAPWjHfWhDVB/2tGkSfnOnrxYqZn1iSeeQKPR0KtXL7Zv3050dDRPP/10uW4YGBhISkqK+3FqamqxndcKCgo4ceIE//jHPwDXEUD33Xcf7733Xo18OJddaGPDkVQAbu8VfFX3OHr0MG3i4wg8FY3Tzw/zbXeUWs8SOAm1R8PaFF+SpLKVmoRPnTrFunXrAJg6dSrTpk0r9w27du1KTEwM8fHxBAYGsmHDBt544w13ube3N3v27HE/njFjBvPnz6+x2RGrDyZjsTsZ2Mqflo0q/utQRkYGqakpTPrjdwDMsx8AL69S6xaF3ItBzhOWJOkipSbhi4ceKjoModFoWLRoEffccw8Oh4MpU6YQFhbG0qVL6dKlCyNHjry2iCuR3eFk5QHX7I/pPa/uANOjRw/TIimRNjFncHr7YJ5Z+nQ8AByFYL98sSRJDU+pGfbYsWP07NkTcK0is1gs9OzZEyEEiqLw119/lXnToUOHMnTo0GLPXe5DvvMLQWrC1pNnScu30srfSN8WfhW+XgjBsWNHmPj7LgDMM2chfE2Xre+7f6qcJyxJUjGlJuGjR49Wdxw14pu/EgG4uWdzlKs4962gIJ+gtFQ6nT6F8DBintUwjn2SJKnyNNiJrIeTczmYnIe3XsO4q1iiDK4PFUf+vhsA8z/uQjRuXJkhSpLUADTYJPz1uV7wjV2b4lGBAzzPE0LgE3OGbieP49TpMN//UGWHKElSA9Agk3BGgZUtJ86iUmBa+NV9IJeUlMgtJ1wnaBTNuAtns6u7jyRJDVvlrMCoY9YdSsHuFAxp04hm5Vgdd6nExASyft1CRF4uFpWK0zfdQnkGNIqa34bRKOcJS5J0QYNLwg6nYHVUMgCTu1d8ox6A48eP8lBmBgCbWrYCS1G5krCl+W14yHnCkiRdpMENR/wem0VSroXmPnr6t6z4tDQA55/7GGe1UoTCzz3CCQkJLdd1ijUDLGev6jUlSaqfGlxP+PtzizMmdWuG6iqmpeXm5tB99UoANrVpwz3/9zxBQeVb7uwTNUPOE5YkqZgG1RNOyS1i55lMNCqFCV2bXvmCUkR//gmdok9SpFazf/TYcidgSZKk0jSoJLz6YApOASPCGuN/FR+Q5WRnEfbfDwH4tUdP8j09KztESZIamAaThO1OwdqDrt3dpvS4ug/kYj94j9bxcRR5erFrwMDKDE+SpAaqwYwJ7z6TydkCK6F+HoQH+Vb4+uzMDDp8thyA3PsfokiIyg5RkqQGqMH0hNcecvWCJ3RpelX7RCQvfZOgtFQK/PxQHvznVcVQFDwTZ5s5V3WtJEn1U4PoCWcWWtlxOhO1ApGdAip8fXZqMp2/du32lvfIfNQXnbdXEZamU+Q8YUmSimkQPeEfj6ThcAr6t/KnsZe+wtdnv/g8jbKzyW4ehHrm7KuOQ1WUAIXxV329JEn1T73vCQshWHNuKGJil4pPS8uJPkGXH1zzgvOffRH9NZy1531olpwnLElSMfW+J3wkJY8zGYX4G7UMau1f4ettTy/Aw2IhsWs39DdOqYIIJUlqyOp9El57yHWI59iOgWjUFWtu3p7faffrFpyKgv3lN658gSRJUgXV6yRssTvZdDwNgPFdKrhxuxDoFj6GWgiih4/C2LtvFUQoSVJDV6+T8I5TGeRbHHQM9KJN44qtbrN8+xXBh6Iw6/Vol7xSRRFKktTQ1esP5n486uoFX9+xgtPSCgowPfsUAIem30bL1m0rJR5zi4fw9JT7CUuSdEG9TcLZhTZ2nnHNDR7ToWJJWHnxWbwyM0kIbIr/k4sqLSZrk7EY5TxhSZIuUm+HIzafSMfhFPRp4UejCvQ+1SeO4//JRziBw/c/hI9fxWdUXPbeBSch73il3U+SpLqv3vaEfzziGooYW5EVckKgf/hBVA4He7r3IOz2Oyo1Jq+j8+Q8YUmSiqmXPeGEbDMHk3Px0KoY1rb8x9AbPv8Ez717yPfwIOH+uXh7+1RhlJIkSfW0J3y+Fzw8rHG5j7NXJSZgfOZJANaOvp7+EWPKrD9//pPXFqQkSRL1sCcshOCnY64kPK5jOecGC4H3o3NRFxRwqG0Y6um34eXlXYVRSpIkuVRJEt6+fTtjxowhIiKCDz74oET58uXLGTduHOPHj+eOO+4gMTGx0l77eFo+cVlm/I1args1lesa/bdfodv6C4V6PWuvH0efvv0rLR5JkqSyVHoSdjgcPP/883z00Uds2LCB9evXEx0dXaxOx44d+f7771m3bh1jxozhtddeq7TX33zcdZrxiLDGaFRX3jdYFR+H19MLAFg7fBRtBg2tsl5wYavHcXaUwxiSJF1Q6Uk4KiqKFi1aEBISgk6nIzIyki1bthSr069fPzzO7cnbo0cPUlJSKuW1hRD8cm6ZckSHJle+wG7H5/57UeXmcLhNW6J6hNOnT79KiaU0tkbDEYGjquz+kiTVPZX+wVxqaipNm17YMjIwMJCoqKjL1l+5ciVDhgy54n0VBUwmY5l1/k7IJinXQqC3nmGdmqG6Qk9YtfhF1Ht2U2gysXbiRAYNHkBwcMU3fS+37AOoc1WYTN2q7jWqgVqtuuJ7URfUh3bUhzZA/WnH1ajR2RFr1qzh0KFDfPHFF1esKwRkX2Gl2ap9rg3Th4c1JjfXXGZdzd49mF58AYAvRo0hT2+kU6fwK77GtfDd90/QqMiu4/OETSZjlf6cqkt9aEd9aAPUn3Y0aVLxocxKT8KBgYHFhhdSU1MJDCw5S2HXrl0sW7aML774Ap3u2vdTcArBlhOu8eBR7cqeG6xkZuBz3z0oDgcHRo7mZMtW9Am/Dk95hL0kSdWs0seEu3btSkxMDPHx8VitVjZs2MCIESOK1Tly5AiLFi3ivffeo1GjRpXyugeTcknNcw1FdG1exiILhwOf++5BHRdLYacufNO1GzqdrkrHgiVJki6n0nvCGo2GRYsWcc899+BwOJgyZQphYWEsXbqULl26MHLkSF599VUKCwuZN28eAM2aNWPZsmXX9Lqbj6cDMKpdE1RlnKZsfHUxul+34GzUiJW33IYjJ4fePXthNDbM8ShJkmpWlYwJDx06lKFDhxZ77nzCBfjkk08q9fWcQrD1pGsooqxZEbofN+D51usIlYrTi1/lYPRJ9Ho9vXr1qdR4JEmSyqteLFs+nJxHer6Vpt56OgV6lVpHfTAK7/vvBaDgqWf5yWYDoGc19oIL2j6D91Wc9ixJUv1VL5Yt/3quFzw8rDFKKUMRqqREfG+bhqogn6JJUzh542RiY2OqvRdsN/VFNB5Qba8nSVLtV+eTsBCCX6NdSXhYWMkP+ZS8XHxvnYY6JRlrvwHkvb2M33buAOC663q7F41UB032HpSzu6rt9SRJqv3q/HDEqbOFJGQX4eehpXtz3+KFRUX43D0DzZFD2NuGkfvpV8SlphAXF4vBYKj2sWDP6OdQyf2EJUm6SJ3vCZ8fihjathHqi1fIWa343HsHum2/4mzchJyvViL8/Nm16zfA1Qs2GAw1EbIkSZJb3U/C7qGIixZo2O34zJmJ/ucfcfr5kf3dGpwtWxEbG+PuBV93Xe8ailiSJOmCOp2EE7LNnEwvwFOnpneIyfWkzYb3g7PQr1+D08eXnBU/4OjcBSGEuxfcu3df2QuWJKlWqNNjwueHIga19kenUYHZjM89/0C/+Wecnl7kfPM99u7hAMTGxhAfH4fB4EHPnr1qMmxJqlccDjtZWenY7darvkdqqoIQohKjqloajQ4/vyao1deeQut0Ev5fdAZwbmpabg4+t9+M7vddOP39yfn6e+zh1wGuGRQ7z82I6NOnL3p9zczVzW//Mt7esgcu1S9ZWekYDEY8PZuWOkW0PNRqFQ6Hs5IjqxpCCAoKcsnKSqdx42bXfL86m4QzC60cTMpFp1YYqCnANGEimiOHcDRr7hqCaN/BXTcm5gyJiQl4eBgJP5eYa4LDuxuYjFAPdouSpPPsdus1JeC6RlEUPD19yM/PrpT71dkk/NvpTATQ2weCbhiB6uxZ7G3akrPiB5whoe56F/eCe/euuV4wgDbjVxSLHvRywYZUv1Q0Ab/66hKg7h6YW5n/4dTdJHzKNRQxZuUyVGfPYh06nNwPP0GY/IrVO3PmFElJiXh4GOnZs+Z6wQDGM6/JecKSJBVTJ5OwLTubPSdSQKVl5PHfKZx9PwXPvAia4s1x9YJdMyL69OlXKfsWS5JUO2VknOXtt9/g6NEjeHl54+npyZEjhwgJaUFqagpeXl54enrh62tiwYKnue22aYSGtsBut9G+fUcWLlyERqPhr7/28c03X/Dqq/+qlrjrXBLW/LWPQ8+8RuGwB+iYHoP3Sy9ScNMtxeqc/1VnypRpJCcnYTR6Eh7esybClSSpFLm5ueTm5pCYmEBQUPA1308IwZNPPs7YsZE899xLAJw8eYLCwgK6dw9n8eJnGTBgEMOHu854TE5OIigoiE8++QqHw8HDDz/A1q2bGT167DXHUlF1JwkXFuL53P/h8d6/+XXELAAGDOyCZXzpyfXiXnDfvrIXLEnVYeXKbzl9+lSZdXJzc/nhh5UIITh06CADBgzCZPLF6Sx9ilrr1m2YOvXmMu/511/70Gg03HjjVPdzYWHtyhWzWq2mY8fOpKenlat+ZaszizWU8B4Y312KAH4JHwnAoF5tL1s/KyuTlJRkPD296NFD9oIlqbbIzc1BCIFeb8DpdJKbm3PN9zx9+hTtL5oRVREWi4UjRw7Rt2/NfGBeZ3rCSnQ09o6d+OuFt0naW4S/UUunpqUfqpeTk0NU1N/07NmLESNGodVqqzna0uV3XIq3jwEcNR2JJFWNK/VYARITEzh06CBOp5Pw8J489tgCQkNDq32ecGJiInfeeSvJyYn07z+Itm3DqvX1z6szPWHnf94ja/N2/qdvCsDAVv6lHmOUmJjA//63ldjYM+zd+weNG1fhEfYV5PAMA+/2NR2GJNWooKBgBgwYRLdu3XnssQWVMibcqlVrjh8/VsE4XGPC3367huPHj/Lbb9uuOY6rUWeSMLNng07HjlOZAAxuU/oBofHxcTgcDvR6A02aNCYlJak6oyyTLv1HlKR1NR2GJNU4Hx8fgoNDKiUBg2tXRKvVypo1q9zPRUef5O+/91/xWpPJxJw5D/H5559USiwVVXeSMJBttnEoORetWqFvC79S64SEhGI0GvH29qZJkwBCLlq4UdM8Yv+N6sRbNR2GJNU7iqLw0kuvs2/fH9x000Ruv/0m3n//Hfz9y3ea+5AhwygqKnIn7X379jJp0jj316FDUVUXu6gju2Y4nYKvd57mqQ3H6B1q4j/Tul227tNPLyA3N4cnnniq0v6nrQy++8ah0ajIqOOLNUwmI9n1YOl1fWhHbWhDSkosTZu2uKZ71KW9I84rrd1NmpT+OVVZ6swHcwC7YrIA6N+y9F7weT4+Pvj4+NSqBCxJklSaOjMc4XQKdp9xjQf3b+Vfw9FIkiRVjjrTEz6akktmoY0ALx1tGpV9RH1d3RREkuoqIUSD2UUNqNS9j+tMEt52Ih2A/i396+ybndflA3x8PODq976WpFpHo9FRUJCLp6dPnf23WRHn9xPWaCpnFW6dScLbT7hO0RjQquzx4NrMaQgGoxGsdfvDIEm6mJ9fE7Ky0q9pf11FqZsna1TKvSrlLtXgr7gs1Ar0Dq27SVif8j1Knh68b6jpUCSp0qjVmms+YaI2zPKoKVXywdz27dsZM2YMERERfPDBByXKrVYr//znP4mIiGDatGkkJCRc8Z52p6Brcx+8DXXm/40SDAkfozq1rKbDkCSpFqn0JOxwOHj++ef56KOP2LBhA+vXryc6OrpYne+++w4fHx82b97MnXfeyeuvv16ue/dvKWdFSJJUv1R6Eo6KiqJFixaEhISg0+mIjIxky5Ytxeps3bqVSZMmATBmzBh2795drvGg/nV4PFiSJKk0lf67fWpqKk2bNnU/DgwMJCoqqkSdZs1cY0gajQZvb2+ysrLw9798Tzfm5cjKDrX6jXWddVc5w/k162pWBtVG9aEd9aENUH/aUVF1ZrGGJElSfVTpSTgwMJCUlBT349TUVAIDA0vUSU5OBsBut5OXl4efnxxqkCSp4an0JNy1a1diYmKIj4/HarWyYcMGRowYUazOiBEjWL16NQA///wz/fr1axCTvCVJki5VJbuobdu2jSVLluBwOJgyZQr33XcfS5cupUuXLowcORKLxcLjjz/O0aNH8fX15a233iIkJKSyw5AkSar16sxWlpIkSfWR/GBOkiSpBskkLEmSVINqXRKuiiXP1e1KbVi1ahX9+vVj4sSJTJw4ke+++64GoizbwoUL6d+/PzfcUPo+F0IIXnzxRSIiIhg/fjyHDx+u5gjL50rt2LNnD9ddd537vXjnnXeqOcIrS05OZsaMGYwbN47IyEg+/fTTEnXqwvtRnnbUhffDYrEwdepUJkyYQGRkJG+//XaJOhXKU6IWsdvtYuTIkSIuLk5YLBYxfvx4cfLkyWJ1vvjiC/F///d/Qggh1q9fL+bNm1cDkV5eedrw/fffi+eee66GIiyfP/74Qxw6dEhERkaWWv6///1PzJw5UzidTrF//34xderUao6wfK7Ujt9//13MmjWrmqOqmNTUVHHo0CEhhBB5eXli9OjRJf5O1YX3ozztqAvvh9PpFPn5+UIIIaxWq5g6darYv39/sToVyVO1qidclUueq0t52lAX9O7dG19f38uWb9myhRtvvBFFUejRowe5ubmkpaVVY4Tlc6V21AUBAQF07twZAC8vL1q3bk1qamqxOnXh/ShPO+oCRVHw9PQEXOsc7HZ7iSm2FclTtSoJl7bk+dI36XJLnmuL8rQBYNOmTYwfP565c+e6F67UJZe2s2nTpnXyHxTAgQMHmDBhAvfccw8nT56s6XDKlJCQwNGjR+nevXux5+va+3G5dkDdeD8cDgcTJ05kwIABDBgwoNT3o7x5qlYl4YZi+PDhbN26lXXr1jFgwACeeOKJmg6pwercuTNbt25l7dq1zJgxgwceeKCmQ7qsgoIC5s6dy5NPPomXl1dNh3PVympHXXk/1Go1a9asYdu2bURFRXHixImrvletSsL1Yclzedrg5+eHTuc6GmXatGm18kOUK7m0nSkpKSXaWRd4eXm5f7UcOnQodrudzMzMGo6qJJvNxty5cxk/fjyjR48uUV5X3o8rtaOuvB/n+fj40LdvX3bs2FHs+YrkqVqVhOvDkufytOHisbqtW7fSpk2b6g7zmo0YMYIffvgBIQQHDhzA29ubgICAmg6rwtLT091jdVFRUTidzlr1nzq4Zj489dRTtG7dmrvuuqvUOnXh/ShPO+rC+5GZmUlubi4ARUVF7Nq1i9atWxerU5E8VauOqdBoNCxatIh77rnHveQ5LCys2JLnqVOn8vjjjxMREeFe8lyblKcNn3/+OVu3bkWtVuPr68tLL71U02GX8Mgjj/DHH3+QlZXFkCFDeOihh7Db7QDccsstDB06lG3bthEREYGHhwdLliyp4YhLd6V2/Pzzz3z99deo1WoMBgNvvvlmrfpPHeDPP/9kzZo1tGvXjokTJwKudiUlJQF15/0oTzvqwvuRlpbGggULcDgcCCG4/vrrGT58+FXnKblsWZIkqQbVquEISZKkhkYmYUmSpBokk7AkSVINkklYkiSpBskkLEmSVINkEpauSnp6Og8//DCjRo1i8uTJ3HvvvZw5c6bMa8LDwwHXktXL7WpWEddyn6eeeoro6GgAli1bVuHX7datm3uaVVXYsGED7733XqllcXFxTJw40f3zlOo2mYSlChNC8OCDD9KnTx9++eUXVq1axaOPPkpGRkZNh1Zuixcvpm3btgC8//77Fb4+NDSUNWvWVHZYbtu3b2fw4ME18tpS9apVizWkuuH3339Ho9Fwyy23uJ/r0KGD+/uPPvqIH3/8EavVSkREBHPnzi3XfR9++GEmTpzIsGHDAFiwYAHDhg0jIiKC119/nT/++AOr1cptt93G9OnTi11rsVh49tlnOXToEGq1mgULFtCvXz8cDgevv/46O3bsQFEUbrrpJmbMmMGMGTOYP38+P//8M0VFRUycOJG2bdsSGhqKr68vd955JwBvvfUW/v7+3HHHHZeNOyEhgXvuuYcePXqwf/9+unTpwpQpU3j77bfJzMzk9ddfp1u3bvz73/8mISGB+Ph4kpOTWbhwIQcOHGDHjh0EBASwbNkytFotQgiOHj1K586d+eOPP1i8eDHg2r3riy++qNP7RkglyZ6wVGEnT550b0l4qd9++43Y2FhWrlzJmjVrOHz4MHv37i3XfceNG8ePP/4IuDbF3r17N8OGDWPlypV4e3vz/fff8/3337NixQri4+OLXfvll18CsG7dOt544w0WLFiAxWLh22+/JTExkR9++IF169Yxfvz4Ytc99thjGAwG1qxZwxtvvMGUKVPcvUyn08mGDRuYMGHCFWOPi4vjrrvu4scff+TMmTOsW7eOr7/+mvnz5xcb7oiLi+PTTz/lvffe4/HHH6dv376sW7cOg8HAtm3bADhy5AgdOnRAURT++9//smjRItasWcOXX36JwWAo189SqjtkT1iqVDt37mTnzp3ceOONABQWFhITE0Pv3r2veO2QIUNYvHgxVquV7du306tXLwwGAzt37uT48eP8/PPPAOTl5REbG0vLli3d1/7555/cfvvtALRp04bmzZtz5swZdu/ezfTp09FoXH/VTSZTmTEEBwdjMpk4cuQIZ8+epVOnTuXauyA4OJj27dsD0LZtW/r374+iKLRv357ExMRibdRqtbRr1w6Hw8GQIUMAaNeunfv0hR07drif79mzJy+//LJ7w5vzm9tI9YdMwlKFhYWFuRPipYQQzJo1q8RwQXno9Xr69OnDjh07+PHHHxk3bpz7nk8//XSJMdKqOtpq2rRprFq1irNnzzJlypRyXXN+VzwAlUrlfqwoCg6Ho0Q9lUqFVqt174ugUqnc9Xbu3Ok+MmfWrFnufSFuueUWPvroozq54ZN0eXI4Qqqwfv36YbVa+fbbb93PHTt2jH379jFo0CC+//57CgoKANdWnhX5wG7cuHGsWrWKffv2uZPuoEGD+Prrr7HZbACcOXOGwsLCYtf16tWLdevWucuTk5Np3bo1AwYM4Ntvv3Vv2pOdnV3iNTUajfveAKNGjWLHjh0cPHiQQYMGlTv2ypCXl4fdbnf3vuPi4mjfvj2zZs2ia9euV5yBItU9sicsVZiiKLzzzjssWbKEDz/8EL1eT1BQEE8++SQtW7bk1KlT7p6w0Wjktddeo1GjRuW698CBA5k/fz4jR44studyYmIikydPRgiBn58f//nPf4pdd+utt/Lss88yfvx41Go1L730EjqdjmnTphETE8OECRPQaDTcdNNN7mGL82666SYmTJhAp06deOONN9DpdPTt2xcfHx/UanUl/MTKb+fOnQwYMMD9+NNPP2XPnj0oikJYWJh7mEKqP+QuapJ0CafTyaRJk1i6dGmxcefzEhISmDNnDuvXr6/0137qqaeYNm0aPXr0uGLd8PBw9u/fX+kxSNVLDkdI0kWio6OJiIigf//+pSZgcB1tk5eXVyWLNRYvXnzFBHx+sUZ5f7uQajfZE5YkSapBsicsSZJUg2QSliRJqkEyCUuSJNUgmYQlSZJq0P8DmsmFH/MupbYAAAAASUVORK5CYII=\n" + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "results_ctrl = np.loadtxt(results_ctrl_file)\n", + "\n", + "v_min = 0.\n", + "v_max = 3.\n", + "\n", + "v_ctrl = results_ctrl[:,0]\n", + "probs_ctrl = results_ctrl[:,3]\n", + "probs_ctrl_err = results_ctrl[:,4]\n", + "\n", + "# fit_bounds=(0, [1, np.inf, np.inf])\n", + "\n", + "ind_vmax = v_ctrl <= v_max\n", + "v_ctrl = v_ctrl[ind_vmax]\n", + "probs_ctrl = probs_ctrl[ind_vmax]\n", + "probs_ctrl_err = probs_ctrl_err[ind_vmax]\n", + "\n", + "popt_ctrl, pcov_ctrl = curve_fit(logistic_growth_simple, v_ctrl, probs_ctrl,\n", + " sigma = probs_ctrl_err, absolute_sigma=False,\n", + " # bounds = fit_bounds\n", + " )\n", + "popt_ctrl_exp, pcov_ctrl_exp = curve_fit(exponential_inverted, v_ctrl, probs_ctrl,\n", + " sigma = probs_ctrl_err, absolute_sigma=False,\n", + " # p0=-15,\n", + " # bounds = fit_bounds\n", + " )\n", + "\n", + "v_plot = np.linspace(v_min, v_max, 100)\n", + "\n", + "with sns.axes_style(\"darkgrid\"):\n", + " plt.figure(0,(5,3.5))\n", + " plt.errorbar(v_ctrl, probs_ctrl, probs_ctrl_err, color='k', marker='.', alpha=.5, label='CTRL', lw=2, zorder=1)\n", + " plt.plot(v_plot, logistic_growth_simple(v_plot, *popt_ctrl), 'r-', alpha=1, lw=2, zorder=100)\n", + " plt.plot(v_plot, exponential_inverted(v_plot, *popt_ctrl_exp), c='C0', alpha=1, lw=2, zorder=100)\n", + " plt.axhline(1, c='C0', ls='--', lw=1.5)\n", + " plt.vlines(0.5, -0.1, 0.5, ls='--', color='orange', lw=1.5)\n", + " plt.xlim(0,3)\n", + " plt.ylim(0,1.05)\n", + " plt.xlabel(\"Cell velocity [mm/s]\")\n", + " plt.ylabel(\"Probability\")\n", + " plt.title(\"Probability for deformed cell\", fontsize=12)\n", + " plt.legend(loc='lower right')\n", + " plt.tight_layout()\n", + "\n", + " # savename = \"fig2A_explain_LogGrowth_FitValues\"\n", + " # savepath = os.path.join(savefolder,savename)\n", + " # plt.savefig(savepath+\".pdf\", dpi=900, format='pdf')" + ], + "metadata": { + "collapsed": false, + "pycharm": { + "name": "#%%\n" + } + } + }, + { + "cell_type": "code", + "execution_count": 19, + "outputs": [ + { + "data": { + "text/plain": "array([2.03423702])" + }, + "execution_count": 19, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "popt_ctrl_exp" + ], + "metadata": { + "collapsed": false, + "pycharm": { + "name": "#%%\n" + } + } + }, + { + "cell_type": "code", + "execution_count": null, + "outputs": [], + "source": [], + "metadata": { + "collapsed": false, + "pycharm": { + "name": "#%%\n" + } + } } ], "metadata": {