Content - 5e9f4a3f5018083acaba6d0b9998f37e55cd4abf - 95b28ff/divergence-tree-analyses/Regression-analysis-on-descendants-in-divergence-tree.ipynb

visit type
https://github.com/blab/mumps-wa-phylodynamics

20 April 2021, 13:07:29 UTC
Tip revision: b8358a0d49d70670dbab9eeffa9972c277b3021b authored by Louise Moncla on 17 March 2021, 19:25:20 UTC
adding in callouts for Wisconsin H and A genomes
Tip revision: b8358a0
Regression-analysis-on-descendants-in-divergence-tree.ipynb
{
 "cells": [
  {
   "attachments": {
    "Screen%20Shot%202021-03-16%20at%202.03.57%20PM.png": {
     "image/png": "iVBORw0KGgoAAAANSUhEUgAABIgAAAK9CAYAAABPS1fnAAABQGlDQ1BJQ0MgUHJvZmlsZQAAKJFjYGDiSSwoyGFhYGDIzSspCnJ3UoiIjFJgf8rAzCDAwM0gzaCUmFxc4BgQ4ANUwgCjUcG3awyMIPqyLsise2neZVfPJkq84A7UXHh/6hVM9SiAKyW1OBlI/wHihOSCohIGBsYYIFu5vKQAxG4AskWKgI4CsqeA2OkQ9goQOwnC3gNWExLkDGRfALIFkjMSU4DsB0C2ThKSeDoSG2ovCLAZmRiZhhNwKKmgJLWiBEQ75xdUFmWmZ5QoOAJDJ1XBMy9ZT0fByMDIkIEBFNYQ1Z9vgMOQUYwDIZbpwsBg9piBgSkVIZa8lYFhpwgDg+BWhJjmXaB3uBgY9gsVJBYlwh3A+I2lOM3YCMLm3s7AwDrt///PQG+yazIw/L3+///v7f///13GwMB8i4HhwDcAYcZd5yDN9FYAAABWZVhJZk1NACoAAAAIAAGHaQAEAAAAAQAAABoAAAAAAAOShgAHAAAAEgAAAESgAgAEAAAAAQAABIigAwAEAAAAAQAAAr0AAAAAQVNDSUkAAABTY3JlZW5zaG9003s/NQAAAddpVFh0WE1MOmNvbS5hZG9iZS54bXAAAAAAADx4OnhtcG1ldGEgeG1sbnM6eD0iYWRvYmU6bnM6bWV0YS8iIHg6eG1wdGs9IlhNUCBDb3JlIDUuNC4wIj4KICAgPHJkZjpSREYgeG1sbnM6cmRmPSJodHRwOi8vd3d3LnczLm9yZy8xOTk5LzAyLzIyLXJkZi1zeW50YXgtbnMjIj4KICAgICAgPHJkZjpEZXNjcmlwdGlvbiByZGY6YWJvdXQ9IiIKICAgICAgICAgICAgeG1sbnM6ZXhpZj0iaHR0cDovL25zLmFkb2JlLmNvbS9leGlmLzEuMC8iPgogICAgICAgICA8ZXhpZjpQaXhlbFhEaW1lbnNpb24+MTE2MDwvZXhpZjpQaXhlbFhEaW1lbnNpb24+CiAgICAgICAgIDxleGlmOlVzZXJDb21tZW50PlNjcmVlbnNob3Q8L2V4aWY6VXNlckNvbW1lbnQ+CiAgICAgICAgIDxleGlmOlBpeGVsWURpbWVuc2lvbj43MDE8L2V4aWY6UGl4ZWxZRGltZW5zaW9uPgogICAgICA8L3JkZjpEZXNjcmlwdGlvbj4KICAgPC9yZGY6UkRGPgo8L3g6eG1wbWV0YT4KhBWDsAAAQABJREFUeAHs3QWclNX+x/HfFkt3SCkKKtgBXBvFFjCwEbuwr4Wd2OK1W7EL+4KK2HGNv91BKChId2zvf78Hz/DMMLs7Ozzj7Ox+zuu1zMxTc573M8MyX37nTFZ5RbMqmlZnZWVVsUXNVy1ZssSys7OtcePGNd+5ij3Kysps3rx51rZt2yq2qvmqgoICmzJlittx/fXXr/kBqtijtKTUFi1eZK1atapiq+RWLVy40Bo2bGj5+fnJHaCSvYqKimzZsmXWsmXLSrZIbrH6O2PGDNffbt26JXeQSvbSNVS/mzdvXskWyS/Wa65F8xaWk5uT/EHi7Ll8+XIrKSmxZs2axVmb/KJZs2bZ/PnznUXHjh2TP1CcPTPtva1TWLx4seXm5lqjRo3inFHyi6ZOnWq6hu3btw/9/Z1p720pLliwwJo0aWJ5eXnJo8bZc8KECaa/+7t27Rr675RMe2+LZ86cOdamTZtQf28XFxfb5MmTnX6P7j1C/bsuVb+31dlUvbdT9Xtbf39OmzbNcnNyrXuP7s47rD9S9Xtb/UvVeztVv7f1vp49e7b7+6hLly5hEbvjpOr3tg6eive2jpuq39vTp09378HWrVtbu3bt9FShtUx7b+vEU/V7+7fffnP/vtW/58L+N26q3tuhvRA4EAIIpEQgOyVH5aAIIIAAAggggAACCCCAAAIIIIAAAhkjQECUMZeKjiKAAAIIIIAAAggggAACCCCAAAKpESAgSo0rR0UAAQQQQAABBBBAAAEEEEAAAQQyRoCAKGMuFR1FAAEEEEAAAQQQQAABBBBAAAEEUiNAQJQaV46KAAIIIIAAAggggAACCCCAAAIIZIwAAVHGXCo6igACCCCAAAIIIIAAAggggAACCKRGgIAoNa4cFQEEEEAAAQQQQAABBBBAAAEEEMgYAQKijLlUdBQBBBBAAAEEEEAAAQQQQAABBBBIjQABUWpcOSoCCCCAAAIIIIAAAggggAACCCCQMQIERBlzqegoAggggAACCCCAAAIIIIAAAgggkBoBAqLUuHJUBBBAAAEEEEAAAQQQQAABBBBAIGMECIgy5lLRUQQQQAABBBBAAAEEEEAAAQQQQCA1AgREqXHlqAgggAACCCCAAAIIIIAAAggggEDGCBAQZcyloqMIIIAAAggggAACCCCAAAIIIIBAagQIiFLjylERQAABBBBAAAEEEEAAAQQQQACBjBHIuIDozjvvtPPOO28V4Pvvv9/uvvvuVZZfcMEF9sgjj6yyvLoFw4cPtyOOOCLuZgsWLLDdd9/dXnjhhbjra7KwrLysJpsntO2PP/5oxxxzjM2dOzdq+7KyMrf87bffjizX/WOPPdYOOOAAu+KKK2zevHmRdf6ODC+88EL/MO7tWWedZffcc09knY4zbNgwGz9+fGSZ7syYMcP1Qf3zP1U5qs/q16233ho5zieffBLZ1x/D306ePDmyne48/fTTbtu//vrLTjnllFX207HV9Lq67rrr3P3gH998843bZ/bs2cHF9tVXX7nXxx9//BG1PPjgtttus6uvvjq4KHJfxz3hhBPs448/jizzd3799Vfr06ePFRcXu0Xz58+3m2++2fVj0KBBdtppp9n777/vN+cWAQQQQAABBBBAAAEEEEAAgdUWyLiAqKioyG644Qb76aefIievZeeee65ddNFFVlBQEFk+YcIE96G/tLQ0sizRO3PmzLGZM2fG3XzhwoUu+Pj555/jrq/Jwnbt2rnQpCb7VLftn3/+aQ899NAqAVF5Wblb7vv98ssv284772w6186dO9ujjz5q+++/f9Thv//+e2d47bXXmu7Ha7NmzXIBxt577+1WT5o4ybbaaiu79957o66TVn7+2ecuWCspKTH/07Zt23iHtSVLltg+++xjl19+ub377ruRbfLy8iL7+mPo9aBzzsrKimynOwpodtttN2fx+OOPO+v8/HzzP40aNXLbv/nmmzZmzJioffXgj6l/uOMqFPTt1VdftS222MIee+wxU3gTr8n2jDPOsJdeemmV1Xr96Lz0+uzbp+8q62+//XZ3fJ3ne+++Zz179nQhmfbTdVJA1q9fP9N2NAQQQAABBBBAAAEEEEAAAQTCEMgN4yBhHaO8aJllNWhc5eH2228/U7WKPsz36tXLbTv+9fGmD89qb731lh155JHuvv/A74MLtzCEP9Zaay1TKKUP8KvTFKz4fq/OcZLdVwHOwIEDTWGGmoK32KojbbPuuutabm6uC3zihRKjRo1yx+nUqZN9/vnnNmDAANt8881N5xfbvv3uW/vXv/7lwqjYdcHHqkBSlZZue/fuHVxlW2655Sr7K+gaPHiwrb322pFtFaSoD1r+7bffuuV67eyyyy6RbWp658EHH7TjjjvODjzwQHv22Wfj7j5lyhT3Gtxmm23c6yR2I1VsqTpo9OjRlpObE7V6+fLlLnh65513XPh0zLHH2Prrr28vvviiBYO09957zzbYYIOofXmAAAIIIIAAAggggAACCCCAQLICaa8gKls0zRaOu8hm3Nnf/rprN3erx1oer3Xr1s023XRTe+WVVyKrRz872gUUO+20U1TFhrbZbrvt3AdrDenZY489rEOHDu6D9dgxY93+GvKkag4t33jjje2+++6LHFcVLCeeeKJ1797dFDJ9+umnkXWqIPF92HXXXc2HJDqOwoqpU6e6bQsLC90xVPmhn6FDh9qpp57qhj4ddthhbpvjjz/e9b+8vNw0lElD23Seet6zzz7bFBqoKcjRULAbb7zRnYO20dCjZFt2drb9/vvvkZBKVTUKeXzzYYUCDfVJFTO+L34b9VnD+zScTK1Fixbu/F4Z+4qr0vHb+dvvvvvOFLCpGkmVS5W1xo0bW//+/Z25HKpqGm6loXKXXHJJ1GYa8qZ+N2jQIGr56jzo0L6DPfH4EzZixIi4h1Hwo/BIwdvhhx++yjZ33XWXPf/88+51quqx2KbQSIGcAjYNjdTr4eGHHo4Kh7SPKoji7R97PB4jgAACCCCAAAIIIIAAAgggkIhAWgOi4pnf26zHj7TsGZ/YGl0bWKeezd1tVsXjWU8cZcUzf4h7Dgp0FApo2I8CmP/+97928MEH20EHHeQ+eCvEWLRokQsN9t13X3cMVcco1FElxkYbbWRDDx9qy5Ytc3PraH4aVWTcdNNNrlLGP+lHH31kHTt2NFWNtGzZ0jQvkW8KOPQcapr3RpUpCgYUMGkIl58nScfUMKNnnnnGBSxjx451+/Xp3ccUaKkpKBp24oqARdU3Gi6l4V6nn366aVjUSSed5LZTtZHChS+++MLuuOMO23PPPd3zas6aZJqeU+ehwEYBj6p/gk3z9+g5VZF11FFHufsKMILtjTfeMA3hU1/UFG4oqImtjPH76NwUKCmMU8ChcC04fMtv17BhQ7v++utXCUb8+uCtnk8B3mabbRZZrH4rYFH4FmxffvmlqzJTpZl+goGXwkJdv+DP6+NfD+5uAwcNtCGHDYlaFnxwzjnnuPMJzsfk12von+ZBUlOApAqg2O302PdZ16ZHjx7WvUfVAZk/PrcIIIAAAggggAACCCCAAAIIJCuQ1iFmi16/0pq3MGvcruHK/jfMseadcix3dqEtemOEtRn69Mp1f99T6HPllVfaa6+9Zqo0URig0EiBj8IUhTCqbFHTkDS1J554wt3qjyZNmrjhQRMnTrTFixe7uW5ULaS5aoJNjzX/jeY1UkBy8sknu6qX4FAfv/35558fGdqm4UF+cmZVEmmokSqZ1FQZ4qqDKj70qwpEbccdd3TD5RRefPbZZ27iYs3hs8MOO7hhRppIORgkKLhRU+Cl5aqOWm+99dyymvyhsEMVPTqGhktpONlll13mzlnHUVCiIGONNdZwh9V9LTv00EMjT6N9ZKNqpESa5vrRtq1bt3aBlIZ76TmDk1Anchy/jUIehYUKzYJN4ZBcY12uueaa4GYm83XWWcct0/A6hY3BNn369ODDKu9rsm1NTK2+NG3adJVt/STYChy33nprFxrq9dqlcxcXPCkQUoA2ZMiKAEqv62bNmkUdR5VS//nPf9wyvUYuvvjiqPU8QAABBBBAAAEEEEAAAQQQQCAZgcQ+1Sdz5Gr2Kfnz/6x48ayKcCg/7pZaXrxwphVN/WSV9QpZVPWiD/MKS1Q9og/SqkjZfvvt3QdvhUSqFPIf/jUvzL///W/bdtttba+99nLH1DxCmsS4TZs27lujFCZoAuJ4TYGGmvaprqkfqmxSU5iiwEgVRvr2KU22rGFi8Zr/Bq7g3DIKAdTU/9imqiY1/1x+vQIwNQUMwbZw0YrHfr3WyUjVSKpuUTChMErDpBRWaA4fDUFTf/Wjb+xSVdUPP6yo7FLFjUIRzcmTaFO45i01t9AhhxwS95u8Ej3epZde6kIshWXBpiolDQ+Mbc8995yr8FHVkn7860Pbbbjhhi5c1GvH/1x5xZWxh4j7uLik2H3LmI738MMPu+ovVY0pIFQlmK6fhijqNaBvXNP8WQofZTBm7IrJsRW2qRLOh0t6Heg14b/NTE+scFFVV5pbSdeHhgACCCCAAAIIIIAAAggggEAYAmkLiMrmTrT8hlU/vdaXzJkQ9zwVCqmCSD8a2uWblmsYlsIjVRWpqTqob9++rlro7rvvdkPP/PYKhT744AP3QVwBga/e8OtX91Yhjn4UQunD/aRJk2yTTTaJe9guXbq45dOmrZx/SRU+apq/KNGm81BTdU2wjXttnHuoUEgtGDxowu0Bew1wyzVkTGFF+/btXUWW5gLSj76WXcseeOABt50qYVRV5PvtFtbwD31TnOZtSqaNGzfOBVaqQAo2hVh+curg8lTez8vNc8MVFXgprNSP5nTS5N7Nmzd39zVszn9rWrAvmiPJz/cUDLVkq5DvlltuiWyuAErfjhYMESMruYMAAggggAACCCCAAAIIIIBAkgLpG2KWk1cxDKzqXrthYtnxu6hhZv4btXwQpKMNGjTIfeW9Plj74WWqglFgoOE8G26wod30n5siT6wKGH2NuD5wq5rjq6++iqwL487XX3/tqkc0DExBiipY1EdVkGhCZ7Uff/zRBQGqglEAcO6557qvVldfVAmjb/PyVSWV9UlD0zTPkqqBFLhoTiENaVJIoWFc+havCy+80N3Xt4Dp28H0XKr+0RxImk/pqquvcuvlrgmpFVZoWF2wyVLDzDS0SX274/Y7gqvdfQVyfoifqps0jE9VSwq+LrjgAjckTcGcJt3W/Ez+28A039KHH35o1157rTuOhvYpxPJfZa/jKGTx3x6n6iHN2RT7LWcaMqcJosOcnNqfZFlZmS1dutT9aJm/36JirGRwjiqtUz80bO2qq67SQ3fdNcRNFUAKCVVppHmfNE/VU0895eZvClZC9enTx71mdVzNk6U5tjQnlqq7NIxNVXQ0BBBAAAEEEEAAAQQQQAABBMIQqLqEJ4xnqOQYOe02sMJlxWZllaREFcsLlxZbXseN4x6h3w79XMDih5f5jVSpo7ln9OFZQYhar5693PxC+sDdqXMnF9ioEiYrK8t9QNccQary0fwx+pYpNa3Tj29+jp3gsuB9v17b5+TkuB/dVxCjb11TCKWwQ2GIwih9g5cqefSNZxq+pcoTzaGk0ERhlvqnYEhBioIEteDzuQV//5Gdle2CLYVP/tvTVHWiShPdatidwiEFaT6M0TAvBWyqtFIoIbM111zTTeysbRSwxRs6pmVapwmkdT57DVgxXC/YHwVhqprReWiybt1//fXXXVilwEgBia6PghNVdPkhdxqKpyF+vrJJw/O0ryrCNORL9/2cO6ocUyimIXHBpmFjCrf8RM/BdbofvE6JrMvKXvka0Paq6FI//GtLrx25LVu+LPZwqzw+88wzXXCl81eIp3mXdG0115TCpHh9lrMqtVQtpXmLulVUoamyqGvXrpG5olZ5IhYggAACCCCAAAIIIIAAAgggUEOBrIpKj0oSmhVH0urKgokaPldkc1WY6IP6/CcOtWa5i61550aRdf7Ooj+X25Ky5tbpuJf8ompvVd2hQCbeJNLaWes0qbWqUIJN56hKj8r2UyWLnwNo/fXXD+5a7X0FT6NGjbLBgwe7bVVJpG/b0vEULJSWlNq06dPc/eDBFKTIyFfLBNdVdl/Bja9K0jZ6rHNV5Y2Cs8quo4IpDYVKtOJG8zBpEmgFQf7byyrrU7zlOjdV3vi5iPw2Cs00z1GrVq1cEOKX1+RWFTkaMqgwJth0DdVvhTthN72uVEFU2Te3xT6fXqe6NjpPNV0fBV0jRoyIGoKmYWeqnvITVeux3jutW7VO+Llin1uPFdzNnz/fWagiKczm39t6n4XZqntvr85zyV+v/3jD/1bnuAprdc0U9vprvTrHC+7r39uqEAyz6T2ivw/091bYTeGtKglr8ndaIn2YMGGC6fWh0DTs111N39uJ9FfbxL63E90vke3092hVf98ncozYbRTa+znyenTvsVp//8QeOxPf2/q9vWjxotDf1/r7U5W2uTm5oX97ZSa+t1P1e1vv69mzZ7u/j1ZniHzsa1mPM+29rT6n6ve2vuRDv1/1bz3NzRlmS9Xv7VS9t3Xuqfq9/dtvv7l/3+rfc2H/GzdVv7fDfC1wLAQQCF8g/vit8J8n7hELO/Sx8qnvWOHkpdaidV7FB4fsiuqRig/P84qtsCTHirr1jrtfsgtjAwl/HAUnlYVDfptkb1UZpEmJn3zySXcIfYuXhgwpHPLNBwD+sW6T+eAXDIeCx6ru3JL5UKV5n5L9EKlzi3d++uCYTF+C56rhZrFDzoLra8N9BX/BwEDXf+TIkdV2TQFG2CFGtU/KBggggAACCCCAAAIIIIAAAvVCIK0BkeU2sHmNe1ijwjlWOHuhZZcVWll2ni3PbmPLG7e1ptkNMv4iaOiQ5gPS16k3zG9omoNGw8ZoCCCAAAIIIIAAAggggAACCCCAQG0RSG9A9LfC8vy2ttza1haT0PuhIWX6oSGAAAIIIIAAAggggAACCCCAAAK1USBtk1TXRgz6hAACCCCAAAIIIIAAAggggAACCNRHAQKi+njVOWcEEEAAAQQQQAABBBBAAAEEEEAgIEBAFMDgLgIIIIAAAggggAACCCCAAAIIIFAfBQiI6uNV55wRQAABBBBAAAEEEEAAAQQQQACBgAABUQCDuwgggAACCCCAAAIIIIAAAggggEB9FCAgqo9XnXNGAAEEEEAAAQQQQAABBBBAAAEEAgIERAEM7iKAAAIIIIAAAggggAACCCCAAAL1UYCAqD5edc4ZAQQQQAABBBBAAAEEEEAAAQQQCAgQEAUwuIsAAggggAACCCCAAAIIIIAAAgjURwECovp41TlnBBBAAAEEEEAAAQQQQAABBBBAICBAQBTA4C4CCCCAAAIIIIAAAggggAACCCBQHwUIiOrjVeecEUAAAQQQQAABBBBAAAEEEEAAgYAAAVEAg7sIIIAAAggggAACCCCAAAIIIIBAfRQgIKqPV51zRgABBBBAAAEEEEAAAQQQQAABBAICBEQBDO4igAACCCCAAAIIIIAAAggggAAC9VGAgKg+XnXOGQEEEEAAAQQQQAABBBBAAAEEEAgIEBAFMLiLAAIIIIAAAggggAACCCCAAAII1EcBAqL6eNU5ZwQQQAABBBBAAAEEEEAAAQQQQCAgQEAUwOAuAggggAACCCCAAAIIIIAAAgggUB8FCIjq41XnnBFAAAEEEEAAAQQQQAABBBBAAIGAAAFRAIO7CCCAAAIIIIAAAggggAACCCCAQH0UICCqj1edc0YAAQQQQAABBBBAAAEEEEAAAQQCAgREAQzuIoAAAggggAACCCCAAAIIIIAAAvVRgICoPl51zhkBBBBAAAEEEEAAAQQQQAABBBAICBAQBTC4iwACCCCAAAIIIIAAAggggAACCNRHAQKi+njVOWcEEEAAAQQQQAABBBBAAAEEEEAgIEBAFMDgLgIIIIAAAggggAACCCCAAAIIIFAfBQiI6uNV55wRQAABBBBAAAEEEEAAAQQQQACBgAABUQCDuwgggAACCCCAAAIIIIAAAggggEB9FCAgqo9XnXNGAAEEEEAAAQQQQAABBBBAAAEEAgIERAEM7iKAAAIIIIAAAggggAACCCCAAAL1UYCAqD5edc4ZAQQQQAABBBBAAAEEEEAAAQQQCAgQEAUwuIsAAggggAACCCCAAAIIIIAAAgjURwECovp41TlnBBBAAAEEEEAAAQQQQAABBBBAICBAQBTA4C4CCCCAAAIIIIAAAggggAACCCBQHwXSGhCVl5VVaV7d+ip3ZiUCCCCAAAIIIIAAAggggAACCCCAQEICaQ2IipfMq7KT1a2vcmdWIoAAAggggAACCCCAAAIIIIAAAggkJJDegKhgaZWdLC6sen2VO7MSAQQQQAABBBBAAAEEEEAAAQQQQCAhgbQGROVlJVV2sry06vVV7sxKBBBAAAEEEEAAAQQQQAABBBBAAIGEBNIaECXUQzZCAAEEEEAAAQQQQAABBBBAAAEEEEipAAFRSnk5OAIIIIAAAggggAACCCCAAAIIIFD7BQiIav81oocIIIAAAggggAACCCCAAAIIIIBASgUIiFLKy8ERQAABBBBAAAEEEEAAAQQQQACB2i9AQFT7rxE9RAABBBBAAAEEEEAAAQQQQAABBFIqQECUUl4OjgACCCCAAAIIIIAAAggggAACCNR+AQKi2n+N6CECCCCAAAIIIIAAAggggAACCCCQUgECopTycnAEEEAAAQQQQAABBBBAAAEEEECg9gsQENX+a0QPEUAAAQQQQAABBBBAAAEEEEAAgZQKEBCllJeDI4AAAggggAACCCCAAAIIIIAAArVfgICo9l8jeogAAggggAACCCCAAAIIIIAAAgikVICAKKW8HBwBBBBAAAEEEEAAAQQQQAABBBCo/QIERLX/GtFDBBBAAAEEEEAAAQQQQAABBBBAIKUCBEQp5eXgCCCAAAIIIIAAAggggAACCCCAQO0XICCq/deIHiKAAAIIIIAAAggggAACCCCAAAIpFSAgSikvB0cAAQQQQAABBBBAAAEEEEAAAQRqvwABUe2/RvQQAQQQQAABBBBAAAEEEEAAAQQQSKkAAVFKeTk4AggggAACCCCAAAIIIIAAAgggUPsFCIhq/zWihwgggAACCCCAAAIIIIAAAggggEBKBQiIUsrLwRFAAAEEEEAAAQQQQAABBBBAAIHaL0BAVPuvET1EAAEEEEAAAQQQQAABBBBAAAEEUipAQJRSXg6OAAIIIIAAAggggAACCCCAAAII1H4BAqLaf43oIQIIIIAAAggggAACCCCAAAIIIJBSAQKilPJycAQQQAABBBBAAAEEEEAAAQQQQKD2CxAQ1f5rRA8RQAABBBBAAAEEEEAAAQQQQACBlAoQEKWUl4MjgAACCCCAAAIIIIAAAggggAACtV+AgKj2XyN6iAACCCCAAAIIIIAAAggggAACCKRUgIAopbwcHAEEEEAAAQQQQAABBBBAAAEEEKj9AgREtf8a0UMEEEAAAQQQQAABBBBAAAEEEEAgpQIERCnl5eAIIIAAAggggAACCCCAAAIIIIBA7RcgIKr914geIoAAAggggAACCCCAAAIIIIAAAikVICBKKS8HRwABBBBAAAEEEEAAAQQQQAABBGq/AAFR7b9G9BABBBBAAAEEEEAAAQQQQAABBBBIqQABUUp5OTgCCCCAAAIIIIAAAggggAACCCBQ+wUIiGr/NaKHCCCAAAIIIIAAAggggAACCCCAQEoFclN69BQc/M4777SpU6fa9ddfH3X0+++/34qKiuzggw+OWn7BBRdYz5497cgjj4xaXt2D4cOH24wZM+zRRx9dZdMFCxa45znxxBNt8ODBq6yvyYKy8rKabJ7Qtj/++KONHDnSLr74YuvcuXNkn7KyMjvuuONs6NCh1r9/f7f87bfftieeeMIWLlxoG2+8sZ122mnWunXryD66I8OsrCy75ppropYHH5x11lm23nrr2bBhw+zzzz83XY9Zs2bZRhttZFrXqlWr4OZR95966il7//337brrrossv+iii+yvv/6KPPZ3/PH8419//dUOO+ww++ijj+ycc86xxYsX+1WR27vuusveeecde/LJJ9311Ln4VlhYaCeddJIdffTRtv322/vF7vabb74xvd60buutt7ZPPvnE7rvvvqht9EDHuOyyy9z5+5VPP/20jR8/3kaNGuUWvfnmm6afX375xfLy8myrrbayY4891lq0aOF34RYBBBBAAAEEEEAAAQQQQACBtAlkXAWRQqAbbrjBfvrppwialp177rmmUKGgoCCyfMKECS50KC0tjSxL9M6cOXNs5syZcTdXmKIP/z///HPc9TVZ2K5dOxdE1WSf6rb9888/7aGHHrK5c+dGbVpeVu6W+36//PLLtvPOO5vOVUGSwrD9998/ap/vv//eGV577bWm+/GagqCbb77Z9t57b/vwww+tT58+LhzacMMNXUCy2267xdvNLbviiitsyJAhppAo2BTilJSURP3onBTaBNvtt99uW2yxhQtdHnnkEfvyyy8tPz8/8tOoUSPLzs42hWaPP/64lZVGB3LFxcXOJPh60vF1jffZZx/Ta6hvn77uKRXsxPZJ+yl4CoZO2vjqq682nbeOf8ghh9iuu+5qCuN0vdWuvPJK22CDDaJer24FfyCAAAIIIIAAAggggAACCCCQBoGMqyDab7/9XEXKmDFjrFevXo5s/Ovj3Qd6PXjrrbci1ULaRk3BRZhtrbXWctVKCgxWpylYURCRrnbvvffawIEDTUGRmoK32FBJ26y77rqWm5truq9AJrapSkbH6dSpk51yyinu/osvvug2U6WSQqjff//dunXrFrXrUUcdZQp1BgwY4IKl4Mqrrroq+NCtf+yxx9y19yuWL19uWqbqIN90LIUzq9tU3aNwZ/To0ZaTm+MOt+WWW65SUaZz07mvs/Y6kadUpZGurarLVBX1zDPPuBDp0EMPjWyjKjRVWjVs2DCyjDsIIIAAAggggAACCCCAAAIIpEugVlQQNSxaYC0Kplrr5ZPcrR5X1hQybLrppvbKK69ENhn97Gj3IX2nnXayl156KbJc22y33XbWtm1bV3myxx57WIcOHVzlxtgxY912GkamShEt1xCr4BCiJUuWmIaRde/e3YVMn376aeTYqlrxfVB1iA9JdByFBhoGp6bhRzqGKnT0o+Fdp556qmkIkoZGqR1//PGu/+Xl5TZ79mw74ogjXJii5z377LNNQYiagpwDDjjAbrzxRncOslDlTrJNlTUKbnxIpcobhTy++QBGYYn6pDDG98Vvoz5rOJmGlqkpZFKlj286f7WmTZv6RZFbeSskGTRoUGRZZXcuueQS23fffd2199sovFF4tfnmm/tFodxqSNrzzz/vXku+4ifegTUsTlVBql4Ltnvuucd5KUC86aab3P1gOKRtW7ZsabvssktwN+4jgAACCCCAAAIIIIAAAgggkDaBtAdELZf/bs1LZ1iHlsXWoX22u9VjLa+sKdDRh3NVYSiA+O9//+vmBDrooIPch3qFGIsWLXIf3hUqqCm4UKijyhbNYzP08KG2bNkyu/DCC91cN++99577MK9KGd80r03Hjh3twQcfdB/oNS+RbxpupedQmzx5sqtsOfDAA13ApCFc5513nlungEChlapIFLCMHTvW7dendx9ToKWmoGjYiSsCFlXAaNiShnudfvrpbliU5shRU5Cj4OKLL76wO+64w/bcc0/3vJqHJ5mm59R5qCLKzx0UPI5CLD2n5m9StY/uK5QJtjfeeMM0hE99UVNgo0BO5/Dvf//b7Sd7LYttCr9UlVNdUwjz7rvv2qWXXhq1qYIYhWvBpsBLVWT+Z/r06cHV9sCDD7hrpCBQP7q2wabheaqCUlNlkIaB6XniNYVWqk7T68k3GakqSv2aMmWKM9t9t939am4RQAABBBBAAAEEEEAAAQQQqJUCK9OQNHSvRcEf1jS/yFqvubK6RIO21mibb/OmLLOyivXxmkIfzeHy2muvWePGjd2HcIVGCnwUpiiEUWWLmoakqWkiZt+aNGlizz77rE2cONFNaqxKIf3EzpWjx5dffrmbJ0YBycknn+zm64kXdpx//vmRoW0a8qQ5itRUSbTNNtu4SiY9VrWLKn+69+hu/fr10yLbcccd3XA5zZ/z2Wef2ccff+wmMd5hhx1s/vz5pnl6giGFghs1BV5arnl5NEF0TdvAQQPtu+++c8eQh4aQabJlnbOaAhSFJGussYZ7rPtaFqyG0T6yUTVSsCkoUUCiAK+yuZyC21d1X33ScwcrhRRsKYTS/EXBpvmAfGWXlt96662R66LHftih7qtpTqFg8xNlKzjSxNQK9vSa6tK5i8nLNwVQCikV1gWbwiFdT12Pr7/+2q1q3rx5cBM3Ubgq19QU9On1QEMAAQQQQAABBBBAAAEEEEAgnQLRn+r/wZ6UzZ1geUWLrHWnRnGftXXnRtagYn3JnFWrYxQUqOpFlUMKS1TF0axZMzcBsL6JSh/qFRKpsmOdddZxx1dYoYqWbbfd1vbaay+3TJNba76aNm3auImV9aH+1Vdfjd+fv7/ZS/tU1zQsyQ+tUpiiwEjf5KVvCFMljIaJxWuqRFJT1Ypv+rYrNfU/tmmYkpp/Lr9eAZiaQppgW7hoxWO/XutkpJBClTMKWxRGae4dBTCaS0cVOeqvfv744w/3bWE//PCDO6xCjhdeeMEFHsHn0X31W9VaqsJSFZUPS2K3q+6xKpQ08bUPrfz2Cqb0jXWxQ9dUEabKMv8T++11L7/0sntt6PWhn+eee84f0t1qGKGu0zHHHONCOz1v7969bczYFfNZ+Y1VzaTQSiFdsGm4nYYUqvk5sn74cYWX367/Tv1d5ZSCrNg5n/w23CKAAAIIIIAAAggggAACCCDwTwqkLSAqmfGt5TepqBfKzYp/vhXL8xvnWfG0FVUYsRspFFIFkX40tMs3LdcwLIVHqipSU3VQ3759XbXQ3Xff7Yae+e0VCn3wwQdumJi+dSu2IsVvl+ytQhz9KIRSpcikSZNsk002iXu4Ll26uOXTpk2LrFeFj1rw6+ojKyu5o/NQU5VLsI17bZx76IdEKQjyTfPlDNhrgHuoIWMKYNq3b++qZzTRtH5OOOEEt+yBBx5w26nKRiGJ77eG9u2+++4uWPLHXXvttd3dZIfBKYjRsLvgUDQ/N5IPYvxzhXGrSaP1zWexrUGDBpFF48aNc8GXKpuCTWGYn5xayzWnk4be3XbbbZHhiFo+5LAh7ivudZ+GAAIIIIAAAggggAACCCCAQG0QSFtAZCUFFcOSKgmH/pbRqKXyiu3iNQ0zU4WMfnwQpO38hMda7oeXqQpGH9w1VGjDDTaMfGuXtlcFjL4CXUGGKkV8VY7WhdFUOaPqHw0DU/CjihUNjVJr0aKFu9Xzq1pH1SiqeNKkx5qsWsPUVJGi0CW2UsbtGPhDQ9MUlGk4l85BcwppuJS+nl4VOJp4+eRTTnYTIytsmTdvnqu4Ouecc+yrr74yzcF01dVXufUanqf5kjR0TMPqgj+qyNFwPYU0rlrmhBXVMuqKghUN2dKE1jofVUSpEkfnqcquYB+1vQKqxYsXRya+Xrp0aVQ1lKq5dJzY6qGnnnrKzXUUW72jY65u0+tH4dibb77pXjOaP0kTaQdDSB9a6fUSbBrud/jhh1swTJK7zkvb6tpPmjjJVWdpqB4NAQQQQAABBBBAAAEEEEAAgdoikLaAKKdlNyssKK3SobCgzHJar6hAid2w3w79XPDgh5f59arU0RwwGoLmq0569ezl5hfq06ePdercyQU2qo7JyspyH/41R5BCFVV66AO9mtbpxzc/x05wWfC+X6/tc3Jy3I/uK4jRt6717NnTBSKa60dDyObMmeOGd+kbzzR865BDDnFzKOmbyhRmqX8KhlTh9PDDD+tQUf1xC/7+Izsr24U8CiD8t6fdcsstdsYZZ5huFc5o6JWCND2/WuuKIXP6ynpVWilokdmaa67pJljWNgrYjjvuuL+fYeWNlmnd9ddf785nrwErhuv5LTQHT9euXd38PfoWNs3Ro4mtNdm3gqhgH++8807T/DxnnnmmO6aqpBQu+aZgSDbxgpjYyan9PsHr4JcFb7OqCSXVF4U8+mY6fSOd5jCSv+aDUlPFmoIuDcULtvkL5rtQLbZfqhrTcD1Vjcmux7o93LflKeTSPFq+mit4LO4jgAACCCCAAAIIIIAAAggg8E8LZFVUi6yYzbmSZ9bqYBBSyWY1WqwhX/ogP//BAda2XY7lt85fZf/CeYUVIUqpdT71nVXWVbagrKzMBTLxJpHWPgprNKm1hhEFm85Rc8FUtl9BQUFkDqD1118/uGu19xU8jRo1ygYPHuy2VSXRZptt5o6nQKa0pNSmTZ/mwpngwTSvkIw09CvRpuDGVyVpHz3WuapKR8FZZddRk3vr29uClS9VPafmYVIQpKFl/tvLYrdX/+UW7I/vU+wyv6/6q3mNNDRLwUq8pnNRODNixIi4Q8Hi7aNl6ov6HTthdGXb67Wk/rRq1aqyTSLL9bqaOGGi/e+j/7mwK7Iizh1VhqnSqrqKMO2qKi1VZGl+rTCbAkhNfi4LBXdhNv/e1vsszFbde3t1nkuvKb3+4w0tXJ3jKqzVNVTYm8jrqCbP5d/beq+E2fQe0d8HYVdSqo+aF0zzn9Xk77REzm3ChAmm14eC6bBfd3pvt2jewnJycxLpSsLbpOq9rQ7oPx+q+vs+4U4GNlS1p58jr0f3HqF6ZOJ7W7+3Fy1eFPr7Wn9/qtI4NyfXfZFF4BKs9t1MfG/X9Pd2okh6X+t3sf4+8kPkE923uu0y7b2t80nV7219i6x+v+o/JDU3Z5gtVb+3U/Xe1rmn6vf2b7/95v59q3/PJfpv3ESvRap+byf6/GyHAALpEchNz9OueNZF+R0ta8Yf1rKk3Bq3rpjjJbeioKmkzJbNK7IFc4psUeM1rXOIHdQvqXhNwUll4VC87WuyTJVBGmalb9dS09Cl4cOHRwVC8QKAZD74VRa8VHduyXyo0pCrqj5Eqv/xzqGyPiZqKquRI0cmunnS2ymcq8mHelWr9f1X32qfL+x/JFX7hGyAAAIIIIAAAggggAACCCCAQAICaQ2IinOb2YIma1vJwpmWO2exxlBVTDpUbiV5TWxJk85WkrPqZMEJnFOt2kTz0mg+IH2FfcP8hnbNNdck9ZX0teqk6AwCCCCAAAIIIIAAAggggAACCNQpgbQGRJJUCLSgUTfLblhmWeXFVp6VZ2UVc+rUpaYhZfqhIYAAAggggAACCCCAAAIIIIAAArVRIO0BkUdxoVBWuPNp+GNziwACCCCAAAIIIIAAAggggAACCCBQuUDdKtWp/DxZgwACCCCAAAIIIIAAAggggAACCCBQiQABUSUwLEYAAQQQQAABBBBAAAEEEEAAAQTqiwABUX250pwnAggggAACCCCAAAIIIIAAAgggUIkAAVElMCxGAAEEEEAAAQQQQAABBBBAAAEE6osAAVF9udKcJwIIIIAAAggggAACCCCAAAIIIFCJAAFRJTAsRgABBBBAAAEEEEAAAQQQQAABBOqLAAFRfbnSnCcCCCCAAAIIIIAAAggggAACCCBQiQABUSUwLEYAAQQQQAABBBBAAAEEEEAAAQTqiwABUX250pwnAggggAACCCCAAAIIIIAAAgggUIkAAVElMCxGAAEEEEAAAQQQQAABBBBAAAEE6osAAVF9udKcJwIIIIAAAggggAACCCCAAAIIIFCJAAFRJTAsRgABBBBAAAEEEEAAAQQQQAABBOqLAAFRfbnSnCcCCCCAAAIIIIAAAggggAACCCBQiQABUSUwLEYAAQQQQAABBBBAAAEEEEAAAQTqi0DaAqKSJXMTMi4tWpbQdmyEAAIIIIAAAggggAACCCCAAAIIIJCcQNoComWzfkmoxwVzpia0HRshgAACCCCAAAIIIIAAAggggAACCCQnkLaAKK9p2+R6zF4IIIAAAggggAACCCCAAAIIIIAAAqEKpC0gCvUsOBgCCCCAAAIIIIAAAggggAACCCCAQNICBERJ07EjAggggAACCCCAAAIIIIAAAgggUDcECIjqxnXkLBBAAAEEEEAAAQQQQAABBBBAAIGkBQiIkqZjRwQQQAABBBBAAAEEEEAAAQQQQKBuCBAQ1Y3ryFkggAACCCCAAAIIIIAAAggggAACSQsQECVNx44IIIAAAggggAACCCCAAAIIIIBA3RAgIKob15GzQAABBBBAAAEEEEAAAQQQQAABBJIWICBKmo4dEUAAAQQQQAABBBBAAAEEEEAAgbohQEBUN64jZ4EAAggggAACCCCAAAIIIIAAAggkLUBAlDQdOyKAAAIIIIAAAggggAACCCCAAAJ1Q4CAqG5cR84CAQQQQAABBBBAAAEEEEAAAQQQSFqAgChpOnZEAAEEEEAAAQQQQAABBBBAAAEE6oYAAVHduI6cBQIIIIAAAggggAACCCCAAAIIIJC0AAFR0nTsiAACCCCAAAIIIIAAAggggAACCNQNAQKiunEdOQsEEEAAAQQQQAABBBBAAAEEEEAgaQECoqTp2BEBBBBAAAEEEEAAAQQQQAABBBCoGwIERHXjOnIWCCCAAAIIIIAAAggggAACCCCAQNICBERJ07EjAggggAACCCCAAAIIIIAAAgggUDcECIjqxnXkLBBAAAEEEEAAAQQQQAABBBBAAIGkBQiIkqZjRwQQQAABBBBAAAEEEEAAAQQQQKBuCBAQ1Y3ryFkggAACCCCAAAIIIIAAAggggAACSQsQECVNx44IIIAAAggggAACCCCAAAIIIIBA3RAgIKob15GzQAABBBBAAAEEEEAAAQQQQAABBJIWICBKmo4dEUAAAQQQQAABBBBAAAEEEEAAgbohQEBUN64jZ4EAAggggAACCCCAAAIIIIAAAggkLUBAlDQdOyKAAAIIIIAAAggggAACCCCAAAJ1Q4CAqG5cR84CAQQQQAABBBBAAAEEEEAAAQQQSFqAgChpOnZEAAEEEEAAAQQQQAABBBBAAAEE6oYAAVHduI6cBQIIIIAAAggggAACCCCAAAIIIJC0AAFR0nTsiAACCCCAAAIIIIAAAggggAACCNQNAQKiunEdOQsEEEAAAQQQQAABBBBAAAEEEEAgaQECoqTp2BEBBBBAAAEEEEAAAQQQQAABBBCoGwIERHXjOnIWCCCAAAIIIIAAAggggAACCCCAQNICBERJ07EjAggggAACCCCAAAIIIIAAAgggUDcECIjqxnXkLBBAAAEEEEAAAQQQQAABBBBAAIGkBQiIkqZjRwQQQAABBBBAAAEEEEAAAQQQQKBuCBAQ1Y3ryFkggAACCCCAAAIIIIAAAggggAACSQsQECVNx44IIIAAAggggAACCCCAAAIIIIBA3RAgIKob15GzQAABBBBAAAEEEEAAAQQQQAABBJIWICBKmo4dEUAAAQQQQAABBBBAAAEEEEAAgbohQEBUN64jZ4EAAggggAACCCCAAAIIIIAAAggkLUBAlDQdOyKAAAIIIIAAAggggAACCCCAAAJ1Q4CAqG5cR84CAQQQQAABBBBAAAEEEEAAAQQQSFqAgChpOnZEAAEEEEAAAQQQQAABBBBAAAEE6oYAAVHduI6cBQIIIIAAAggggAACCCCAAAIIIJC0AAFR0nTsiAACCCCAAAIIIIAAAggggAACCNQNAQKiunEdOQsEEEAAAQQQQAABBBBAAAEEEEAgaQECoqTp2BEBBBBAAAEEEEAAAQQQQAABBBCoGwIERHXjOnIWCCCAAAIIIIAAAggggAACCCCAQNICBERJ07EjAggggAACCCCAAAIIIIAAAgggUDcECIjqxnXkLBBAAAEEEEAAAQQQQAABBBBAAIGkBQiIkqZjRwQQQAABBBBAAAEEEEAAAQQQQKBuCBAQ1Y3ryFkggAACCCCAAAIIIIAAAggggAACSQsQECVNx44IIIAAAggggAACCCCAAAIIIIBA3RAgIKob15GzQAABBBBAAAEEEEAAAQQQQAABBJIWICBKmo4dEUAAAQQQQAABBBBAAAEEEEAAgbohQEBUN64jZ4EAAggggAACCCCAAAIIIIAAAggkLZCb9J5p2vHOO++0qVOn2vXXXx/Vg/vvv9+Kiors4IMPjlp+wQUXWM+ePe3II4+MWl7dg+HDh9uMGTPs0UcfXWXTBQsWuOc58cQTbfDgwausr8mCsvKymmye0LY//vijjRw50i6++GLr3LlzZJ+ysjI77rjjbOjQoda/f3+3/O2337YnnnjCFi5caBtvvLGddtpp1rp168g+uiPDrKwsu+aaa6KWBx+cddZZtt5669mwYcPszTfftMcff9yWLFliffv2dcds1KhRcPOo+0899ZS9//77dt1110WWX3TRRfbXX39FHvs7G220kem5fPv111/tsMMOs48++sjOOeccW7x4sV8Vub3rrrvsnXfesSeffNJdT52Lb4WFhXbSSSfZ0Ucfbdtvv71b/MUXX9h9991ns2bNsk022cTOPPNMa9mypd8lcjt+/Hh75pln4ro8/fTTpvWjRo1y28tEP7/88ovl5eXZVlttZccee6y1aNEicjzuIIAAAggggAACCCCAAAIIIJAugYyrIFIIdMMNN9hPP/0UMdOyc8891xQqFBQURJZPmDDBhQ6lpaWRZYnemTNnjs2cOTPu5gpT9OH/559/jru+JgvbtWvngqia7FPdtn/++ac99NBDNnfu3KhNy8vK3XLf75dfftl23nln07kqSFIYtv/++0ft8/333zvDa6+91nQ/XlOQcvPNN9vee+9tCkZ23XVXUxjVvXt3F54MGTIk3m5u2RVXXGFar5Ao2BTilJSURP3onL755pvgZnb77bfbFlts4UKXRx55xL788kvLz8+P/CiYys7ONoVmCq3KSqMDueLiYmfiX0+ff/659e7d2z2vzkNGO+64o2m7YJsyZYoddNBBLgCK9/q6+uqrbbfddnP7HXLIIc5EYZyut9qVV15pG2ywQdTrNXh87iOAAAIIIIAAAggggAACCCDwTwrUigqi7PISa1C82HRblpVrRXnN3G08iP32289VkIwZM8Z69erlNhn/+nhXAaMHb731VqRaSNuoKbgIs6211lquWkmVIKvTFKwobEpXu/fee23gwIEuBFEfFLzFhkraZt1117Xc3FzTfQUysU1VMjpOp06dXJh06qmnRrZTCHLUUUfZokWLrHnz5lG7arlCnQEDBtiHH34Yte6qq66Keqz1jz32WFT10PLly90yVQf5pmMpnEm2KahSQPTggw+6QyhAUwXad99954IoLVRYdOCBB9qGG27oKpdin+uTTz5x1UeqLlNVlKqMVL106KGHRjZVFZrCqIYNG0aWcQcBBBBAAAEEEEAAAQQQQACBdAmkvYKoUeFsa734F2tts61t/gJ3q8daHq9169bNNt10U3vllVciq0c/O9oFFDvttJO99NJLkeXaZrvttrO2bdu6ypM99tjDOnTo4Co3xo4Z67bTMLJ99tnHLdcQKw0t8k1DpDSMTJUwCpk+/fRTv8qFBb4PqjTxIYmOr1BBw+DUNIRJx1CFjn40vEsBiiptNDRK7fjjj3f9Ly8vt9mzZ9sRRxxhOk8979lnn20KQtRUzXLAAQfYjTfe6M5B26hyJ9mmyprff/89ElKp8kYhj28+gNFQKPVJAY3vi99GfdbwPg0tU1OVjqqNfGvWrJm7q+1im7wVkgwaNCh21SqPL7nkEtt3333dtfcrR48e7cKrzTff3C9a7duOHTvavHnzXAWUDqb7akEXVasp4Lntttvcutg/7rnnHuelAPGmm25y94PhkLbXkLVddtkldlceI4AAAggggAACCCCAAAIIIJAWgbRWECkEalY6xzqs1cSymqzsSvnSEps5bU5FuhLfRIGOhujoQ7qGEP33v/+1O+64w815o/lkFGKoykNDejQXj5qqYzQU6dJLL7VbbrnFhh4+1KZPn24XXnihm+vmvffec6GOhmf5pnltFP6omuSyyy4zzUukoERNw61UFaM2efJkV9ly6623urDn5JNPtvPOO88Nm1JAoNDq+eefd1VHqipR2NSndx+3n+alUVDUtUtXdyxVwChM0XCvr776yg3RUlXPww8/7IIcHUfVPDrfZ5991j2v9tH8PzVtw04cZoP2HmSqiNIwKM1PpOoZ3xRiqcLJz9+kuYgUygTDjjfeeMM0xGrPPfd0u2mOIN90Hqo40hCteHPtKPxSU0hUVdN1fPfdd93wseB2CmIUrgWbAi9VkfmmKrNguPPAgw+4+ZT8+tjAS+cqY80RpOBJ105zHq2xxhpulxdeeMF0nb/++ms3jM0fx9/KS1VRmmtIw9D0ePfddveruUUAAQQQQAABBBBAAAEEEECgVgqsTGX+4e6VFy22JoWzVgmH1A2FRR06N7KyKbOtvHDJKj1TJYkCotdee80aN27sPoQrNFq2bJmbcHjs2LEuZNGOGpKmpomYfWvSpIkLVyZOnOgmNValkH40Z0yw6fHll1/u5onRcCgFP5qvRxVJse3888+PBCka8qQ5itRUSbTNNtu4SiY9VujgqoN6dLd+/fppkQtQFGRo/pzPPvvMPv74YxdQ7LDDDjZ//nzTPD0KQ3xTcKOmwEvLNS9PMgHRwEED3dApHUNhk4aQKQjTOaupmkpDx3w4ovtaFgyItI9sVI0U22SikEvntDpNfdJzByuFFNBp3qDY+Y00lMtXduk5Feb4gEuP/bBD3VfTPEfBpnmrdM3WXHPNyPxACuj02tKcVMccc4xponRVsfm5nIL7KxxSIKbroRBJLXZonYI4Va6pKejT64GGAAIIIIAAAggggAACCCCAQDoF0hcQzfjO8hrmRlUOBSEUEuU1zLGSv74yW3tlVYu2UVCgqhdVDqmpIkdDmRT86JuoNOeLgiNVs6yzzjpuG1VzaDiWwgpVmahpcmvNV6MhVH369HHDlVRdtNdee7n1wT/8N3tpn+qaJiLW0DI1hSkKr1SFoqomVcIotIjXVImkpnl7fFMli5r6H9v8N2v55/Lr5aAWO7/RwkUr5jvy67WNjBRSyEZBj8IoTfatChjNpaP1Gtam9scff7gw6ocffrC1117bhRy+osZtEPhDoY7Co3HjxjnXwKoa3VWFkuYfiq0y0rH1jXVNmzaNOp4qwqqag+jll162nNycyD4KBv0wOC3UN7UpePNBkkIhBWS77767qwySqYKp008/3VWwaR99W5yCKg0Z03A7+an5ObJ++PEHUxjnW/+d+tsvv/7iQk55ExB5GW4RQAABBBBAAAEEEEAAAQTSJbBq2cc/1JPygoWWm1f102t9ecGCuD1SKKQKIv1owmDftFxDhBQeKZhRUwigr1vXV6DffffdbuiZ316VHh988IEb7qVJh2MrUvx2yd4qxNFPmzZtXBAwadIk99Xp8Y7XpUsXt3jatGmR1ZocWS34dfWRlZXc0XmoBYda6fG418bpxoU+ug1+M5fmyxmw1wAtdkPGFMC0b9/eVWT179/f9HPCCSe4ZQ888IDbTkPvVNnj++0WVvyhCiRVeClk8QGXX1fTWw0J1BC6LbfcMrKrnxtJczuF3VTxFHwuzSmlkFHVXaoeUgCliiCFSgoh1XRfk01rSKImHtcwQjXN6aShd5qryA9H1PIhhw1xoaTu0xBAAAEEEEAAAQQQQAABBBCoDQJpqyDKatTSSoqjv3I8FqSkqMzyGrWKXewea5iZ/0YtHwRphSY81iTCqvTww8tUMaQP7pqfaMMNNrSb/nNT5JiqgNG3VIfbcYgAAEAASURBVKlqR/PvKCAIs2mYkap/NAxMQcpzzz3n+qjqEj8vj76CXf1V5YrCCPVfX+muvqgiRdUrsZUysX1UZZTmWVI1kEINTRqtb9BSSKHKlm+//daFG7qvAESTL+u5NNxJcyApwLjq6qvctpo7SBNSK4DRsLpgk6WGmalqRn274/Y7gqsjVTFaX1BQYKoAUtt6661d5U2wjwqotI2fB2jp0qWu8kp9Vnv11VddFVPsEDV905i+WU1eYbdtt93WzfekoFGvC80Fpcqufjv0c1VA+mp73zTETEGa5pvq0L6DHX3M0Xb44YdbgwYN/CZ211132WabbeZeW6pO2nyzzW15wXI3SXlkI+4ggAACCCCAAAIIIIAAAgggkGaBqkt4Uti5rDU2teLCUitfEj0HjH/K0orlxUWllts5fgigD+wKWPzwMr+fKnU0B4yGoPlKkF49e7n5hTSMrFPnTi6wUXVMVlaWG7qkOYJU5aNKD32gV9M6/fjm59gJLgve9+u1fU5OjvvRfQUxmq9GYYMCEc31ozBKcxlp+Ja+8UxDuDRJtOa50TeVKcxS/xQMqcJJE1SrBZ/PLfj7j+ysbBcmKXzy356moXJnnHGGm5Bbw+5U+aIgTc+vpiFzCthUaaWgRWaad0dz6GgbBVYKj2Kblmnd9ddf785nrwErh+Npgm8NLVPT19RrDif/ozmZFHgF+6i5fFSNc+aZZ7pjqkrKTwKuY6gSSTbBibO1PN7k1FquFrwOK5ZE/5mVvfKaRq9Z8UjXX0MY9drRUDzNo6RQKzhELN5+8xfMd6Fa7KTZ3SrmF9J8SZtssonz7LFuD9O3tynkUpWVXgM0BBBAAAEEEEAAAQQQQAABBNItkFVRLbLq948HeqXVlQUTgc1qdFdDvgrnTLR5r42wpiWzrX3FhNQ5TVcWMykcmjVtuS3JbWed9rnamnTqmdDxy8rKXCATbxJpHUBhjYYFaThQsOkcNRFxZfupysXPAbT++usHd632voKnUaNGRYYdqZJIFSU6ngKZ0pJSmzZ9mrsfPJjmFVLYoaFfiTYFN74qSfvosc5VQ+sUnFV2HRVM6ZvRgpUvVT2n5mFSkKSKKP/tZVVtH1wX28fYdZq8WRVEClbiNZ2L5u0ZMWKE+wa7eNvEW6ZrqH7HThgdb1stU5in12mrVvEr2IL76XU1ccJE+99H/3NhV3Bd7P3Zs2e7fldXEab9VFmlSbSDcyTFHi+ZxwogNfm5LDp27JjMISrdR2Z63frhd5VuWMMV1b23a3i4qM31mtLrX9+IGGZTWKtrqLA3kddRTZ7bv7d9tV1N9q1qW71H9PeB/t4Ku+lbJxW61uTvtET6oInl9fro2rVr6K87vbdbNG8RNW9ZIn2qbptUvbf1vPrPh6r+vq+ub/HW6+9DP0dej+49QvXIxPe2fm8vWrwo9Pe1/v7UEPPcnFzrXvFFFmG2THxv1/T3dqJeel/rd7H+PoodIp/oMSrbLtPe2zqPVP3e1jcF6/er/kNSc3OG2VL1eztV722de6p+b//222/u37f691yi/8ZN9Fqk6vd2os/PdgggkB6BlalMGp5/ef6KbwMr/2OW5TXIqZiTKKti2Fm5qxxa2qC9+fVhdc1PNB17PAUnlYVDsdvW9LEqgzR3jb5dS01faz98+PCoQCheAJDMB79gOBTsZ3XnlsyHec37lMyHyMr6GOxvVfdlNXLkyKo2CWWdPsTW5EO9Ko76/qtvtc8d9j+Sqn1CNkAAAQQQQAABBBBAAAEEEEAgAYG0BkTqn0KgwgatLa9ksWVX/C9lWXaeFTdtZmUVw6bqQtNwKM0HpEmOG+Y3dN+SlcxX0tcFC84BAQQQQAABBBBAAAEEEEAAAQRqp0DaAyKxKAwqzGtRO4VC6JWGlOmHhgACCCCAAAIIIIAAAggggAACCNRGgbpRplMbZekTAggggAACCCCAAAIIIIAAAgggkCECBEQZcqHoJgIIIIAAAggggAACCCCAAAIIIJAqgbQFRHO/HpvQOc2f8L+EtmMjBBBAAAEEEEAAAQQQQAABBBBAAIHkBNIWEBUvnJFQj5fP+zOh7dgIAQQQQAABBBBAAAEEEEAAAQQQQCA5gbQFRKWlxQn1ePn86Qltx0YIIIAAAggggAACCCCAAAIIIIAAAskJpC0gSq677IUAAggggAACCCCAAAIIIIAAAgggELYAAVHYohwPAQQQQAABBBBAAAEEEEAAAQQQyDABAqIMu2B0FwEEEEAAAQQQQAABBBBAAAEEEAhbgIAobFGOhwACCCCAAAIIIIAAAggggAACCGSYAAFRhl0wuosAAggggAACCCCAAAIIIIAAAgiELUBAFLYox0MAAQQQQAABBBBAAAEEEEAAAQQyTICAKMMuGN1FAAEEEEAAAQQQQAABBBBAAAEEwhYgIApblOMhgAACCCCAAAIIIIAAAggggAACGSZAQJRhF4zuIoAAAggggAACCCCAAAIIIIAAAmELEBCFLcrxEEAAAQQQQAABBBBAAAEEEEAAgQwTICDKsAtGdxFAAAEEEEAAAQQQQAABBBBAAIGwBQiIwhbleAgggAACCCCAAAIIIIAAAggggECGCRAQZdgFo7sIIIAAAggggAACCCCAAAIIIIBA2AIERGGLcjwEEEAAAQQQQAABBBBAAAEEEEAgwwQIiDLsgtFdBBBAAAEEEEAAAQQQQAABBBBAIGwBAqKwRTkeAggggAACCCCAAAIIIIAAAgggkGECBEQZdsHoLgIIIIAAAggggAACCCCAAAIIIBC2AAFR2KIcDwEEEEAAAQQQQAABBBBAAAEEEMgwAQKiDLtgdBcBBBBAAAEEEEAAAQQQQAABBBAIW4CAKGxRjocAAggggAACCCCAAAIIIIAAAghkmAABUYZdMLqLAAIIIIAAAggggAACCCCAAAIIhC1AQBS2KMdDAAEEEEAAAQQQQAABBBBAAAEEMkyAgCjDLhjdRQABBBBAAAEEEEAAAQQQQAABBMIWICAKW5TjIYAAAggggAACCCCAAAIIIIAAAhkmQECUYReM7iKAAAIIIIAAAggggAACCCCAAAJhCxAQhS3K8RBAAAEEEEAAAQQQQAABBBBAAIEMEyAgyrALRncRQAABBBBAAAEEEEAAAQQQQACBsAUIiMIW5XgIIIAAAggggAACCCCAAAIIIIBAhgkQEGXYBaO7CCCAAAIIIIAAAggggAACCCCAQNgCBERhi3I8BBBAAAEEEEAAAQQQQAABBBBAIMMECIgy7ILRXQQQQAABBBBAAAEEEEAAAQQQQCBsAQKisEU5HgIIIIAAAggggAACCCCAAAIIIJBhAgREGXbB6C4CCCCAAAIIIIAAAggggAACCCAQtgABUdiiHA8BBBBAAAEEEEAAAQQQQAABBBDIMAECogy7YHQXAQQQQAABBBBAAAEEEEAAAQQQCFuAgChsUY6HAAIIIIAAAggggAACCCCAAAIIZJgAAVGGXTC6iwACCCCAAAIIIIAAAggggAACCIQtQEAUtijHQwABBBBAAAEEEEAAAQQQQAABBDJMgIAowy4Y3UUAAQQQQAABBBBAAAEEEEAAAQTCFiAgCluU4yGAAAIIIIAAAggggAACCCCAAAIZJkBAlGEXjO4igAACCCCAAAIIIIAAAggggAACYQsQEIUtyvEQQAABBBBAAAEEEEAAAQQQQACBDBMgIMqwC0Z3EUAAAQQQQAABBBBAAAEEEEAAgbAFCIjCFuV4CCCAAAIIIIAAAggggAACCCCAQIYJEBBl2AWjuwgggAACCCCAAAIIIIAAAggggEDYAgREYYtyPAQQQAABBBBAAAEEEEAAAQQQQCDDBAiIMuyC0V0EEEAAAQQQQAABBBBAAAEEEEAgbAECorBFOR4CCCCAAAIIIIAAAggggAACCCCQYQIERBl2weguAggggAACCCCAAAIIIIAAAgggELYAAVHYohwPAQQQQAABBBBAAAEEEEAAAQQQyDABAqIMu2B0FwEEEEAAAQQQQAABBBBAAAEEEAhbgIAobFGOhwACCCCAAAIIIIAAAggggAACCGSYAAFRhl0wuosAAggggAACCCCAAAIIIIAAAgiELUBAFLYox0MAAQQQQAABBBBAAAEEEEAAAQQyTICAKMMuGN1FAAEEEEAAAQQQQAABBBBAAAEEwhYgIApblOMhgAACCCCAAAIIIIAAAggggAACGSZAQJRhF4zuIoAAAggggAACCCCAAAIIIIAAAmELEBCFLcrxEEAAAQQQQAABBBBAAAEEEEAAgQwTICDKsAtGdxFAAAEEEEAAAQQQQAABBBBAAIGwBQiIwhbleAgggAACCCCAAAIIIIAAAggggECGCRAQZdgFo7sIIIAAAggggAACCCCAAAIIIIBA2AIERGGLcjwEEEAAAQQQQAABBBBAAAEEEEAgwwQIiDLsgtFdBBBAAAEEEEAAAQQQQAABBBBAIGwBAqKwRTkeAggggAACCCCAAAIIIIAAAgggkGECBEQZdsHoLgIIIIAAAggggAACCCCAAAIIIBC2AAFR2KIcDwEEEEAAAQQQQAABBBBAAAEEEMgwAQKiDLtgdBcBBBBAAAEEEEAAAQQQQAABBBAIW4CAKGxRjocAAggggAACCCCAAAIIIIAAAghkmAABUYZdMLqLAAIIIIAAAggggAACCCCAAAIIhC1AQBS2KMdDAIF6J/Dmm2/aCSecYLNnz4469y8+/8IuvPBCKy8vj1r+yCOP2AUXXBC1LJEHL7zwgq233nqVbnriiSfaJZdcUun6RFfE9jfR/ararqSkxI4//nj74IMPVtnsuuuus3vuuSeyfPLkyfbvf//b9tlnHzv77LPtxx9/jKzzd2R4zDHH2JIlS/yiVW6ffvppt41fUVZWZldccYXdeuutflHU7bx582zYsGE2fvz4qOXxHjz11FN20kknRVYtXLjQPZf6FPvz6quvRrbTnQkTJli/fv2suLjYzjjjjFW21/4FBQWm15Xux16PwsJCtzzWUuc3cuTISs9Pz61zO/bYY91z67Ganuumm26yAw44wI444ggbN27cihWBP0tLS+3www+3b779xi0tLSk1XYPTTjvNBg4caEcddZTJRMtpCCCAAAIIIIAAApkpQECUmdeNXiOAQC0SaNu2rd1///327LPPRvXqvvvvc8HH//73v6jlCo2mTZsWtSyRBwpDFC5U1hQY/N///V9lqxNePmjQIHvyiScT3j6RDRVyKFD47rvvVtl8zJgx9tZbb7nlU6dOtS222MJk1qNHD/vkk09sm222sWXLlkX2U7AyfPhwe+ihh+yll16KLI+9c/XVV9tuu+3mFstOgdPll19u7777buymNmniJNtqq63s3nvvtZ9++mmV9cEFCpmGDBniAhG/PCcnxxSiKAjzPwqN1McZM2b4zdzt3XffbZtttpnl5eU5ky+//NLy8/MjP40aNbLs7Gz74Ycf3P5lpWVR++v8ddxgP/35XX/99fbee+9Fbe8fTJkyxQ466CAbNWqU66Nfvu+++9oNN9zgwkeFRXvuuaeNHTPWr3a377zzjs2fP9823WRT03H69O3jQqFff/3VunTp4s5RJgcfcnDUfjxAAAEEEEAAAQQQyByB3MzpKj1FAAEE/lmB8pJCK10wxbLyGllOi66VPrk+7K+11lo2duxYO/nkk912ixYtinxQf+aZZ2y77bZzy7/66iubPn266UN52O27b7+z7JzVz/1VwZOupqofhSUffvihu1U/fvvtN2vcuHGkSy+//LLNmjXLVa7IVsFEbFOwpG0GDx5sqgzafffd3W3v3r1jN7XPP//cBgwYYJtvvrnbZ5UNAgtUKaOgS9urj741bdrULfePdTtixAgXcB122GGRxcuXL3fB0iuvvBJZpmMpzEq2Bc9Pr8V4TaHSgQceaBtuuKF99NFHkU2++eYbe/31153Blltu6ZZvv/32NvrZ0TZw0MDIdrou2l/tuOOOcwGnXsvB51NYpPOjIYAAAggggAACCGSmQEYEREun/5yQbnl5mRUtWWxLi+YktH2iGxUVF1r5wllu86XTo4eKJHqMyrYrKyuxwor/GV+6vHllmyS9vGDpUivPz7OS3AZJHyPejsWlxVZY8b/MS5c1i7c66WWuvwvnWWleri1tUJD0ceLtqGtYUlxiS5c0ibd6tZYVLl5UUd3QuOJ//MN9OxUWFVpZRUVC9uKVH4xXq6N/71y8cL6VL15S4dHYlpYvDOOQkWMUFlR8OMzOsvIGDSPLwriTqve2+la4fJkVV1R/lDXIj3S1vHCxFX79mBVO+dKyc3Pddcht2sqa/esoa7TR4Mh2wTuqTrntttvcB2RVgGiYjiqLNETo8ccfd8N+VBWiahm1PfbYw+bMmeMqYd5++2035EfD1C677DK3XlUeGnqlKhR9YL/55pvdcv2hIVIazqMP4xqidOaZZ7p1x59wvKvmuPLKK+2iiy6y3Iq+q3pFAUCzZs3cfv3793fb6hjqrypP+vTpY926dTMFBAouVJly7XXX2qOPPeqGLK2//vou7Hjuuedcf7beemv7z3/+Y2uuuaY7vs5FQ+ZuueUWU7ikwENDxho0qPnfPTJSsKPj9OrVy/V17bXXdrf+D1X5KPDROaovCjm6do0O8PT8Gi6lPihc0nmfe+65duqpp5qGaAVbixYt3PILL7jQOnXuFFy1yv2NN97YhSkKlYIBUeyGum4atqVqIwVevo0ePdq6d+9um2yyiV+02rfB89NrKHZImp5A575gwQL3ugmGZOuuu66r6tpoo40i/ZBPy5YtI49///13d643jbzJvvzqSzf0Ta/PYDikjasa/hg5GHcQQAABBBBAAAEEaq1AuJ9oU3Sak8at/GCUyFOEGw9FP+Ok76Mfh/VoRfwU1tH+mePMTNHTFFUcd1KKjp2qw0bPPJOqZwn3uHKeH+4hU360VL63feezy0us1bLJ1qRJlnXqURGCNlhRkVO8sMDmvne7lS6eZU23HuY3j9yqUkWBy1tvvuUqLxQE7bjjjm5eFw3f0fAvzTuj5ZqzRR/q77vvPhfyPPbYY269Ao+dd97ZhR2aJ+aBBx5w4ZAqZrKysiLP9e2337owRsOrzjrrLBfI6MO5hv5o2JLaX3/95YYhKXB68MEHTaGR5pDR0LaPP/7YzfHz8MMP26abburmBlJVk8IFVcI8//zzLsDqt0M/W2ONNVwQo6qZO++80x37qquusl122cUNgVJligKaa6+91v2oD5qbR8O1dLyatoMPPthuvPFG22CDDVyVlUIeVVv581c1kebm0XA+PYeGoamKSKa+KZxRf3/55Re3qGHDhqahV5U1hSSJzt2kOZHUFBBV1RSgKRjSvFDBpuBKoWGwKYDxQ+y0XMFYp04rg6oHHnwgcv5aH1ulU935ae4qBYJff/11VFilY+l16MMhXV9tq+MrUPJNr0P5tmjZIjKHlEJBGgIIIIAAAggggEDdEsiIgKhukXM2CCBQmwWaFv7lwqGWXRpFdTOvRZ6tkZ9t0z993PLX2c7yOmwUtX67bbczVaKMGTvGttp6K3v//fddMKMAQx/AVTnSs2dPFywosFHTRMy+aQiagqQvvvjC2rRp4xYrzFFAE/ywrhXaX6HAv/71L1epo4AmXvWG5t/RnDtqqhRS0KKqJc3zo74qfFHwooofzTGj/qkqSG3zzTaPDDFSkKUKIU1irNauXTsXfmkfX+Wj4EOBjZruq0/JNFUCff/996Z5ejSsScGbwihN9KzwS31p3769m09Ix9ewJ4VTd9xxR2QYmsIhhXPxTJLpU0330Vw9CmQuvfRSd538/jovVWdpHqBge/LJJ13lll+mfY888kj/MFJ15hdojqNEmwI1BVIKfxQG/vxz5RW5Cvfmzp3r5lJSFZeuhQJABYma/FpN4ZuaKtJ802tKQ+/UFDDqutEQQAABBBBAAAEEMk9gxX+NZ16/6TECCCAQukC2lVmDokXWsv3KIUFRT9Iwx5q1yreCn8dFLdaDnNwc23vvvd08RKrA0fAyH5goxNDwLFUCqWk7NU28rGFku+66q2kYlz58FxUVudBFc9KoCsUHOcFJmt3OFX/4YUCxQ6b8+uCtQh01HV9VSk2aNHFBhUIMDQ0LzpMT3E/z26hfCo98UzClpvAhXtNzxfZJw93UNDdTbNMy9cc37a9+6dvLFGyoYkjz5iis0NAmtUMPPdQFVn5CZg25800ThsdW7vh1/8StKqBUPRT8ljM9r4bGqUJKIUqwadJyDf/yP8FwSNu9/NLL7nWlOa70o9dSok190PVTMHX66afbNddc43ZV5ZmGNgabAkpVGal6SxVnanrN6jXYt29f99hXGyns8k3nozBS28X7BjS/HbcIIIAAAggggAACtVuAgKh2Xx96hwAC/6BAdlmRZVXMoeSHlcV76vyG2VYya8XQpdj1qtDRUC19uFbVi77ZSk0VI6rI0Hw0qhRSeKSm6hh9ANfwJn372DrrrOOW6w+FBtpHlUcagqZhPmE1DStSIKPqJgUZCh0qGw7WqlUr97TBb13z4UDnTp0T7pIqlTTvjuZDCjZVreh4PnhQCBRs++23n3uowElWMtGcS5pTSD8KJjSHkqpc1BQkaRvZpqOpmuauu+5yVV+ai8o3DdvSdfyngytVD+m11Lx5c1f1o2uvpseqQlO1Vey11+vCf1ueQq3jjz/en4arWlMFl4YZ+qbjKHzy3xjnl3OLAAIIIIAAAgggkFkCDDHLrOtFbxFAIIUC5VYxwXY189CXlVVs0DD+xOF+XhaFRIMGDor0VNU3CkAUhPgJpbVS37SleYd22GEHVyWjiZnVZs+e7YYVac4gVfvoA3ls1YnbMMk/FGIoANA8OuqbKnQUwGgybDU93w8//uCqSVTZpH5owmUN21JTdY+26bdjP1f14hbG+UNDrRSIaGJohUPnn3+++8YxDW3TkCRVJ/lj6Zux1GShIXb6GnuFGArV9FyqWlLQpjBIx/FN56H5ehRyyFPD29Tf2Amy9fXtCp/8V9AvXrzYBSR+ziYNwfOTO8tC6/W19Z99/pk7dw1h69ChgzuGjuXnAdJ2OoZCEjVV4Sh0899m5/upCifNdbTFFltE9vXrwrgNnp/Ow5+fzIJD2jTETKGPhh4qwNI56lp07tzZjj76aPv2m2+doYIlvR712lCwFRla1rSZq+pSVZxe77pOmnRbAaLmxKIhgAACCCCAAAIIZK4AFUSZe+3oOQIIhCxQml0xtCw7x4oXRlexBJ9m6RKz/LW2Di6K3Fd1hibzVaCx9TbR2/j5e4KVLQqL9AFbwYPm0VE4okqbmTNnugmvFXIoLFHV0dChQ6MmKo48acWd7KwVf5XrG8B803H045tfp2UaGqZqJQ2BU7ihoEgVT/6r14cNG+aGIikQ0jpNvq1vONNXpOvnjz/+cNU8wXlofLWUnk/PpR+dhyaT/vTTT103ZKBKH83PpOBLIYOGk7377ruR4XK33367mwdHE3rra+cnTpxoL774ogsoxo8fb6eccoo/pcitfPRNZxqOpjAjWPHiN9KQNAVOGv6niindv/jii/1q9+1vWqbqo/POO8+t1xxLGnKlIV2at0lNz6HtdO0Umui+rpuagjfNkaTj+kodt6LiDwVX8fql9f7a+G1jb11VW+zCmMf+/DQJerzzi9k88lChoM5J35K31lpr2aC9B7nqKw1FU5B0yCGHRObE8jvpOqpSSyGUrmO3bt1s2223tT///NPNHeW34xYBBBBAAAEEEEAgswSyKv6nscr/L9fq4IeMME5P/1M7efSKb4IJ43gcAwEEEAhLoGHRPGtWPNM6rlVRJVQx51CwLZtVYMuKm1vbo18MLo57X4GCqkwUFvlhWvE21N+HquLQ/C2xLbY6xa9XMKGKFVWqJNP22WcfNxm1whg1VcxocmqFVAqwNLeR5jzKb5Dv5lbyz+ErcIJDp/y6ym41v5BCFM2vo2FtvmJH56DjxFb6+ONoriRV+8QGLX69v1UFUVlZmZtQWXMW6dviglVafrtkblXh1KJ5C1tSkQrGuz6JHlPXUZVQI0aMcOes14XOLRiwJXqs6rZTSKVQsaa/t/W7XpVrbVq3iVxz9VnBk8I3X93Wo3uPyHr1pbSk1ObMneNeM8m8HnXt5OyHXVZ3fjVZL3fNfVWT12six9c5L1q8qMr3dSLHid1GfxeoEis3J9e69+geu3q1Huv9pPe0n7dstQ4Ws3PseztmddIPVRWnfuvvjzCbXm96revvoy5duoR5aPd3fm17b1d3gnrdKaSu7u/a6o4Tu16VtHoPtm7d2v1HQOz61Xmcae9tnevq/t6uzEv/4aP3SceOHUN/r6TqvV3ZubAcAQRqh0BGDDFrukZi/1BS1FVaXGK5DcI9LQ0pKSgodFesceMVwwjCunzqc0nFPzbz8qI/iIZxfH3ozKqoLMjWnCohtvIKj7LSMssJuc/qb2FhsfuHSsOGDULssVV8gCx3w0dyclZWWIT1BCUVr7mcig8hgWKNUA6tD056fYTd56KiEvfhVBUf+fl5ofTVH6S04nWhV1t2yM6pem+r367PFZ1eWcXR3Yrn/mTTJ0+wZi3zLb9RxdTVFa8fVQ6V5za3loMq/7p071CT26qGjqUiQFDfVLmjYVi//vqr+8e75j/SELj9998/Mswq3ocFhTs+4En0HCv7cFdd4KLgqLLwqLLn7t27t+kn7FZdX6t7Pl1H/y1g1W2brvUKlBRmBttll13mHsbOCxXcRpOzqwKOhgACCCCAAAIIIJD5AuEmKSnyWGev4QkdOVX/E6n/xdJEqmrrVMzHEWZL1f9Eqo+p+t+KVP1PpPo7Y8YMVxWhIQthtlT9T6T66KsM9EEpzJaqKgMNodHcMPrgrv9xCrOl6n8iU/Xe1rlX9j+RxTO/r/i2stdsqSakbtjEGvXa2hpvdkiYXGk7loaraQjXB+9/YAWFBXbGGWe44Wb+vZ22jvHECCCAAAIIIIAAAgggkDaBjAiI0qbDEyOAQL0VyOuwkemnrrauXbvakMOG1NXT47wQQAABBBBAAAEEEECghgLhj7epYQfYHAEEEEAAAQQQQAABBBBAAAEEEEAgvQIEROn159kRQAABBBBAAAEEEEAAAQQQQACBtAsQEKX9EtABBBBAAAEEEEAAAQQQQAABBBBAIL0CBETp9efZEUAAAQQQQAABBBBAAAEEEEAAgbQLEBCl/RLQAQQQQAABBBBAAAEEEEAAAQQQQCC9AgRE6fXn2RFAAAEEEEAAAQQQQAABBBBAAIG0CxAQpf0S0AEEEEAAAQQQQAABBBBAAAEEEEAgvQIEROn159kRQAABBBBAAAEEEEAAAQQQQACBtAsQEKX9EtABBBBAAAEEEEAAAQQQQAABBBBAIL0CBETp9efZEUAAAQQQQAABBBBAAAEEEEAAgbQLEBCl/RLQAQQQQAABBBBAAAEEEEAAAQQQQCC9AgRE6fXn2RFAAAEEEEAAAQQQQAABBBBAAIG0CxAQpf0S0AEEEEAAAQQQQAABBBBAAAEEEEAgvQIEROn159kRQAABBBBAAAEEEEAAAQQQQACBtAsQEKX9EtABBBBAAAEEEEAAAQQQQAABBBBAIL0CBETp9efZEUAAAQQQQAABBBBAAAEEEEAAgbQLEBCl/RLQAQQQQAABBBBAAAEEEEAAAQQQQCC9AgRE6fXn2RFAAAEEEEAAAQQQQAABBBBAAIG0CxAQpf0S0AEEEEAAAQQQQAABBBBAAAEEEEAgvQIEROn159kRQAABBBBAAAEEEEAAAQQQQACBtAsQEKX9EtABBBBAAAEEEEAAAQQQQAABBBBAIL0CBETp9efZEUAAAQQQQAABBBBAAAEEEEAAgbQLEBCl/RLQAQQQQAABBBBAAAEEEEAAAQQQQCC9AgRE6fXn2RFAAAEEEEAAAQQQQAABBBBAAIG0CxAQpf0S0AEEEEAAAQQQQAABBBBAAAEEEEAgvQIEROn159kRQAABBBBAAAEEEEAAAQQQQACBtAsQEKX9EtABBBBAAAEEEEAAAQQQQAABBBBAIL0CBETp9efZEUAAAQQQQAABBBBAAAEEEEAAgbQLEBCl/RLQAQQQQAABBBBAAAEEEEAAAQQQQCC9AgRE6fXn2RFAAAEEEEAAAQQQQAABBBBAAIG0CxAQpf0S0AEEEEAAAQQQQAABBBBAAAEEEEAgvQIEROn159kRQAABBBBAAAEEEEAAAQQQQACBtAsQEKX9EtABBBBAAAEEEEAAAQQQQAABBBBAIL0CBETp9efZEUAAAQQQQAABBBBAAAEEEEAAgbQLEBCl/RLQAQQQQAABBBBAAAEEEEAAAQQQQCC9AgRE6fXn2RFAAAEEEEAAAQQQQAABBBBAAIG0CxAQpf0S0AEEEEAAAQQQQAABBBBAAAEEEEAgvQIEROn159kRQAABBBBAAAEEEEAAAQQQ+H/27gLOrSp94/g73nbq0BaKFSguRXdxFli8xd3dF3cW18XdpUjRRYu7LfDHobgUL6XUdXz+9znDCTeZZCaTuSG96e/wGSZzLSffm3QmT95zLgIIFFyAgKjgp4AOIIAAAggggAACCCCAAAIIIIAAAoUVICAqrD/3jgACCCCAAAIIIIAAAggggAACCBRcgICo4KeADiCAAAIIIIAAAggggAACCCCAAAKFFSAgKqw/944AAggggAACCCCAAAIIIIAAAggUXICAqOCngA4ggAACCCCAAAIIIIAAAggggAAChRUgICqsP/eOAAIIIIAAAggggAACCCCAAAIIFFyAgKjgp4AOIIAAAggggAACCCCAAAIIIIAAAoUVICAqrD/3jgACCCCAAAIIIIAAAggggAACCBRcgICo4KeADiCAAAIIIIAAAggggAACCCCAAAKFFSAgKqw/944AAggggAACCCCAAAIIIIAAAggUXICAqOCngA4ggAACCCCAAAIIIIAAAggggAAChRUgICqsP/eOAAIIIIAAAggggAACCCCAAAIIFFyAgKjgp4AOIIAAAggggAACCCCAAAIIIIAAAoUVICAqrD/3jgACCCCAAAIIIIAAAggggAACCBRcgICo4KeADiCAAAIIIIAAAggggAACCCCAAAKFFSAgKqw/944AAggggAACCCCAAAIIIIAAAggUXICAqOCngA4ggAACCCCAAAIIIIAAAggggAAChRUgICqsP/eOAAIIIIAAAggggAACCCCAAAIIFFyAgKjgp4AOIIAAAggggAACCCCAAAIIIIAAAoUVICAqrD/3jgACCCCAAAIIIIAAAggggAACCBRcgICo4KeADiCAAAIIIIAAAggggAACCCCAAAKFFSAgKqw/944AAggggAACCCCAAAIIIIAAAggUXICAqOCngA4ggAACCCCAAAIIIIAAAggggAAChRUgICqsP/eOAAIIIIAAAggggAACCCCAAAIIFFyAgKjgp4AOIIAAAggggAACCCCAAAIIIIAAAoUVICAqrD/3jgACCCCAAAIIIIAAAggggAACCBRcgICo4KeADiCAAAIIIIAAAggggAACCCCAAAKFFSAgKqw/944AAggggAACCCCAAAIIIIAAAggUXICAqOCngA4ggAACCCCAAAIIIIAAAggggAAChRUgICqsP/eOAAIIIIAAAggggAACCCCAAAIIFFyAgKjgp4AOIIAAAggggAACCCCAAAIIIIAAAoUVICAqrD/3jgACCCCAAAIIIIAAAggggAACCBRcgICo4KeADiCAAAIIIIAAAggggAACCCCAAAKFFSAgKqw/944AAggggAACCCCAAAIIIIAAAggUXICAqOCngA4ggAACCCCAAAIIIIAAAggggAAChRUgICqsP/eOAAIIIIAAAggggAACCCCAAAIIFFygvOA9oAMIIIAAAggggAACkQr89NNPdvrpp9u+++5rq666auLYY8eOdcv//e9/2wILLJBY/tZbb9mNN95ol156qfXu3TuxvL0bv//+uw0ZMsQeeeQRW3HFFVttrmM+9dRTdv/991tFRUWr9R1Z0NzcbCUlJR3Zpd1tzzvvPOvXr5/tv//+Sds+9NBD9uqrr9rll1/ulk+fPt0uu+wye//9961///625ZZb2mabbZa0zzvvvGM33XSTHXroobbaaqslrfM/fPXVV7brrrvaIYccYq+99ppfnPguoxtuuCHxc/jGhAkTTOdtm222sY022sit8uctvJ2/rW0XWWQR/6MdffTRtvjii7uvu+66yy1vaGhw38vLy23jjTe2HXfc0f71r3/Zeuut5+4nsXNwQ+dQj/Giiy4KL7YPPvjA2Zx77rlJz6n77rvPvvjiC+vSpYvbXufu+uuvT3oe1NfX2/rrr29XXXWVrbDCCjZp0iQbPny4jRo1yvTcGjRokG2//fa2zjrrJN0nPyCAAAII5EeACqL8uHJUBBBAAAEEEECgYAID5x1oTzzxhF1xxRVJfbjjjjtcEKQ34eGm7V588cUOhUPav66uzn788Uerra0NHy5xW4HKSy+9lHF9YsN2blx88cV24IEHtrNVx1ePHDnSnnvuuVY7vvfee3bPPfe45Y0NjS4wUbix4IILmoKazTff3AVI4R0Vril4UV8zNQUhK620kgtJFM6Ev1555RV75pln0u76/fff2xprrOHCo88//zyxjQKl8DF0W+tvu+22pDBt3LhxLsTZYostXPii9dq3srLSqqqq3Jduq915552mx5/a3n77bbv77ruTFj/55JPu8WgfhTvh9t///te+++67RP/69OmTFA5pW/nLU+HQ888/b0suuaSdeeaZNmXKFJtvvvlMAdi6667rAqTwsbmNAAIIIJAfASqI8uPKURFAAAEEEEAAgYIJlJWX2SabbGKPPvqoe4PuO/LAAw+4m6oGOfXUU91tBSCq8tl99939ZpF9V6iiUKSz1UPffPNNZH3q6IE++PADe/fdd11oonBHTYHNQgstlDiUKrNkPXToUHvwwQftt99+swEDBiTW68asWbNc+KLATNVWu+22W2L91KlTXfh01llnJZb5G7pvVSvpvhX0hNvKK69sCv3CbYMNNnDVPwsvvHBi8a233ur6NnDgwMSya6+51mbVzLLS0lLr1q1bYnm2N2655Rbbb7/9XIWPf175fRsbG10QpeeVnoeZmqqldIyJEye6yipVOD388MM299xzJ3ZRcLb00ksnfuYGAggggED+BKggyp8tR0YAAQQQQAABBCITqPvpbZsy8mgbd9Pm7ku36376v4zH32qrrVwlxuuvv+62GT16tAs6rrzySvvkk09cpYlWvPHmG267rbfe2jTkRyHFcsst5wKOnXfe2TS8Sk0hg97AK/jYbrvtTMPYfFMlyd/+9je3bp999nGVRVqngGjYsGFuMwUo2k9DlPSGX8OHNGzLN1WQaDichritvfbaduyxx9oZZ5zhhsQpNFBFlMKGu0e0VLHce++9bnv1R1UmGhLm24Ybbuj6q8BG6xWaqNIpl+aHtX344YeJ3dV3v1wLFZb06tXLrrvuOvddP6c2hXKLLbZY2qF4cqiurrYDDjggdTd3PIUoqrZRtU9bTQaqBPPhn7bV0DwNfTvooIPa2rXD6wb0H2Aj7hphZ599dqt9vbUekyqa9LxKbQrZnn32Wdtrr73cuVL4Nfy24UnhkPbRudUwQBoCCCCAQP4FYlFBNPrJC7OSCH7/WWN9g02tjPZhNTU1W31NS+n06G9bxlFn1aEsNlKfG4JP7iZVlGWxdcc20ac3JSWlwSdDJR3bsZ2tmwOPpsYmmxhxn9Xf+tp6aww+yRr9WUuZcztdyXq1zqH+QBpfFn0m2hA85yYFY/dLomW2pqamoM9mv0fc57q6lpL2yWVlNquqc/NBpJ6AxuB5IYbSiPucr9e2+u/6HHRan6BG2Wpq6tw5HBeU8Ef9+o7ba9s51zfahPLSpDdUUXjXzapxr5MxX1VG/ryL22tbng3B63tqRfA7MMJ/j/RvZ/2slt+B33/TJdJ/6+L42laf8/F7W/8W1dfWWUPwy2T0V61DgPKJX1n5+C+s99xV1nPulr9zZk39yMY/+Lb1WnMfq151n1YvKc0po6YqDoU3ClQ0J81hhx1mF1xwgWmOGAUwCh4UbiiUUbWKKlzOP+98q6mtcdUdCi/22GMPN5/RzTff7LZT2BMOSB577DE7//zz7YcffrCDDz7YzcGjsEPzyCiYUtOwIVXXaL6bq6++2lR1onlxNFxLgYvCI1UxqSJG/VOo8fLLL7tQRHMCqfrlsEMPs8GLDTYFUgqvFDapskZVSgoSNMePQhjdp46toXOaX0hz/pxwwgl2zTXXuL505H8rDFnB1lprLff41Sf1UR7du3d3h9FrROHZTjvt5Bw1x5ACopNOOinJSGFZ6lxHOsDkyZNdUKY5o/x8PeH+6fHoscitvaZgSMPINGTLNw2h0++NTTfd1C9y31948QUX3Og89uzZ0z1Gv4HmFdL8UeH20UcfhX+0ocOGup+//PLLpOX6wS/zcwfp+aWKpV123SWxraqHZNa3b1/7+OOP3XlbdPCiifXcQAABBBD46wXa/03z1/ep1T1OH/ttq2WFWjB9an7uueVP7/wcO25HbQo63PJZZdx6Hq/+yrn153nxegxx6K1e27y+83+mZub/LriHQGBGnn4HxhE3n6/r6VOSRSrqp1mXWT/awEFBINHtzw+UulWXW7fqChvzxm1WPvdiVrXw2kk7atiQKmhUeaPwQRUsO+ywgwstNNGxfj4jCIgef/xxFypoGNjqq6/uAiJ/oBEjRriJmRXeqP3yyy82zzzz2HHHHed+VuijpvDDT8ys26lhgtvoj/8pqFJTsOO3VdWQjqVKmaWWWsoNv1JQpHl31C/dp4Z0+VBC968ATFVGamuuuaYLVxTUKKhSO/HEE23PPfd0txV6qVoll6bheqrK0bw9mpdIk1BfeOGFbpJpTfStAEaBlL8vTQx+7bXXujl1VMmk5iu2dtnlz4DE90VzFqkySMFaZ9oLL7zgqqhS5w9SEKMqndQPQvRc8E3h22effeZ/dAGPQr9w+/rrr8M/tnlbj1suqkTr0aOHHXXUUbbrbrvaBv/cwFV0qaJI50rPQTWde4VI4SZzzeukpueWJt2mIYAAAgjkVyDaj8zz21eOjgACCCCAAAIIzHEC1Q0TrG//oII5FA4lEIJlfftV2awPWoZdJZb/cUPDxvTGXuGIQhtdpUpN4YCG/jw+8nH3XVUnvmnCYe2nN/e6OpnezCu00VWq9IZdb+RVQTNzZvpoVsOBMk1a7e9D3/3V0rStrgymK3MpIDrnnHNcoLDttttmnLtIw5M0obFvClg0PE2TIqdrmfqkKiAfcoX305xAXbt2TSxSSKWKKAVNqrJSpZSqqdT8VccUGqm6SFdGCy/3t2Xvq47cBsH/NPeOqqmOP/74tNVDfrtsvp922mkuEPTzJGkfzY2k6iu5prYJ4yfYzz//bGPGjEkKh7SdKnsUHIa/fEiYepx0P5cFVcoKz/z8Rn4YmiYtV1MFmp5HqvpSW2aZZdzzNDwUTVVlCpo09EyTVdMQQAABBPIvQECUf2PuAQEEEEAAAQQQyFmgvHGWdemeuehb62rGpq/u8MHPkUce6YaX+aFHGjKlqpEDDzrQ9ctPJKyQQ1U5GialOYH23nvvRL9PPvlk92ZdVR8KkXxAktigkzc0MbEmVtYQ61NOOaXN46vvqmbyTRNta3jZ/PPP7xdl9V3BxP/93/+5YV5+Bx1L1TjLLrusW6T+aJlvmhhaw750BTdNRq0ARk4KOzSkSpdt188aTqdww09One4qbP/5z39c9ZCGwHWmPf300/bGG2+4SrHwcTTUTVVkHXUJHyOK27/++qs7jCrB1BSqhT30PFVQd/nll7v1+p+GQx5xxBFMUJ0Q4QYCCCCQfwECovwbcw8IIIAAAggggEDnBNqaV0rrNDlSmqbQRXML6Q26rx7SZpp3RsPMVD2iOYB8ZYvCElXiqFpkxowZicoNzSWkIUHaTxM+q+LH75PmbnNapCoRVTWpT5oDRwGLryhRtYkqoVQhpAmQ9Vh0GXUNn1OVzEknn+QChmyqXDS/kb/yl4Y+KZhQ9ZKqqTQsaqutt3L90BA1NYVVyy2/nAuCNAxLw/LUFwUvfjJqVVYp8PDzHV1yySVuX5lpWJoCpXBlj1bK1Ady4WolrdOcQxqm5puqtaZNm+Z+VMWVbiu48k3VQzqPq6yyil/k5l7UnEkHHtASAiZWRHRD969+6Hmipu9+QnMFQBo+qHOlIW9yGTJkiC2/3PJuOJ7CRz8kT/tqCNnhhx/uKqn02PVcUJWWzm/qkDltT0MAAQQQyI8AAVF+XDkqAggggAACCCAQiUBTWRerndGQ8Vi1Mxqtau6FMq73VxELzzmjjbfffnu3j4aT+aahSLrqmYZ/KQDxgYMqZXT1s8rKSptrrrnchMbhy7RrSJFvGlbk57vx37UuPKm131bfS4MLaigU0pxJqr7RhMm6QtoxxxzjJnrWNqrIUWikqhJVMOlnzaukSqd5553Xbr/9dhdg+XmQ3HFDFx9Q/3wfNWm3KoTUNK/R22+/7cKUYVsMc0GRJlge+djIxKTNqr7SVdd0dTZVHOlKZZrAWfMeKSBS0OGHy7mDBv/r06ePG4angEZBSbrJqTX3kO5fcxqlNk0gruFdvi2//PJuSJYqkjTZtoyeeeYZt1qP55133rEzzzzTb+6+a70Cts023yxpuf+hpDRz6hg+b4ntg3Aw3EaNGuX6oYoqNc0XpfmGVDGl54gmIdfzR1+aiFpzGmk+J4VHGsIWvpS99teE4vJUNZTmwho0aJB7DmqomkI5GgIIIIBA/gVKgisvpP/I6Y/71upMv9Bz7Z4+XRh9/zG57s5+CCCAAAIIIIDAHCNQVT/Fetb+YvMuEkxSXZny2V5dk439ocb6/PMEq1pik7QmGgqlCpTUECPtxsFCba85eFLfwGt7VYxoPh5/tS1dgUuXMteyzjRdKU0hj7+UfE1NjQuwFDQMHz7cHVphh67G5e/b35/+rsy2mklz/qiqRRYKIMJN96n7UMiRrulvYrmkTqbst1Uw0tDQkLS/vBTcaA6e1Cohv1+67zpfMvWu48ePd6FLR/4mv3vE3danb59WVy8L35/sFAYpnIuyKfDTY1fTsEEFi77JQ1egW3zxxf2iVt9lqb717dPXhUrhDXRcXdGtI57h/TPd1jDCqdOmunAv0za5LleVmp63misryqaKOr1eFZIqNIyyRfXajrJPHAsBBPIvkHlAe/7vm3tAAAEEEEAAAQQQaEegtqKXzWqcab9+N9nm6l9llX/MR1Q3vcEmT2i06uWGZgyH2jl02tV6M58uHNLGmcKTtAfqwELNgaTgQMPbllhiCXcVLVUt+UofHSocmIQPnW045PdRsJDukvF6A58aPvl99F3hTKZwKLxd+La8VCnU0RZFYBO+pHxH7z+q7VU5FA6HdFxVfrXXdI6iDoDau0/WI4AAAgiYxSIgWn6fG7I6VxoLrU+GMv1Rk9VB0mykT5Q0DlpNf7RE2eL4aUVHP4nM1kufrmgeAX26kvqpXrbHyLSdzqH6HfWnK7o/Ped69ezV6hOuTH3Jdnm6TyKz3bet7VSePmnSJGehT5yibPn6JDJfr2099nx9EqlPqHUONU+HhhpE2fL1SWS+Xtt67Pn6JFLzgOj5oSEIUbyhCp+nuL221fdcqgzCjzndbVVU6BLeaoMXHRzpv3VxfG3n6/e2/v3UpMvlZeW26OBF050Kq/l8pE3/4F6r/eZHU/l3l34LWs9/7GRdlh6Wdvs4LdRV0nRFLc0BpLl5dtpxJ1tt9dUi/5suTib0FQEEEEAAgb9aIBYB0V+Nwv0hgAACCCCAAAKzm0CXpYaZvoq1qRLIz4fkP9gp1sfK40IAAQQQQGB2FEgZyD47dpE+IYAAAggggAACCCCAAAIIIIAAAgjkU4CAKJ+6HBsBBBBAAAEEEEAAAQQQQAABBBCIgQABUQxOEl1EAAEEEEAAAQQQQAABBBBAAAEE8ilAQJRPXY6NAAIIIIAAAggggAACCCCAAAIIxECAgCgGJ4kuIoAAAggggAACCCCAAAIIIIAAAvkUICDKpy7HRgABBBBAAAEEEEAAAQQQQAABBGIgQEAUg5NEFxFAAAEEEEAAAQQQQAABBBBAAIF8ChAQ5VOXYyOAAAIIIIAAAggggAACCCCAAAIxECAgisFJoosIIIAAAggggAACCCCAAAIIIIBAPgUIiPKpy7ERQAABBBBAAAEEEEAAAQQQQACBGAgQEMXgJNFFBBBAAAEEEEAAAQQQQAABBBBAIJ8CBET51OXYCCCAAAIIIIAAAggggAACCCCAQAwECIhicJLoIgIIIIAAAggggAACCCCAAAIIIJBPAQKifOpybAQQQAABBBBAAAEEEEAAAQQQQCAGAgREMThJdBEBBBBAAAEEEEAAAQQQQAABBBDIpwABUT51OTYCCCCAAAIIIIAAAggggAACCCAQAwECohicJLqIAAIIIIAAAggggAACCCCAAAII5FOAgCifuhwbAQQQQAABBBBAAAEEEEAAAQQQiIEAAVEMThJdRAABBBBAAAEEEEAAAQQQQAABBPIpQECUT12OjQACCCCAAAIIIIAAAggggAACCMRAgIAoBieJLiKAAAIIIIAAAggggAACCCCAAAL5FCAgyqcux0YAAQQQQAABBBBAAAEEEEAAAQRiIEBAFIOTRBcRQAABBBBAAAEEEEAAAQQQQACBfAoQEOVTl2MjgAACCCCAAAIIIIAAAggggAACMRAgIIrBSaKLCCCAAAIIIIAAAggggAACCCCAQD4FCIjyqcuxEUAAAQQQQAABBBBAAAEEEEAAgRgIEBDF4CTRRQQQQAABBBBAAAEEEEAAAQQQQCCfAgRE+dTl2AgggAACCCCAAAIIIIAAAggggEAMBAiIYnCS6CICCCCAAAIIIIAAAggggAACCCCQTwEConzqcmwEEEAAAQQQQAABBBBAAAEEEEAgBgIERDE4SXQRAQQQQAABBBBAAAEEEEAAAQQQyKcAAVE+dTk2AggggAACCCCAAAIIIIAAAgggEAMBAqIYnCS6iAACCCCAAAIIIIAAAggggAACCORTgIAon7ocGwEEEEAAAQQQQAABBBBAAAEEEIiBAAFRDE4SXUQAAQQQQAABBBBAAAEEEEAAAQTyKUBAlE9djo0AAggggAACCCCAAAIIIIAAAgjEQICAKAYniS4igAACCCCAAAIIIIAAAggggAAC+RQgIMqnLsdGAAEEEEAAAQQQQAABBBBAAAEEYiBAQBSDk0QXEUAAAQQQQAABBBBAAAEEEEAAgXwKEBDlU5djI4AAAggggAACCCCAAAIIIIAAAjEQICCKwUmiiwgggAACCCCAAAIIIIAAAggggEA+BQiI8qnLsRFAAAEEEEAAAQQQQAABBBBAAIEYCBAQxeAk0UUEEEAAAQQQQAABBBBAAAEEEEAgnwIERPnU5dgIIIAAAggggAACCCCAAAIIIIBADAQIiGJwkugiAggggAACCCCAAAIIIIAAAgggkE8BAqJ86nJsBBBAAAEEEEAAAQQQQAABBBBAIAYCBEQxOEl0EQEEEEAAAQQQQAABBBBAAAEEEMinAAFRPnU5NgIIIIAAAggggAACCCCAAAIIIBADAQKiGJwkuogAAggggAACCCCAAAIIIIAAAgjkU4CAKJ+6HBsBBBBAAAEEEEAAAQQQQAABBBCIgQABUQxOEl1EAAEEEEAAAQQQQAABBBBAAAEE8ilAQJRPXY6NAAIIIIAAAggggAACCCCAAAIIxECAgCgGJ4kuIoAAAggggAACCCCAAAKayk3QAABAAElEQVQIIIAAAvkUICDKpy7HRgABBBBAAAEEEEAAAQQQQAABBGIgQEAUg5NEFxFAAAEEEEAAAQQQQAABBBBAAIF8ChAQ5VOXYyOAAAIIIIAAAggggAACCCCAAAIxECAgisFJoosIIIAAAggggAACCCCAAAIIIIBAPgUIiPKpy7ERQAABBBBAAAEEEEAAAQQQQACBGAgQEMXgJNFFBBBAAAEEEEAAAQQQQAABBBBAIJ8CBET51OXYCCCAAAIIIIAAAggggAACCCCAQAwECIhicJLoIgIIIIAAAggggAACCCCAAAIIIJBPgfL2Dv7b+49YqZW0t1mH1tc11HVoezZGAAEEEEAAAQQQQAABBBBAAAEEEMifQLsB0biPns7fvXNkBBBAAAEEEEAAAQQQQAABBBBAAIGCCzDErOCngA4ggAACCCCAAAIIIIAAAggggAAChRUgICqsP/eOAAIIIIAAAggggAACCCCAAAIIFFyAgKjgp4AOIIAAAggggAACCCCAAAIIIIAAAoUVICAqrD/3jgACCCCAAAIIIIAAAggggAACCBRcgICo4KeADiCAAAIIIIAAAggggAACCCCAAAKFFSAgKqw/944AAggggAACCCCAAAIIIIAAAggUXICAqOCngA4ggAACCCCAAAIIIIAAAggggAAChRUgICqsP/eOAAIIIIAAAggggAACCCCAAAIIFFyAgKjgp4AOIIAAAggggAACCCCAAAIIIIAAAoUVKG/v7vsP2cRKraS9zTq0vq6hziZ++kKH9mFjBBBAAAEEEEAAAQQQQAABBBBAAIH8CLQbEA1YaSsrKYk2IJo+fToBUX7OJ0dFAAEEEEAAAQQQQAABBBBAAAEEOizAELMOk7EDAggggAACCCCAAAIIIIAAAgggUFwCBETFdT55NAgggAACCCCAgP3000+2zz772P/+978kjbFjx9qBBx7o1odXvPXWW277yZMnhxe3e/v333+3wYMH2wcffJB22xtvvNG23nprq6+vT7u+Iwubm5s7snlW25533nl20003tdr2oYcesiOPPDKxXNXvZ599tnss8nvyyScT6/yNd955xw455BCTZab21Vdf2aqrrprkcc8999jBBx/cahdte/jhh7v7PPPMM+2XX35ptU14gc7BHnvskXRuTznlFHde9VwIf1166aXhXV1/1l57bfvwww/t2muvTdrW76d1P//8s1v3+eefJ+2vH84666y0lpkenz/Ajz/+6I75ySef+EVJ35999lnbd999k8z8BkcffbRdf/31/kd7/vnn7cQTT3RmO+ywg+lxTpkyJbGeGwgggAACbQsQELXtw1oEEEAAAQQQQCB2AgPnHWhPPPGEXXHFFUl9v+OOO0yhzfDhw5OWa7sXX3zRevfunbS8vR/q6upMb/Bra2vTbvr+++/bSy+9lHF92p3SLLz44otdsJVmVacWjRw50p577rlWx3jvvfdMwYZaY0Ojrbfeei6IWHDBBW3ChAm2+eab26uvvpq0n8KIu+66y9TXTO2qq66ylVZaySoqKtwmCn522WWXxH35/RQO/e1vf7NRo0bZ8ssvbwpJ1lprLZs2bZrfJOm7Aisd984777RJkyYl1mmaiIaGhqSv2267zT766KPENrohBz2uFVZYwXnouVNVVZX0VVlZaRMnTjTtny6s0jH0HAo3WaR7fH4bBYfbbbedO6aeR6nthx9+MAU9t956q3sM4fXjxo2zyy67zLbYYgsXHu2000624YYbuj7069fPbarQaumll7aamprwrtxGAAEEEMgg0O4cRBn2YzECCCCAAAIIIIDAXyxQ8+kjVjf2M7PmJquYZynrusw2FkwW2aoXZeVltskmm9ijjz6a9Mb6gQcecNvef//9duqpp7rbCkCeeuop23333Vsdp7MLVN2hUMQHIrke75tvvsl1107v98GHH9i7775rCo0Uwqh9//33ttBCCyWOrcosWQ8dOtQefPBB++2332zAgAGJ9boxa9YsF+AoMFPba6+97Pbbb3dh0+uvv+6W+f9dc801NmjQIHvhhRestLTUttlmGxfevPLyKzZ02FC/mft+yy232H777Wfbb7+9+fPrNzjnnHP8Tfdd96MQSZU34XbDDTe4Y/hlqgq77rrr/I+J7x9//HHidns3Tj75ZFMllsK01Mfn9z3hhBNM1Vm9evXyixLfFR7pMS2zzDL2xhtvJJb7GwqN5D1w4EBX3XXffffZ3XffbTvvvLPfxFQRp3PXpUuXxDJuIIAAAghkFohFQDTuvUcyP4LQmiZrdr98m7p2Cy3t/M2GxkZrnNpSnjpueuuS2s7cQ2NQLl0ffOpWn4dfXLV1NTajvMLKS8s608VW+zY2NVp9Q73VVUb7y1b9bZwx0+qC/o6b8GGr++3MAncOGxusprKqM4dJu+/M4A++2uD8laX5Az3tDlkurA8+8WsK3gDMqqjMco/sNps+c6Y11tbYzOCTwHFjume3U5Zb6QqFatPLo+1zvl7b6mttfZ2VlpRaRXm0/xzWTpsafOrcYNPGd4v89R2317aca4LX96zgeVEWvNGJsjUEn5Q3B/9NnNw98udd3F7bcp0xa6Y1dO0a6dVHG5uarHHKZHfaxk37LNJ/6+L42s7X7239+9k4bbo1B/8ejZsyqvXLZOqP1vTlE1Ze0mBdq1tWz/zuFZv6xs3Wd6NTrWLQGq322WqrrUwVQ3pzrmqU0aNHuzfLV155pRu6pGFCSy21lL3x5htuGI4fCnb++ee7oEEVGuuvv74bNtS9e3d78803XbCg4yy++OKm7fr06ePuVxUsGg6lig8FAgqGVHGi74888og9/fTTLkBROPH3v//dVYzMDH4fHXHEEXbUUUe5Y2h40EknnWRff/21LbbYYrbiiiu6N/8aWvbwww9befDvtEKvPXbfw3bZdRe799577ZJLLnEVTEsuuaQLCdZZZx13LFWSKChQQKGhX8suu6y7Tz2OjjZ/sRYNsfIBkcKbcFNIo5BDocprr71m+lkBSbgplPOPS8uXW245dz4UYKQGKMccc4wpPFE4pNa9uqXf+vsrtQ3oP8BG3DXCVl5l5VYBUeq2CgX1vBgyZEhilc6ZKpR8xVRiRSdv6Dmi6qxvv/221ePToXVOVQGkKqkNNtig1b0dd9xxLuBRv1ZZZZWk9XpOaGignsu6reeBhteFwyHtoIq4f/7zn0n78gMCCCCAQGaBaN8RZb6fTq0Z+9FTHdo/ffFthw6RceOxratfM27LitwEGoLdxua2K3t1QEDOMzqw/eywaT5f2/l6fC1vq/N1dI7rBSb6G3y3fL5Ofs/T78B89jmOT4mxPyT3uqy53vrM+Mb69qu0LnP9+eGMcqLaiXU2fuTJ1n/X4VbWd1DSjhtvvLH7WdVBCogUqCyyyCJ22GGH2QUXXGCquDjjjDPc8CKFG5qDRmGFKlzOP+98qwk+TFBlit7EK1jYbbfd3JttVWloeFVj8AGab4899pgLjBQ2aD6d1VZbzQ444ADTHEUKlNQ0F4yqaxT0XH311S7MUCWLAiUFLhpqpComhVrqnwKAl19+2Q1zUtCjSpHDDj3MBi822M0BpDDgoosucqGNqpTWXXdd09AshTC6Tx1bQ+f2339/NzeQAhdV5nS0rTBkBTe8S/PgqE/qo8IIHzYpoFA1i4Y4yXHXXXd1AZHCLh8u6T4VlqkvvikEUpN5atNQtnC78aYb3bEV2KU2X1H05Zdfpq5K+lnDv+SpYX/hpjBLfe/bt29i8dSpU131kl+gaiiFbL4pEPTn1S8bM2aMm4/K/6wKKR1TAVFqUwXW3nvv7UzCx/Xb6Xzr3CmU01C31PbiCy+659+mm27qQkk9tzbeqOX5nrotPyOAAAIIZC/Q8rFE9tuzJQIIIIAAAggggMBfKNC1drx1614ehEOt3yhX9a207j3LbPo7N7fqUbdu3dwQHM0no6YKFs3notBCQ5b0s9rjjz/u5nHRMLDVV1/dBUQKHRTYaO4dHyhoKND48eOtrKzM9txzTxfIuAME/1P4oTfrBx10kKtOSZ3jxm+n7wqqFHSoAklN22ookN7kK5BSVZPCKLU11ljDhU3zzDOPG9KlfqlaSJU6CsCOPfZYdywdU01BjW+arFj93HLLLd1Qpbb65PdJ913D9RSuaBiWhiodeuihLizRROBqqpJRWKL7UlOQpJ9VEeWbJmBWxZbm4+loe3zk43bhhRe6uaPSDcXK9ninn366ez6oMss3zU+kQE7nLdzU32233TbxpUqdcFM1mULB8NeMGdl97KR5qzR0TFVemvA7tX333Xdu0mqFeeFKp/B2N99ysxuipworP7F6z549w5u455KGoOlLgRQNAQQQQKB9gVhUELX/MNgCAQQQQAABBBAoToEujVOtR+/W4ZB/tD16Vdhv371jrWdxMXc1JwUWGkKkgMRPTq2gSFU8Ch8UXGgyX980DEwVHJr3R0GBKnzUFL5oGJnm3tFQrnRX/9J2miA406TVWu+bnxBb2/bv39822mgj96ZegY4CHwUUmeYu0hv+8LAkVZno6mAKF9K1TH1SFVC6q1ypgqZrMFzTN/VDFVH60lxEGvJ08803myaZVnCkphBH1USqkFLTcoUg/vaOO+6YqDpyC7P4n4bmDdtimKuoUqiSa1OIpWFsqdVKCgd1HlR9FW4K5lKvgBderyu6pQ7dkn827fLLL3f9UBCo55OahjPq6mm6Spqeezonel5qvQ+AVBEmg4UWXMgNV9TzV03HUfv0s0+T5mdaf7317cuvvnTPbZ0nVanREEAAAQTaFqCCqG0f1iKAAAIIIIAAAoUVaG6ysoqSzH0I1jXVzUq7Xld4UtMl2zW8TFepUtMVsTRk68CDWio4NLePmkIOVeVomJQqYDQMyDcFRRoupOoRDR3TEKoo29xzz20LL7ywNQXzXuny7ApgMjX1PXwlLU20reFl888/f6Zd0i7XBMj/93//lwghtJGOpcmh/dAn9UfLfFt55ZXdMDZVwmgyagUaclLIouBM1VH6WcPpFHz4yanTVcv4Y6b7/swzz7iqLM2vo5CvM+20005zQZ/6Hm4K/XQZ+7+yyUfzMy2wwALWo0cP96XKKIV1mmtJ/dF6VQRpvSrh1PSzKrjuuPMO5+LPtcJBVa+pyknBnm+ap6qzbv5YfEcAAQTmFAEqiOaUM83jRAABBBBAAIFYCjSXVlh9TXDVsqqytP1vqGm0iu5/zh8T3kihi+YW0sTJ4UDHDzNTFYaCHz+fjsISVYJoeJmGSb311lsuWFJljCpitFwhg+bIqa7+Y7bs8B124rbuS2GA5q3RcDZNYq25flS9owBBk1erQkhD3FSNozf/Gj6n/ugS86o6Uf/aa5rfSPMUKYTQBNma60bVS6edeppVdalyQYOqV2688UZ3KIVVukKZrgim4W0amqe+aOiS5u9R0/0rqNCQLYUautS8LgevAEbVUbo/P8G12yH4n67SpcuvK0BS0yXs9VgVgqjiS6Gdzo0ms37llVfccVUto3OlsElzSSmQUoCl4V1+iJe+y8+fU80XJFtN1h1uOr+ab0pXUstH0xxV6kvq41P1Veqk0woDda4322wz1xVVuPn2xRdfuOee5suSjeaekne4qfpI4aeOe95559mKK6xos2pmuUq08HbcRgABBBBoW4AKorZ9WIsAAggggAACCBRUYGZZb5s8vi5jH6YGs+F3WWrTjOuHDRvm1oXfdGuBH7Kkq5f5pjmANBRJw44UgPg38gocFMZoqJauTqbLumvYjm8KbXzTvDD+6lv+u9aFJ2z22+q7riSpCY4VDqm6RJUimt9Hkzj7UEsVOR988IELqxTQ6GfNqaNKp3nnndeFHApjNDm2b+H7Vv98HzVptyqE1DRc7u2333ZDwzSUS0GRJnse+dhIV2WlbVR9tfTSS7vKFlUcaf4jhUdrrrmmC4g095AfLqft1XR1NwUeGoaXOjl1yxbmJszWY1VIpXBLt/1V2DTRt5rM1ScFX5pzSRVeqlrS5ewV5qnpKmDa11cHaXiYQiofzChY0f7+XLqdgv/5wG+uuebyixLfw3aJhaEb7a3XpiNGjHAVYekeX+hQHbqpqiqFcH4Cdr+zho9pOOTyyy/vhilqInMFa7oCmoZP+mowvz3fEUAAAQTSC5QEnwg1p1/VslSrM/1Cb2u/ttbpj4zR97dcuaGt7ViHAAIIIIAAAgggYNZ31nfWrbLe+s4bXMWs8o/P9+qbbPK4Bmuq6m99dxmRkUlDoVSNkxpiZNpB22uojqqPUpvWKXjwkyVrfhhVEmWaKyh1/0w/60ppqmTRVavUVFmjAEvhxfDhw90yVdyoKkVVJOEWrpYJL093e+LEifbjjz86i9Q5aXSfug+FK+ma/iaWi3/sqdvIxVcQ+XWqClKQpjl7wnMa+fXZftfk4LLwf5OrHwqFOtPUL833pMvR+2FcnTleeF8FfnrsqgZTqBhVu3vE3dalaxc3pKwtTw2B1HpfRZXN/WsY4dRpU124l832HdlGAaCet+muyNaR46Ruq4o6vSYVknb2+ZB67Khe26nH5WcEEJi9BWIxxGzRTY7KSrE5GKM/dfo069Uj3TSNWR0i7UZ19bXBpzXj3LoF5l8g7Ta5LmxqarAZwR9tPbp37pd8uvufHpT1VlVVWEV5ZbrVOS+rb6y32uCPqO7V6f+AyvXA6u+kSRODPzLLbZ4B8+Z6mLT76Rw21DcEfwBFWw6vO9MfE92ruwWflkb7cqqtq7Wm4A/hrl1bxt6nfWA5LJw8ZVLwR9t0q3Zl/K0/NczhkIldaoJy7uCjYOtSmfzHe2KDHG/k67Wt7syaNdNKg0+WqyozTwCbS7fH/T7OTdLau3evyF/fcXtty2/6jGnuj+PysopcODPu8/MvvwSf/DdZv/79In/exe21LaQp04IqhB49rST4L6rWEPyb/+uvY93h5hs4b6T/1sXxtZ2v39v69/N3XSEsqL4ZOHC+NKev2WqCK5WN+eb/rCJ4g6wzXDerxqqXXMf6bHBqmu1zX6QKoXThkI6odfqKumk4lQILDW9bYoklXBWPKmV8pY/uTyFUuiCqIyGAjqPgwE8krZ990xv41PDJr9N3hTOZwqHwduHbCpsuvvji8KJIbkcRBqgCS+FanJrmFVLw1F6LMpRq775YjwACCBSTQLTvaPMkUz1wyayOrDHYNcEnQ9VpPvHK6gAZNioLwpCSmpY3vNUDl8iwVW6L9WlFQxAwVAelyFG3hjx9WqFPKoKPIq06KD+Psqm/JU1jrSwYw189cFCUhzadQ/W7upOftqXrVG3wnOvWs5fpMrhRttI/PomszvBJZq73NaN8nJWUTrLywKI6+MQpytYc/KGpsvOoP4nM12tbj70p+ENTbxTa+iQyF6Oyhm5WEpzDqn79I399x+21Lb/6oMqgWwRVBqnnonRGmZt/o2v/BSJ/3sXttS2bWUHAUB2qMkj1yuVnVVSUzGwJBLrNMzjSf+vi+NrO1+9t/ftZUh982FBWHvwOXDTtqare8hJrnjXR6n77zIJk1Cr6L2ml1dFVZ6S9079ooebW0RWsdLUtVX/stONOttrqq2UMqv6ibnE3CCCAAAIIzFECsQiI5qgzwoNFAAEEEEAAAQQyCJR07WtVg9bKsDbei1UJ5OdD0nAv94FUvB8SvUcAAQQQQCBWAkxSHavTRWcRQAABBBBAAAEEEEAAAQQQQACB6AUIiKI35YgIIIAAAggggAACCCCAAAIIIIBArAQIiGJ1uugsAggggAACCCCAAAIIIIAAAgggEL0AAVH0phwRAQQQQAABBBBAAAEEEEAAAQQQiJUAAVGsThedRQABBBBAAAEEEEAAAQQQQAABBKIXICCK3pQjIoAAAggggAACCCCAAAIIIIAAArESICCK1emiswgggAACCCCAAAIIIIAAAggggED0AgRE0ZtyRAQQQAABBBBAAAEEEEAAAQQQQCBWAgREsTpddBYBBBBAAAEEEEAAAQQQQAABBBCIXoCAKHpTjogAAggggAACCCCAAAIIIIAAAgjESoCAKFani84igAACCCCAAAIIIIAAAggggAAC0QsQEEVvyhERQAABBBBAAAEEEEAAAQQQQACBWAkQEMXqdNFZBBBAAAEEEEAAAQQQQAABBBBAIHoBAqLoTTkiAggggAACCCCAAAIIIIAAAgggECsBAqJYnS46iwACCCCAAAIIIIAAAggggAACCEQvQEAUvSlHRAABBBBAAAEEEEAAAQQQQAABBGIlQEAUq9NFZxFAAAEEEEAAAQQQQAABBBBAAIHoBQiIojfliAgggAACCCCAAAIIIIAAAggggECsBAiIYnW66CwCCCCAAAIIIIAAAggggAACCCAQvQABUfSmHBEBBBBAAAEEEEAAAQQQQAABBBCIlQABUaxOF51FAAEEEEAAAQQQQAABBBBAAAEEohcgIIrelCMigAACCCCAAAIIIIAAAggggAACsRIgIIrV6aKzCCCAAAIIIIAAAggggAACCCCAQPQCBETRm3JEBBBAAAEEEEAAAQQQQAABBBBAIFYCBESxOl10FgEEEEAAAQQQQAABBBBAAAEEEIhegIAoelOOiAACCCCAAAIIIIAAAggggAACCMRKgIAoVqeLziKAAAIIIIAAAggggAACCCCAAALRCxAQRW/KERFAAAEEEEAAAQQQQAABBBBAAIFYCRAQxep00VkEEEAAAQQQQAABBBBAAAEEEEAgegECouhNOSICCCCAAAIIIIAAAggggAACCCAQKwEColidLjqLAAIIIIAAAggggAACCCCAAAIIRC9AQBS9KUdEAAEEEEAAAQQQQAABBBBAAAEEYiVAQBSr00VnEUAAAQQQQAABBBBAAAEEEEAAgegFCIiiN+WICCCAAAIIIIAAAggggAACCCCAQKwECIhidbroLAIIIIAAAggggAACCCCAAAIIIBC9AAFR9KYcEQEEEEAAAQQQQAABBBBAAAEEEIiVAAFRrE4XnUUAAQQQQAABBBBAAAEEEEAAAQSiFyAgit6UIyKAAAIIIIAAAggggAACCCCAAAKxEiAgitXporMIIIAAAggggAACCCCAAAIIIIBA9AIERNGbckQEEEAAAQQQQAABBBBAAAEEEEAgVgIERLE6XXQWAQQQQAABBBBAAAEEEEAAAQQQiF6AgCh6U46IAAIIIIAAAggggAACCCCAAAIIxEqAgChWp4vOIoAAAggggAACCCCAAAIIIIAAAtELEBBFb8oREUAAAQQQQAABBBBAAAEEEEAAgVgJEBDF6nTRWQQQQAABBBBAAAEEEEAAAQQQQCB6AQKi6E05IgIIIIAAAggggAACCCCAAAIIIBArAQKiWJ0uOosAAggggAACCCCAAAIIIIAAAghEL0BAFL0pR0QAAQQQQAABBBBAAAEEEEAAAQRiJUBAFKvTRWcRQAABBBBAAAEEEEAAAQQQQACB6AUIiKI35YgIIIAAAggggAACCCCAAAIIIIBArAQIiGJ1uugsAggggAACCCCAAAIIIIAAAgggEL0AAVH0phwRAQQQQAABBBBAAAEEEEAAAQQQiJUAAVGsThedRQABBBBAAAEEEEAAAQQQQAABBKIXICCK3pQjIoAAAggggAACCCCAAAIIIIAAArESICCK1emiswgggAACCCCAAAIIIIAAAggggED0AgRE0ZtyRAQQQAABBBBAAAEEEEAAAQQQQCBWAgREsTpddBYBBBBAAAEEEEAAAQQQQAABBBCIXoCAKHpTjogAAggggAACCCCAAAIIIIAAAgjESoCAKFani84igAACCCCAAAIIIIAAAggggAAC0QsQEEVvyhERQAABBBBAAAEEEEAAAQQQQACBWAkQEMXqdNFZBBBAAAEEEEAAAQQQQAABBBBAIHoBAqLoTTkiAggggAACCCCAAAIIIIAAAgggECsBAqJYnS46iwACCCCAAAIIIIAAAggggAACCEQvQEAUvSlHRAABBBBAAAEEEEAAAQQQQAABBGIlQEAUq9NFZxFAAAEEEEAAAQQQQAABBBBAAIHoBQiIojfliAgggAACCCCAAAIIIIAAAggggECsBAiIYnW66CwCCCCAAAIIIIAAAggggAACCCAQvQABUfSmHBEBBBBAAAEEEEAAAQQQQAABBBCIlQABUaxOF51FAAEEEEAAAQQQQAABBBBAAAEEohcgIIrelCMigAACCCCAAAIIIIAAAggggAACsRIgIIrV6aKzCCCAAAIIIIAAAggggAACCCCAQPQCBETRm3JEBBBAAAEEEEAAAQQQQAABBBBAIFYCBESxOl10FgEEEEAAAQQQQAABBBBAAAEEEIhegIAoelOOiAACCCCAAAIIIIAAAggggAACCMRKgIAoVqeLziKAwOwo8Pzzz9sBBxxgv//+e1L33nv3PTv55JOtubk5afntt99uJ510UtKybH546KGHbPHFF8+46YEHHminnnpqxvXZrkjtb7b7tbVdQ0OD7b///vbaa6+12uyCCy6w66+/PrF89OjRduSRR9qWW25pxxxzjH322WeJdf6GDPfZZx+bPn26X9Tq+7333uu20Ypff/3VTjjhBNt6663t0EMPtc8//7zV9uEFH3zwgR188MH2088/JRbfcccd7ni63/CX3MOtvr7e1l57bfvwww/t2muvTdpW+x100EE2atQo+/nnn926dH0566yz7Kabbgof1t2eMmWKHX744XbbbbclrdNzcK+99nLHu/DCC23WrFlJ6/XDV199Zauuuqqpf2q633POOcd22GEH53LKKafYt99869bxPwQQQAABBBBAAIE5T4CAaM475zxiBBCIWGDuued2b+YfeOCBpCPfeNONLvj43//+l7RcodEvv/yStCybHxSGfP311xk3Vfjy9ttvZ1yf7Yphw4bZ3SPuznbzrLZT6KRQR8FIahs5cqS98MILbvGPP/5oK620ksls8ODB9tZbb9kaa6xhM2fOTOymgOP44493IckjjzySWJ5649xzz7WNNtrIxo8fb6ussoo999xztvzyy9v7779vq6++ugtoUvfRz08++aTrw/3332+TJk1KbFJaWmoKusJfr7zyij3zzDOJbXRDj2fChAm2wgoruPt84oknrKqqKvFVWVlp+po4caJ7DOmeCzrGiy++mHRc/aCA6b777rONN944sU5B2IYbbmhNTU02aNAgO//8822XXXZJrPc3rrrqKve4Kioq7OKLL7all17a3X+XLl2sT58+NmLECBu82GB77733/C58RwABBBBAAAEEEJiDBMrnoMfKQ0UAAQQ6JNAw4RtrmPSDlVZWW8X8q1hJafp/MhUELLTQQvb444/bIYcc4u5j6tSppvBATW/o11prLXdblSljxoyxrbbayv0c5f9GfTzKSss6n/urgqdQTWGHwpTXX3/dfVc/vvvuO+vWrVuiS48++qiNGzfOhg4d6mzThSEKlrTNNtts46p4FOromDrOcccdZz169LBnn33WBS6JAwc3brnlFttvv/1s++23t9TAb7fddjN9+aZzvOCCC5qqfcLthhtucMfwyxR0XXfddf5HV92j/uhxdaQp4FEVmZ5XAwcOTOyqCqzDDjvMtF5hmAK2vffe29S/nj17uu1UUXTnnXfaSy+95Kq4ZKB9Lr30UlNgpNbY0GhPPfWUrTBkBfcz/0MAAQQQQAABBBCYswTSv9uZzQxmjPkiqx41NzdZ3fRpNqNufFbbZ7tRXX2tNU8Z5zafMSZ5qEi2x8i0XVNTg9UGn4zPmNXyR3ym7XJZXjNjhjVXVVhDeWUuu2fcp76x3mpramzGzB4Zt8llhevvlInWWFFuMyprcjlExn10DhvqG2zG9OqM2+S6onba1KC6oZuVZggPcj5uXa01NTZa6bQ/3xjneqzwfvVTJlnztOmBRzeb0TwlvKrTt2trgmEtpSXWXNml08cKHyBfr23dR+2smVZfVmZNlVWJu2yaNNpq3r8jCId+sYouVdbc1BxUjTRZz1V2turVDkhsF76h4VBXXnmle/PftWtXe/rpp02VRar4uOuuu+yKK64IniOlrrpE+22yySbuzbwqYVQpoqoYDVM7/fTT3Rv1o4852v773/+aqjt0bFWF+KZj3XPPPe6+jjjiCDvqqKPcqv0P2N/mn39+F1houFB5ebmNHTvWVbgoENF+66+/vttWt9VfVSVp2JEqT1ZeeWVTtYuGHp1/wfl2x513uEqTJZZYws4++2zXHw1xUvWNggWFIzq+HouGzF1++eWmcGnzzTd3lVOqkulok5GCHR1nqaWWcrsvvPDCSYdRAKMKGj1G9eWjjz6yBRZYIGkbDVnbY489XKWOAiQFcj5kUjij1rt376R99MOA/gNsxF0jbOVVVm4VEKVufNlll1l1dbU7b37d999/74InnZ8om4JFDS2T8TrrrJN0aD2/dP5807lWCw8VVDXUYostZiuuuKIL1hZZZJGkcEjbl5WX2dBhQ3WThgACCCCAAAIIIDAHCsQiIPr26cs6dGqijYeS7/rbT5J/juqnlvgpqqP9Ncf5LU93UxccN26zYCTPPJMnmIgPK+c/B89EfPA8HS6fr23f5fKmGus9Y7T1nrvSui0RBLclLWuaZzTY+I/us8baKdZz3eP85onvqlRR4PLC8y+4N9kaIvSPf/zDtttuO9OcMBr+te6667qASJUvCituvPFGF/KoskPrFXhssMEG9u2337rhWBrqpIBIYVN52Z//XH/88ccuKNDwqqOPPtoFMpqb6IcffkhUg2jOHc1To8BJVTGqctl9993d0LY333zTzfEzfPhwGzJkiJsbSFVNCqi6d+9uDz74oAt91l1nXZtnnnlcvzQ87JprrnGPV/PW/POf/7RPP/3UBVsKaBRg6Ut90Nw9q622WlJwkoBq58aOO+5oF110kRv+pFBHIY++l5S0nAhV3Wi+HVX36D5UnaMKLZn6phBL/f3yyy/dov79+7vvCsP+85//2N133+0qhDSULrX5gMTvm7re/zx58mRTQCRfnSPfFF7ttNNO1rdvX7/IVfL4IXRa2KtXL1Po5pvOswKxcNP50GNT09Cxfffd193Wfepriy22cI9FYeSyyy7r1ul/CoWuvvpq99zT/fimwExzQKl98sknbuidrxzy2/AdAQQQQAABBBBAYM4W+PMdx5ztwKNHAAEEnED32rHWa64gHOr355t+rSipLrd+85Xarx89Zl0Gr2+V862cJLbWmmu5N/4jHx9pq62+mr366qsumNGbfL2BVwXHkksuae+++65brp01EbNvGoKmIEnzv+iNu0IOVecoxNCQoXBT4KNQ4u9//7ur1FFAk27yas2/c8YZZ7hdFY4oaNEQJM3zo/BA4YuCF1X8aOiR+qeqILUVVwgqTf6oJlGQpQohhV1q/fr1cwGE9vFVPgogFNio6bb6lEtTJZACDA3J0nAzBW8KoxSiyEV9UeCjqio1DQVTMKVQxFcIKRxSOJdqosqhb775xmqCCkjNaSRjVXnl0jSHj4bCKQzzTVVgt956qzvXfpm+6/Fsu+22iUV6TKrg8k2BnSaQDrcZQQWobxpWpgqiY4891vbcc083j5PuVxVa6ke4qdJL277zzjuJxbp/VYX5oXgKt/zQM7+R+qM5mtQUJHlfv57vCCCAAAIIIIAAAsUvUFr8D5FHiAACCGQnUNpcb+X1M6x6rqr0O1SWWo9elVbzVcuEyuGNNDxHVR2ah0gVOAoefGCiEEPDxTR3jpq2U1NIofljNMGwKkoUWNTV1bmrUenNvEIFhTGqeknX/BCp2tradKuTluk4ajq+qpQ0NEpXrzrttNNcNdKuu+6atL3/QRMpq18Kj3xTMKWWaQ4d3VdqnzTcTU3z4qQ2LVN/fNP+6peuXqaqJVUMvfHGG65aSQGM2s477+wCKz/PU3hIl67+lXplMe0jL22nkEihiSq3cmkyUSCl4YHh6iGdXwVvqhQLN02yrVDGf/lKLL+NQh09b8JfGg7mm+ZOUqWXKqsUNuqxa56kxx57zG/ivquaSVda0zxC4f1V1aTKLFWHqS2zzDLu8bsf/vif5tHS81D31dZE6OF9uI0AAggggAACCCBQXAIERMV1Pnk0CCDQCYHSpnor0STPZX+MK0tzrIoupdY0+fs0a8xV6GhokIZzqeqlLJjbSE1BjObVOfPMM91k1b5qRZUkChV0aXpdfUzzwqipEkZXlNIVtFRldOKJJ7Z6Q+82zPF/Or4CGVU3qQpGwYSGl6VrurqVWvhKWwpX1OYbOJ/7ns3/VKmkK4ilXvFLQ9J0PD9Myl+C3R9Tl6VXU+DkJ6fWnEuaS0lfqpLSHEoaLqemIMlPTu0WBP/TPE263LxvGpY177zztqra8evb+67ATm5+QnK/vYKYdMGUX5/rdwU7qfM56bmlPvimSjE97/S88cGk1vnJqcP9UnWQzoMm8vZNoZacfOjol/MdAQQQQAABBBBAYM4RICCac841jxQBBNoRaC4pcxNSt7VZYzBZdUnXltAkdTtN1qymkGjY0D/nt1H1jQIQLQ8P3dEbdIVEmnRYQ8/8PDS6rSoQvVlXAKIWrrBxCzrxPw0zU5WI7kfD2FShE54jR0O4Pv3sU/vwww/dcCzNXXTJJZe4ih5V9ai6R9us+4/kSpnULingUjim4XZqCrpefvllN7RNE3OrqmrTTTd1x1JFjJosNJ+Q9tH9K4TRfalqSZVBCoN0HC3XlyqtNCRq1KhRLvDQ8Db1NxyoKJA79NBD3f1peJ0ei+x1vlL7qPl+pk0LLnbwxxAvVXlpeJ5vv//+uxsKqOFeCpp807mTo4aARd3koSFjCrl0/4+PfNwURqkyTU3BkMLHf//73+58aaiYvtRvVUypmig8TFEB2yqrrOIev+bNUkD37TffujmrVOVEQwABBBBAAAEEEJgzBVpq/ufMx86jRgABBJIEGkuDK5aVVVjdpDqr7JP+ClwzppVYj5VbLlmftHPwgypz9GZe1UCrr7G6G87lt9H8PXojrkDIN71R1zAlzSOjAEnhiCptvv/+e1fNoZBCw7s0KbSu5KXhP+laaUlL1q8rgPmm4+jLN79OyzQ0TNVKqjTRfSgoUsCguW7WW289O+igg1zooKFQmmxbIYIqjDQ0SU19VTWPrpblAwVfLaX1ui99/fbbb24yaYU6qh6SgYa4aSiUH96lx6wheb5yRZdq/9e//pUYpqX7evjhh52DLkvvK4V0P75p/ib5qL+6otcXXyRf+VKB0s8//5wIVDQMTCGRwjptqwmv1Uf1RUGThlv5pmOrKShSIKRztdBCC7nAyW+j735yal8dFl7n7cPLwrfbW6/zdPPNN9vJJ5/s7lf9VyimeaH0uOSppudJuMk1PDm1X6e5nDQpuoI+XY1OlUNqOq4Cvd12281vyncEEEAAAQQQQACBOUigJLjiSXNbj1erw28y2to223X6VHP0/cdkuznbIYAAAn+ZQFX9FOtR84sNWKCblXVPztCn/lJrjd0WtD47tsyD01anVKmi4T2qfvHDtNJtr38PGxsb3Zvz1PWa66ZXz17u8uPhdQqNNPdNeIhReH17txWMaDJqhTFqGr6lnzXZswIshTUKRKoqq5LuW8O/NNFzuHKmvfvS/EKaEFlBkqqg/JWz9Bh0nHClT/hYCpJ0X37i6fC68G1VQqnqR5Nbq7pJwYeCt3RN/Vc//HxMfhvfR/+z/57J368Pf1fAprmBUifGDm/jb+t5ocfmL0fvl2fzXX3S8ynT72VVh80111xuvSqh1C/NcdTWOdO5kGGm56ncfHXb4EUHJz0nsulzW9vofvWY0gVrbe2XzTo9fs191dZjz+Y4qds0NjTa1GlTM3qlbp/tz/q3QEM5dcXCRQcvmu1uWW2n15Ne0z6IzWqnLDdKfW1nuVu7m2kyefU7dUL1dndsZwM931SJp3+P5p9//na27tjqzry227un8Gu7vW07sl7PO4XU7f1b25FjaltVzOo1qCs6pv6b29FjpW4ft9e2+t/Z39upBv5nfeCj14mGTUf9WsnXa9v3ne8IIDB7CiS/+5k9+2jzDNk0q541WbN7Q1bdtVtW22e7UUPw5m3y1Clu87n79M12t6y2awwCuPrgzVl4otOsdsxio9q6muAP+QorL22ZByWLXbLapLGp0eob6q1LZfJVnrLauY2N1N9pM2a6/vYOPsmOsukcNjY2uDe8UR5Xx5oZvNmrCt6sl4WqNaK4j/rgDWRTc5NVVaSvZMn1PqYHbxJqamsCi0rrUd0yaW2ux0rdr66hzi2qLI+2z/l6bauztfV1pgqcij8mUday5jHv2m+jX7Eu1VXBa7M5eO4Ec7nMaLLyAYtZ76GXaJPImp84ON0Bw5dKT7c+12UamqRhWLpylu5DFU+q1NGk2P4zg3RvFhTu+IAn2/vO9AerqlXaagqOMoVHmfbTsCl9ZWrqe7o3Kpn6mOk46Zb7Kp5066Jc1pHnhAKo1KucpetLe+ci3T4sQwABBBBAAAEEECg+gVgERP1X3ior+Xx9EqlPsaYFE6mq9Q+uNBRly9cnkepjvj6tyNcnkervzGA+lMpg4tX+gwZFyezm5cjHJ5HqZEeqDDryoPL2SWQwWXJ9MDdMt6Cqo3/wiVOULV+fRObrta3HnvaTyODfnKbpvwdXK3vaGiZ+byVV3a3X/Ctb1cLrRMlVsGNpCJEmJX7t1ddcWKghRhrG5F/bBesYd4wAAggggAACCCCAAAIFE4hFQFQwHe4YAQTmWIHS7v2s20q7F+3j15CsXXbdpWgfHw8MAQQQQAABBBBAAAEEOibw54ymHduPrRFAAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQQAABBIpEgICoSE4kDwMBBBBAAAEEEEAAAQQQQAABBBDIVYCAKFc59kMAAQQQQAABBBBAAAEEEEAAAQSKRICAqEhOJA8DAQQQQAABBBBAAAEEEEAAAQQQyFWAgChXOfZDAAEEEEAAAQQQQAABBBBAAAEEikSAgKhITiQPAwEEEEAAAQQQQAABBBBAAAEEEMhVgIAoVzn2QwABBBBAAAEEEEAAAQQQyFrgX//6l916662ttj/uuOPsxRdfTFo+fvx4O+CAA+ytt95KWp7ND8cff7ztscceaTfV8dZbbz37+uuv067vyMKmpqaObJ7Vtg899JDJI7V9/vnnts8++9ivv/7qVjU3NzvLnXfe2XbddVe7/vrrrbGhMWm3yZMnu31uv/32pOXhH+rr623ttde2Dz/8MLH4gw8+cH4//fRTYplujB071k444QTbcsstTefyyy+/TFqf+sPEiRPtoIMOsmeffTax6o477nB90mMJfx144IGJbfyNo48+2j0uPTfC2/rb9913n9tUfZFbarv//vvTWmZ6fH5/mRx55JH23//+1y2aMmVK2vtXP5588km/m/t+7733um39Qp23c845x3bYYQfbeuut7ZRTTrFvv/nWr57tvhMQzXanhA4hgAACCCCAAAIIIIAAAsUn8OOPP9pJJ51kCjd8e/fdd+3iiy+2c8891y9y3/WG/6abbrIBAwYkLc/mB4VLv/32W9pNf/jhB3v55Zczrk+7U5qF7733ni255JJp1nRu0TvvvGMKUVLbL7/8Yrfddpsp9FE76qijbN9997Xq6mr3deKJJ9o5556TtJuOo30UmCn0SNdGjhxpEyZMsBVWWMGtVuCx0kor2Z133mmTJk1K7DJt2jRbbbXVXGC36aabuqDq73//u33//feJbcI3FIJo+xtuuMEUkvhWWlpqDQ0NSV+vvPKKPfPMM34T933cuHF22WWX2RZbbGGjRo1yj6OiosKqqqoSX5WVlW5b9VXnI7W9/fbbdvfddyctzvT4whvJ8oorrjDtr1ZWVmaNjY1JfVZoJFuFZuGm5/FGG23kFul5vfTSS7vtunTpYn369LERI0bY4MUGp+1v+DiFul1eqDvmfhFAAAEEEEAAAQQQQAABBIpEoCmoXikta/PBbLXVVvbYY4+5N8errLKK21YVF2qqEvn999+tX79+7mcFF8suu6wtvPDC7ueo/rfjjjvaNttsYwobOtO++eabzuzeqX0VsA0fPtwFa8ccc4w71nnnnWcKIcJNAdvQoUPt8ccfN3nqcac2BTj77befW3zLLbe429tvv7098MADSZu+8vIrpnDt448/tp49e7rqrgUWWMAefPBB833wOyj023zzzW3FFVc0BT3htttuu5m+fJs6daotuOCCdtZZZ/lF7rsqzdT3gQMHJpZfe821Vlbe9nMssXGaG209Pr/5o48+apdeeqkNGTLEL7Lu3btbahXW2Wef7cIyVW/5puo0PV45v/baa6566bDDDnPH8883VXk99dRTtsKQlkDO7zu7fKeCaHY5E/QDAQQQQAABBBBAAAEEEIiRQP2YD23yI/+yX69ez8Zcua6Nu3moTX/juoyPYNiwYW6dQiI1BR0axqOqi169eiVCiZqaGhdqKFBS0zarr7669e7d2w2H8kObVL2h4U6qMlpuueXsxhtvdNvrf9OnTzcNW5pvvvncOj+E7f3333dv/lVlpP1VOaOhSjq+jqNhQ3V1de443333nW2yySZu+eKLL+6Opzf8Cmc0bEjD1LRew67UdOwNN9zQba/jatiXbxpadPrpp7tjDBo0KKlPfpuOfldY44e5zT333C7I8Md444037JNPPnFDmlTRoiAotan6R8O/9tprL7dqQP8BNuKuEabwI7XNO3Bet8hXMOkczZgxw+aff/7UTd25lNMTjz/hqn1abRBaoCohVUFpOKFvel4o3NLwtChbW49P96MAbM8997Sbb77ZdL4zNVUPXXLJJa4ySxVNvul8a2ijKpv+85//2CKLLJIUDmk7BVxDhw3tVNDl7y8f32NRQfTxra3HI7aFMaatlZ1c9/H/OnmADLv/lGH57Lz4xzx1rjY47sd5Ona+Dvtzvg6cx+P+HhxbX3Fq+Xxt58tBr+24vb7z9drOl7GOm6/P8OL42s7n6+TTPP0OzGef8/W8y9fr2v0OfDU/vY7jazs/EmZynpivg+fpuHF8nfwaWOgrTi1fr+18Gnyfz4O3ceyq+snWY9Yv1nvuKuu9UPAmubyrNc6stymfP2QTf3jD+u44vFVFkUKMddZZx4U/qhhRxYXelO++++6uMkVVK4ccckhiPiLN2aKgR2/G9957bxs8eLBpXhpt88ILL9jJJ5/shjppiJKGr/3885+/uRWQKKzRMCvdl+5Dw7RmzpzphjxpmJOGXX300Ud2/vnnuy/15eCDD3ZDoxRYKKBQEKL7eumll+zwww+3iy66yIVJqm7SHDuHHXqY9enbxw1ZW3/99d2QqDPPPNMefvhhd6yuXbu60EFzB2lIkkIiVbKEW+EAmAAAI4VJREFU+9QGc9pVJSUlrm8KIV5//XVXkaMqIFX0+KZASH3UMC+t0xw4CoQUTvmmbXbaaSfr27evW6TgQs0HcO6HP/6nYWcy0fnbeOONTb6qslGlUGpbbLHF7NRTT01d3OpnhU0KiGQSrn567rnn3JAuDWULtxdefMEN99IyBTNrrbVWYrXmFQoHhFqhcxtubT0+hYLbbbed6T41dC91yFv4OKow0v2H501SaKQqI2+ncE7BnK8cCu8/O9+ORUA0OwPSNwQQQAABBBBAAAEEEEBgThIobW6wHjVjrN983ayiV0XioZf1KLe+Pcwm/zLWpr5wrvXc8LTEOn9DFT8akqTARJU7CjAUbCjA2HbbbV3QouFQCy20kJsLR/u9+eabfndX7XLGGWe4nzUvjgIkffl5X/yG+tlvp/WqRlLVULqmyg/1Q023fbCgSZoViihoWWKJJVxApKBEt5dZZhk3VM6HDgqx1FT5ovBgjTXWMA210vFUlaKWqU8KzjraLrjgAvvHP/6RCJsUOGnokiqaFLwoGNMcOGoy79+/vwtQNBRNTeGYhnFpIudsmuxUsaT5eGpra913hSoaUqUhWLk09U9WCuXCTcGVqpo0X1G46Tnim4aeffbZZ/5HF8z4yjS/sCMTkR977LEu8EsNmfyx/HfNy6T5iU477bSkUEvhkM6HrzzSOdBQvHDTfgq/1Pbff393XsLrZ4fbyeKzQ4/oAwIIIIAAAggggAACCCCAwGwrUFU3ybp0K08Kh8Kd7T13uU3/9FlrrpsZXuxu+3lw9GZeFUP+Tb8qNzTMTMPJFHT44WjaSVUyGrqz6qqrukoff1ANTZtrrrnccr0xT72ilN/Oz2vkh4755em+a1sFIGq6T1UXKWhSvzVkyAdJqft+/0d1jgIP3/72t7/Z6NGj/Y9J3zP1ScOtFLqEJ/LWjn5ol9b7pjBIhtpeIZb6quYnuZaHqmJ22WUXt1yVS36yas21I+91113XrWvvfxpWpyothWe6reF0mrhZVUy5NFVfXX311W6YVrh6SMP+NEG5nxcpfOwJ4yc4B1mEwyFto0oozbUU/tJjz6apku2qq65yw+U0FFCVYrqq2/PPP+/6Fz6GKsjShVoKBsMVRQoQVUUUbhp2qKo2PZ87El6Fj5Hv2wRE+Rbm+AgggAACCCCAAAIIIIBAEQlUNs+ybt3LMj+iqlIrD+ZhaRj/ZattNMRJQ5NU8TJmzBjThMhqGoqlChtVuGiol59/SEOZdBl2vblWGOJDEO2jUEiTASuE0RtyH4RoXRRNIY4mUFbVjAIrDWPq1q1b2kNrriNVRfk5gbTRp59+6vZPu0OGhZpLSe3VV19N2kJzBSnQ8cPIfNCjjdRPzaHkAzCFFRp+pYosDXvTl656piBJ1VlqqtIJBxpuYRv/UyAkb18tpOqe5Zdf3gVFbeyWcZWCJQUtGi4YbgqxNDl1urmNwttFeVtDF//97387sx49epi+ysvLXf/CVUCqorr22mvd5NN6vvqm56hsffip5ara0jA1hU++KYw74ogj3Fxaftns9r1wQ8xKgmyquWl286A/CCCAAAIIIIAAAggggAACbQg0W0nwXq6NDdyqYIOS9CGSLl2uiZA1DCscBKiaSCGQghDNdaOmYVoaHqVhSJo3KFwlpEoTXWpelxLXVdEU4ETZ/ve//7mqGU1WreFsd911l6t40pAw9VGhgKpqNOxI8yWddNJJ7nGpAkV9e+KJJ7KqsNFQL11SXcPRNttsMxegyULDyBZddFF3HIU+MtP8QwoedDU2hWka1qShXxoqdtxxx7nheKpckVPqHD664piCIYVtqo655557krgUbmniaX2p6buG5ykUUkh36KGHugovzTv0zjvvuNu+0kdzQylo8YGP9vNVUKrIkp+qnxQs6Wp1F154ofsKBy1+cuqrr7o6qV9R/ZDp8el8pk7MrQofhZkKjnxTn9OFWjpvmuNKk1P7pkBOz2VVeSkMVUjXtUtXe/W1VxPVYH7b2el7wSqIKrr1mp0c6AsCCCCAAAIIIIAAAggggEAWAvWl3Wz6tOCy9pnazMZg+FGzVczTUg2TupnCFLXUIUB+wmNVkPjJfbWtQoR5553XBRCq6PFN4ZFCJl3d7Morr3TVHVqnEEVfvimUUNOy8G2/XhVCvmm9vhRWTJgwwVWVzDPPPO5NvQICH7qowknD2zRXkuaj0bxECmA0EbUmfdY8SwoJ9KXWVp80L829997rhn/pcWuInR6XrqimqhMN6VLl1EknnuSOtd5667lJqjUZtiqndt55ZxdQHH/88W4OJA2F00TSqU0hmyqRNPePhmSlzn00atQoN2/Oyiuv7HZVH1RNM2vWLBfQ6fgKhFRFpX7JwAcrmk9Kw7t8U/Cn6huFaLrKm277iZ91/3JT4BRuWq/KqM023yy8OHG7pPTPc5pY+McNf17Dy8PPAS1v6/GF90t3W9VDmp9IgVG4ikzD3e688043p1B4P51HVbdpDitNaq3KsMGLDXbPB4V/u+22W3jz2eZ2SfDEbzP71epU2M72XmnijO/etN/e+W9nD8X+CCCAAAIIIIAAAggggAACf6lAs8094yubK5hrqGquP+fccV1oCoKVn+qtaqlh1n2twyPpld6T6g163z59W10eXOsU5KSGHZ29Y81Bs+KKK7qJi/v06eMOd/nll5u+NN+Qb6qMUYgSbnq/q2qZbN9HKxRRlY0fvuWP1djQaFOmTklcZcwvD3/X1bMUvmR7X9pXwZJCJT+hcvh47d2WtyqmFMqFQxldHU6hiA/22jtOuvV3j7jbXRHOh3DptpndlimkVBDkg8BM/dN5UgWTfy5l2q7QywsWEOnJNO7VG2zqmC8yGlRUdbeqPgMyrk9dETxXrbG+IRjvGu3IuabgH7mampaJyrp165J6t536WX1uCF74FRV/ptadOmBoZ00aVhIM5SttI2kNbZ71zebAo6mxycoi7rP6W1tb7/6h6dLlz/K8rDvWxoY6h/rHrKws+qK5huA5VxaMUQ19SNFGT7JfpX9A9PyIus91dQ2uPFeflFRV/XnViex7lnnLxuB5oVy/NGLnfL229Uhcn4NOh3/BZX6E2a+pqalzvwRaflFG+/qO22vbOdc3Wml5aYf+eMlGe9asGvc66VJVGfnzLm6vbXk1BK/v8orgd2DmD9iyYU3aRv92zprV8juwa9cukf5bF8fXtvqcj9/b+reotrbOvUa6dk15s5V0Rjr+Q75+b6snjXl6befr97Ze13X6vR38nqoK/t2IsuXr97b6mI/Xto6br9/b+ntOv6s0h0dlxH+X5+v3dr5e2y3O+fmbXO9P9FrRcznqvxezeW2XzhxnVb++a12ry6x7rzIrD15XtTWNNnVyUDk03xDrPazlql4yiGNT1cyaa67pAiLN76NQRPPJ3BoMBwvPNRPHx0afEcgkUNCASKVZU797z35990GrnTYh0ceKbj1trsXWtP4rb5VYls0N/WLWbOhRp8c1NTVuojT1QaWDUTalwlOnTc1LkqiUUjPCh2fSj6LvmvxMCbFS4yjb/7d3HrBy1N4ednoCIaG3BAgk1MCTaAKeAFGF6PAAIRRAIAXRe6969B5Ch9BCCaK3P70pQOjw6BAgEAKE3nsL8/wz8jK7d2fsm929mdn7Wbp3Zmc8ts/ncjxnXJRerVqv9A4bNqyZQVsD329u0bb0ImPNikBlbvCgwR2+ZjQavpSS5lnXfpFoNFwN85SCEwsN1W2m09cSGVrSwy6bEX6r6rbSpq8+6sCm5z83I82aj6481Jz5Zn8pKFvdFk8Nv9WXtEa+KtXLF80PV/nQgo3NLndlq9vio6+rGu7ema+I9bimr+mrpt+BZcTwEU1t68pYt1ult9V+Tp8+3b5g9TbDRwxPZ0HD563S20pYq+p2q/S26rXWvlB7lF73pGHINoBW6W2lrRV1W+G2Sm9r4WHpV02z8Ts1Kb5muFbp7VbVbcncKr09depU179Vf67ZfdzYuv33T1+an//vGvPb1KeM+f1X03vOhcyAkZuZ/iP/mULWjDyflWFIT02cONFMmTLFvWNqXaRmv2vOSvmIGwK1BJo71KY29IjfgxZf2Qwcspz55Yt3zS+fv29mX2gZM2C+YaZnn+aO1IlICl4gAAEIQAACEIAABCAAAQhAIJJAz4HzmTnWPtj9RT5SKm/68KnFhfWHg0B3IDDLDUSC3LPvADNw6H+5v+4AHRkhAAEIQAACEIAABCAAAQhAAAIQgECRCDR/QZYiSUdaIAABCEAAAhCAAAQgAAEIQAACEIAABIIEMBAFEeEBAhCAAAQgAAEIQAACEIAABCAAAQi0NwEMRO2dv0gHAQhAAAIQgAAEIAABCEAAAhCAAASCBDAQBRHhAQIQgAAEIAABCEAAAhCAAAQgAAEItDcBDETtnb9IBwEIQAACEIAABCAAAQhAAAIQgAAEggQwEAUR4QECEIAABCAAAQhAAAIQgAAEIAABCLQ3AQxE7Z2/SAcBCEAAAhCAAAQgAAEIQAACEIAABIIEMBAFEeEBAhCAAAQgAAEIQAACEIAABCAAAQi0NwEMRO2dv0gHAQhAAAIQgAAEIAABCEAAAhCAAASCBDAQBRHhAQIQgAAEIAABCEAAAhCAAAQgAAEItDcBDETtnb9IBwEIQAACEIAABCAAAQhAAAIQgAAEggQwEAUR4QECEIAABCAAAQhAAAIQgAAEIAABCLQ3AQxE7Z2/SAcBCEAAAhCAAAQgAAEIQAACEIAABIIEMBAFEeEBAhCAAAQgAAEIQAACEIAABCAAAQi0NwEMRO2dv0gHAQhAAAIQgAAEIAABCEAAAhCAAASCBDAQBRHhAQIQgAAEIAABCEAAAhCAAAQgAAEItDcBDETtnb9IBwEIQAACEIAABCAAAQhAAAIQgAAEggQwEAUR4QECEIAABCAAAQhAAAIQgAAEIAABCLQ3AQxE7Z2/SAcBCEAAAhCAAAQgAAEIQAACEIAABIIEMBAFEeEBAhCAAAQgAAEIQAACEIAABCAAAQi0NwEMRO2dv0gHAQhAAAIQgAAEIAABCEAAAhCAAASCBDAQBRHhAQIQgAAEIAABCEAAAhCAAAQgAAEItDcBDETtnb9IBwEIQAACEIAABCAAAQhAAAIQgAAEggR6B32UzEPv3q0RqV+/fi0j0aNHj5aFXbaAW8m5VSx69iyfnbWMnFtVt1tVLhRuGTm3kkcrwu7Tp08rgnVhlrFut6qetLIstyrNLSsYNuBW6e1Wcm4lj7KFXUbOZawnZeTcqrrdqjqi9JaRc6t4EC4EINA4gR6JdY0HQwgQgAAEIAABCEAAAhCAAAQgAAEIQAACZSVQvqEPZSVNuiEAAQhAAAIQgAAEIAABCEAAAhCAQEEJYCAqaMaQLAhAAAIQgAAEIAABCEAAAhCAAAQg0FUEMBB1FWnigQAEIAABCEAAAhCAAAQgAAEIQAACBSWAgaigGUOyIAABCEAAAhCAAAQgAAEIQAACEIBAVxHAQNRVpIkHAhCAAAQgAAEIQAACEIAABCAAAQgUlAAGooJmDMmCAAQgAAEIQAACEIAABCAAAQhAAAJdRQADUVeRJh4IQAACEIAABCAAAQhAAAIQgAAEIFBQAhiICpoxJAsCEIAABCAAAQhAAAIQgAAEIAABCHQVAQxEXUWaeCAAAQhAAAIQgAAEIAABCEAAAhCAQEEJ9Ppf6wqatplK1tSpU83ll19ubrvtNvPpp5+a5ZZbzvTq1atDWEmSOD/jx483zz77rBk8eLBZcMEFO/i7//77zY033mjWXnvtqntPP/20i+ehhx4yf/zxh1lyySWr7vsf8nfFFVeYtdZcy/TsWW2Pm/HXDHP22WebX3/91SyxxBL+kVIcY+X/6aefzJVXXmkmTJhgJk+ebBZffHEzcODAKhn//vtvM27cOPPuu++aFVZYoXIvNo/0wA033GCeeeYZs9JKK1We9yfvvPOOueiii8xcc81lFlhgAX+5FMdYzjHl/vvvvzdnnXWW6dOnj1lkkUUq8sdyziqvqgPK4wceeMD88MMPZpllljE9evSohF/0k1j5JUcM57ffftuMHTvWjBw5sqqsh+rCrbfe6tqaiRMnGv099thjlXajM2ksKu+Q/Ol0x5T7LM6hPFIdUFn1nN9///1Ku6G26PbbbzfSC4888oiZMWOGGTFiRDpphT8Pye8FiC1TWTowL56vvvrKnH766RXGYq341P7Xury2u9ZvkX7nyZ9OZwznLB2ocPLqwosvvmguvfTSKs6LLbaY03XpNHQHHRjTvsysDgxxRgdW97GzOMfkkcptVnnNqwvp8l7U85i2wKc9RtYsHRjL+eGHHzbXXnutWXHFFU3//v1d1O2gAz1DjhCAwMwTqLZY1AnnzTffrHO1mJfUWG6//fZOuQwbNsxcffXVZvTo0UYNXq079thjXQd2oYUWMt9++63ZaaedzBNPPFHlTR2vww47zNx5551V12V82n333d21Oeec0xx++OHuhbDKk/0hA9X+++/vDFF/zfir9rYZe+5Y1zi/9tprHe4V+UKs/L/99psZNWqUufvuu91Llgw42267rePi5ZNx7IADDjAXXnihee655/xld4zJI3nUi/Qpp5xipOxqnRTlPvvsYz744AOz/PLL194u9O9YzjHl/qOPPnJ5cdNNNzkWacFjOdcrr9dPuN4cfPDBzuAkg4gMI3oxLJOLlT+G86RJk8x2221n7rrrLmcs8xxi6sLNN99spk2b5owSMkwMGjTI9Or9T8c7No0+vqIdY+T3aY4p91mcQ3n0yy+/mGuuucboKMb6UxvundqR0047zRmSBwwYYA499FCjMl4WF5I/LUdMmcrSgaF43nrrLXPfffdVGNdy9unIa7u9nyIeQ/Kn0xzinKcDQ3VB/GS88GVZhvk55pgjHb3pDjowpn1pRAfmcUYHVvexszjH5JEKblZ5DdWFqkJf0B+htsAnO0bWLB0Yy1nvdgcddJCZd9553QdyH3fZdaCXgyMEINAgAWvRznX2a2vu/SLd3HfffRNrkEhsZ8klyxoFEjsiJbFfQKuSaTuv7ro1CFWuH3nkkckWW2yRWGOSu2YbcufHGheSVVddteLPNr7JmmuumYwZM6ZyzX75d36tQahy7a8//0p22GGHxBqe3D09l3b2i6q7vs0221SFlfZTxPNY+ZV2+2XCyWgNcE6U33//Pdloo42SE044wf3+7rvvEmvQc9d0tIYid13/YvJI/j755BOXP+KsMGqdVYDJOuusk/g01N4v6u/OcA6V+zfeeMOV2T333NOxuv766ytix3LOKq/WGJLYUXCV8G655Zaq+lK5UdCTWPmV/BBnO/LElXdrVHBH+xW0InWoLqjNUltljcWVZ/xJZ9LonynaMSS/T29Muc/jHMqjV155xXH2OsLHq6PaI+WBHTlUuSw9oDJeFheS38sRU6aydKDCCMVjR/Em9iOKj67uMdR2132oIBdD8vtkhjjn6cCYuiDGdvSmj67usTvowFD70qgOzOOMDvy3j53HOZRHvvDWK68xdcE/X9RjqC3w6Y6RNU8HxnD+8ccfXb9bfcK0LmwHHeg5coQABBojEBxBZBuPBk1QXfO4RgnZl1iz9dZbV6Zyaaj1KqusYh5//PGqRGik0DzzzGOsoady3Rpq3PSR6dOnu2uaMqapZda4UPGjE9vIGw2flX/vNttsM3eqETLejTlnjBtBYA1P/lLlqJFF+jJtDSWlm1oWK7+E1Ve3TTbZpPKFvm/fvmarrbYyjz76qGOhIa2rr766mx6mvEq7mDzSlCeNXtlwww2NNe6lH3fnGi2jod8XXHBBJQ0dPBX0QiznmHKvKX0777yz49CvX78qiWM455VXTdlTffDOdjCMRuWVxcXIL1liOKtN0eipvffeu4P4obrw8ccfu2cG9B/g2iGVbe9i0+j9F/EYkt+nOabcZ3GOySNNY9WUMX3lVrm26tNHbWabbTajUVzrrrtu1TWFWwYXI7+XI6ZMZenAmHjEeciQIea9Ke8ZtQm1LtR21/ov0u8Y+X16Q5zzdGBMXdBIJk2dFmeNRKp13UEHSuZQ+9KoDszjjA78t4+dxzmUR8rHrPIaUxf0fJFdqC3waY+RNUsHKowYzscff7yLTqNl00tflF0HeoYcIQCBxgkEDUTrr79+47F0QQhff/21i2WRof+uraILWmvFG318MuyXy6o1WHRdnVk5vTTIacrZsssu687T//Ss3JCF//Gvcxk+Fl544Uo8MoBYK76xo4yMpimk3Z9//mkOPPBAs9566zljSfpeGc5j5Pdy6CXMc/XXhg4ZapRXekGQsULTy9JTPLy/mDw6e8zZzgh3zDHH+Mcqxy+++MKcdNJJ7vcee+xhttxyS2NHt1TuF/0klnNMuV900UXNbrvtVtUR8PKHOIfKq6ZaygCrMn300Ue79TD22msvH3zhjyH5vQAxnNdaay2z8cYb+0eqjqG6MGXKFOd/6//Z2pXVNf57DTdFRxdj01gVYcF+hOT3yY0p91mcY/LovffeM2K9+eabGzua0eWXXv7ktDbX0ksvXVk/S4aNe+65x2ywwQY+eYU+xsjvBYgpU1k6MCYeMZaxTeVZ6/epnbBfrX30Jq/trngq6EmM/D7pIc4hHahwsvoamoqjtBx33HGO82qrrebWmPNxdxcdKHlD7UsjOjDEGR34bx87j3Moj/LKa4xe8OW+qMdQW+DTHSNrlg5UGCHOej/RGnzffPON2XTTTc2OO+7opqnq2bLrQMmAgwAEmkMgaCBaeeWVmxNTi0PxX8/69utbFZMaPK03kXbyW2u48Qu02eGdaa8dzn08fm0Q70EdPT2rxv2II44wMlostdRS/nblKKORXjzUqSujC8mflslOKevA2efPH3/+kfba4TyUR1Jy1113nesQ1+alAtPC4HInnniiW3NEL+4asaWvOGVwsZy9P8/Vy1av3Pt76WOIc6i8ah77zz//7NbAkDFJ7vPPP09HUejzkPw+8Y1yDtUFjVS89957zcsvv+zW4pIBQ2ubqRMXm0af1iIeQ/L7NHvOWe2r91fv6J/NqwtaE+7BBx80r776qluEeu6553aGzdrwlF6tzyAD96677lp7u5C/Y+T3CW+kTMXEo40i1NaKsxahfumll8zFF1/sog+13T6NRT3GyO/T3gzOWXVBozVUlp9//nnXbnhdJ75y3UUHStbY9kV+a10oj0Kc0YEd+9i1jPU7lEd55dXXuay6UC++ol0LlTOf3kZlDXHWiHp90D733HNdGzHffPMZO9XMfPnllz4J7lhGHVglAD8gAIGGCAQNRA2F3oUP+92ptHtK2ukLW+2UF+1WVtsY+t/zzz9/+vEO5z6s2mHzGnmkNGjUigxF2hlH003UUZbTzjn6oqrdvNQ4n3/++e6+OhfareCcc87pEFcRL4TkT6dZfj1Xf135IWNcPaOO96NjXh5peK2McOKoBcTFWQuiaqSYzpUXdq0RtyizRg5p5xyNItKOdpqGWAYXy7kz5b6e3HmctSNZXnnVKDCV6/3228+cd9555owzznC78ulabb7Xi7sI1/LkT7cFjXIO1QUN8x46dKgb5aX6oYXV5TTcPDaNReCZlYaQ/P652HLv/aePMXkkQ754yqljvMsuu7hNDdQZ9k67UsqQpLbqkksucSMd/b0iH2Pk9+lvpEzFxKNRodoZVE7trqYAy/ipNiPUdvs0FvUYI79PeyOcY+qCwleZVvvhdZ2McnLdRQdK1tj2RX5rXUweZXFGBxrXTvqyWss2/TuUR3nl1Yef1e9Ox1PU85hyprQ3KmseZy0ZotGd2nxnjTXWcB+xTz75ZPfOosXuvSurDvTp5wgBCDROoG0MROokaXirjC3eqZF74YUXOkwV0zbcaiQ/++wz79U89dRT7lxh5Lnhw4e72+l4ZOSRUUhT0rTGjqbzzD777G5NCz8ySb9l0JChQltKaq6v/nr16uWmqNXuPJKXhll5LyR/Om3i8eSTT6YvOc4xu4nl5ZEMPprGpOGxnqPyv3fv3m5bcY2eEXddq3W6VwYXy7kz5b6e3HmcNXowr7xO/2R6pdz7sBWe3OTJk/2lQh/z5E+3BY1y7mxd8AY27TASm8Yig46VP7bc15N1ZvLoqy+/qmorpDM07VVrQumLdr3pr/XiLsK1zsjfSJnqTDyei4xtKssaARBqu/0zRT12Rv5GOHe2LshYoRGH4izXXXSgZI1tX+S31nU2j9Kc0YH1+9i1jPU7lEd55bWzdaFe/LP6Wmw5a1TWPM5qu+T8Uefa+VBO/We5MutAJwD/IACB5hAIrXGd3okn5HdW39cOStqFxm6rntjRJImd5uV+awcru5Buot2FrJU80Q5jdmHjZPTo0cnUqVMTa8RwOy+ld2OSHzt1Jplw3QR3T+fWCOREtFM/3PN22Hwy5d0pbpebejtoybPCV5r8s7WMlKb0jmi194v4O0/+119/PTnkkEMS+0JQkV3y2al3iXbPEgu7Zk1FLHERW+1coZ1hdC72MXlUCcSe2NFZVbuYXXXVVS7fnn3m2cR2mhP9VtzK/7K4PM533HFHYqcUOFHyyr2X1U6zdGy1A59YiLN2r+gs53R51fN2apTbrchOK3O7QNmpk46z8r8MLiR/ZziLh7jaUT+OgXbMEnc53w5k1QW79Xpy2WWXuXbLGpzdDoja4XBm8qiI3EPyn3nmmYndFMAlPa/cy0Me57y6oOe0g+J/7vpPovKqdkj14dRTT3Xx2lFEiXZ1UTshP/YjgPuzHxnc/TL8y5O/mTowLx774SU56qijErW9djRnYqeYOaZ2w4C6CGvb7rqeCnYxT/60Dgy1LxIrSwfqXl5dEF87Yjmxxvjkww8/TNT2asdV6Vo5tfP63a46MM051L6Ix8zqwDzOalPQgStU7RKbxTmUR6HymlcXlL9Fd6G2oFk6MMTZGuhdX9mux+faZ/XppAft6KykHXRg0csB6YNAWQhoF5dcZ3cFy71fpJt2d5Fk3LhxrlOkTr4dYZLYYe0uib7RHD9+vPs9bdo09xImf/rTy7a1nFfEkWHI3/NHbVsvJwWoRtVfl6HJLrBXeTZ94uNtJwNRnvy+46zOm5xewmSMEyspId1PO21t7zn6o51u57yE8igdTu1LhpSxXvx8mOrE6aWvTC6PszpL6vxLzrxy7+WVX8/CHydNmuRud4Zz2kCkh2XMkHHUh6k6J4NrmVye/J3hbBc7rnDwPHT0dT+vLtjFkCv1RM/IUKGXa+/y0uj9FP2YJ7/qpzqucnnlXvfzOOfVBdWVsWPHVvSDOKuN8PmjPEjnW/pcxpUyuDz5vS5qhg7Mi0cGeX0k8PzU7qt9znK1bXeWvyJdz5O/VgeG6m6eDsyrC9Kxo0aNqnBWXy1tzGx3HVjLOa99UdmZWR0Y4owO/KeP7etnHue8PAqV17y64OMu+jGvLWiWDhSDPM56V/EfQtRG60OU76+3gw4sehkgfRAoC4EeSmhzxiIVJxSJZDv9Hda50W4UWnAw7bQgnIZbprd6TN/PO9eivIpLu5h1R5clfxbn0LpDWQwbySPbkTdKz6BBg7KCL/z1epw1zF0LfaeZZpX7WAEb4az6pvjT6YmNtyj+6snfKs5ZnLTGgqajZk2FrJfGovCLTYdkqJVf5adP7z5uCpIPp1659/dCx7y6oDZB62tpjRw/vD4UXtnuZ8mf1TbPrA7Mike8NFVBee3XIiobw5j0ZsnfbM55dcG+OFfWLquX5nbVgZI1i3Nt+1KPS71ree1riDM6sB7R+tfq6QDvM1Re8+qCD6Pox3rlrNk6UAzyOKt9tiPgOujiorMjfRCAQNcQaEsDUdegIxYIQAACEIAABCAAAQhAAAIQgAAEINAeBNpmker2yA6kgAAEIAABCEAAAhCAAAQgAAEIQAACXU8AA1HXMydGCEAAAhCAAAQgAAEIQAACEIAABCBQKAIYiAqVHSQGAhCAAAQgAAEIQAACEIAABCAAAQh0PQEMRF3PnBghAAEIQAACEIAABCAAAQhAAAIQgEChCGAgKlR2kBgIQAACEIAABCAAAQhAAAIQgAAEIND1BDAQdT1zYoQABCAAAQhAAAIQgAAEIAABCEAAAoUi8P83c4RyLEDWcgAAAABJRU5ErkJggg=="
    }
   },
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Parse divergence tree and output tips with evidence of onward transmission using nextstrain jsons and run a logistic regression\n",
    "\n",
    "In addition to the beast analyses using the structured coalescent model to estinate differences in transmission patterns among epidemiologically defined groups in Washington, we wanted to use a simpler metric that could be used with any divergence tree. Trevor came up with the following idea: if you sample from a transmission chain densely enough relative to the substitution rate, then you will end up with some tips in the tree that lie on internal nodes, ie., their divergence values are the same as the internal node ancestor. This means that you could have sampled the true ancestor to tips downstream of that internal node. Although this is probably unlikely for the mumps data we have we still expect that tips sampled closer to the internal nodes should be less removed from the true ancestor. You can then start to interrogate the tree using tips that are like this and have descendants vs. tips that don't have descenants. For example, in the tree below, Washington.USA/15.17/FH146/6, Washington.USA/3.17/FH67/G, and Washington.USA/2.17/FH84/G all lie on an internal node, and would be classified as having descenants. All other tips in this tree would not. \n",
    "\n",
    "![Screen%20Shot%202021-03-16%20at%202.03.57%20PM.png](attachment:Screen%20Shot%202021-03-16%20at%202.03.57%20PM.png)\n",
    "\n",
    "The code in this notebook will read in an auspice json produced by the Nextstrain pipeline (a v2 auspice json, although it is easy to adapt to v1 jsons). It will then iterate through the tree, and for each tip, assess whether that tip lies on an internal node, ie., has a branch length ~0. If it does, calculate the number of tips downstream from that internal node, minus the current tip, and store that value. Also calculate the downstream branch length. The output is a dataframe containing each tip, and its number of descendants. \n",
    "\n",
    "In the 2nd half of this notebook, this dataframe is used to look at 3 factors we believe might be important to mumps transmission: community status, vaccination status, and age. All of these factors are evaluted with a logistic regression run in R."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "import sys, subprocess, glob, os, shutil, re, importlib\n",
    "from subprocess import call\n",
    "import imp\n",
    "\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import scipy as scipy\n",
    "from scipy import stats\n",
    "\n",
    "import json\n",
    "import collections\n",
    "from collections import Counter\n",
    "from Bio import SeqIO\n",
    "from Bio import Seq\n",
    "import Bio.Phylo\n",
    "\n",
    "import rpy2\n",
    "%load_ext rpy2.ipython"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Functions from augur to read in nextstrain jsons and convert to trees\n",
    "\n",
    "These functions are copied and pasted directly from nextstrain augur/base/io_util.py"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# function to use the json module to read in a json file and store it as \"data\"                \n",
    "def read_json(file_name):\n",
    "    try:\n",
    "        handle = open(file_name, 'r')\n",
    "    except IOError:\n",
    "        pass\n",
    "    else:\n",
    "        data = json.load(handle)\n",
    "        handle.close()\n",
    "    return data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# original code that Trevor gave me for parsing through tree jsons and returning descendents\n",
    "def all_descendants(node):\n",
    "    \"\"\"Take node, ie. dict, and return a flattened list of all nodes descending from this node\"\"\"\n",
    "    yield node\n",
    "    \n",
    "    # this will recursively return all internal nodes (nodes with children)\n",
    "    if 'children' in node:\n",
    "        for child in node['children']:\n",
    "            for desc in all_descendants(child):\n",
    "                yield desc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Biopython's trees don't store links to node parents, so we need to build\n",
    "# a map of each node to its parent.\n",
    "# Code from the Bio.Phylo cookbook: http://biopython.org/wiki/Phylo_cookbook\n",
    "def all_parents(tree):\n",
    "    parents = {}\n",
    "    for clade in tree.find_clades(order='level'):\n",
    "        for child in clade:\n",
    "            parents[child] = clade\n",
    "    return parents"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "def annotate_parents(tree):\n",
    "    # Get all parent nodes by node.\n",
    "    parents_by_node = all_parents(tree)\n",
    "\n",
    "    # Next, annotate each node with its parent.\n",
    "    for node in tree.find_clades():\n",
    "        if node == tree.root:\n",
    "            node.up = None\n",
    "        else:\n",
    "            node.up = parents_by_node[node]\n",
    "\n",
    "    # Return the tree.\n",
    "    return tree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def annotate_parents_for_tree(tree):\n",
    "    \"\"\"Annotate each node in the given tree with its parent.\n",
    "    >>> import io\n",
    "    >>> tree = Bio.Phylo.read(io.StringIO(\"(A, (B, C))\"), \"newick\")\n",
    "    >>> not any([hasattr(node, \"parent\") for node in tree.find_clades()])\n",
    "    True\n",
    "    >>> tree = annotate_parents_for_tree(tree)\n",
    "    >>> tree.root.parent is None\n",
    "    True\n",
    "    >>> all([hasattr(node, \"parent\") for node in tree.find_clades()])\n",
    "    True\n",
    "    \"\"\"\n",
    "    tree.root.parent = None\n",
    "    for node in tree.find_clades(order=\"level\"):\n",
    "        for child in node.clades:\n",
    "            child.parent = node\n",
    "    # Return the tree.\n",
    "    return tree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "def json_to_tree(json_dict, root=True):\n",
    "    \"\"\"Returns a Bio.Phylo tree corresponding to the given JSON dictionary exported\n",
    "    by `tree_to_json`.\n",
    "\n",
    "    Assigns links back to parent nodes for the root of the tree.\n",
    "\n",
    "    >>> import json\n",
    "    >>> json_fh = open(\"tests/data/json_tree_to_nexus/flu_h3n2_ha_3y_tree.json\", \"r\")\n",
    "    >>> json_dict = json.load(json_fh)\n",
    "    >>> tree = json_to_tree(json_dict)\n",
    "    >>> tree.name\n",
    "    u'NODE_0002020'\n",
    "    >>> len(tree.clades)\n",
    "    2\n",
    "    >>> tree.clades[0].name\n",
    "    u'NODE_0001489'\n",
    "    >>> hasattr(tree, \"attr\")\n",
    "    True\n",
    "    >>> \"dTiter\" in tree.attr\n",
    "    True\n",
    "    \"\"\"\n",
    "    node = Bio.Phylo.Newick.Clade()\n",
    "    node.strain = json_dict[\"name\"]\n",
    "\n",
    "    if \"children\" in json_dict:\n",
    "        # Recursively add children to the current node.\n",
    "        node.clades = [json_to_tree(child, root=False) for child in json_dict[\"children\"]]\n",
    "\n",
    "    # Assign all non-children attributes.\n",
    "    for attr, value in json_dict.items():\n",
    "        if attr != \"children\":\n",
    "            setattr(node, attr, value)\n",
    "\n",
    "    \"\"\"The 2 lines below needed to be altered to work with auspice v2 jsons. The original\n",
    "    code for auspice v1 jsons is listed here: \n",
    "    node.numdate = node.attr.get(\"num_date\")\n",
    "    node.divergence = node.attr.get(\"div\")\"\"\"\n",
    "    \n",
    "    node.numdate = node.node_attrs[\"num_date\"]['value']\n",
    "    node.divergence = node.node_attrs.get(\"div\")\n",
    "    \n",
    "    if \"translations\" in node.node_attrs:\n",
    "        node.translations = node.attr[\"translations\"]\n",
    "\n",
    "    if root:\n",
    "        node = annotate_parents(node)\n",
    "\n",
    "    return node"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Read in metadata, find proper parent node, and add branch lengths"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 198,
   "metadata": {},
   "outputs": [],
   "source": [
    "def read_metadata(metadata_path):\n",
    "    metadata = {}\n",
    "    \n",
    "    with open(metadata_path, \"r\") as infile: \n",
    "        for line in infile: \n",
    "            if \"strain_name\" not in line: # skip first line\n",
    "                new_strain_name = line.split(\"\\t\")[0].replace(\"[G]\",\"/G\").replace(\"MuVs/\",\"\")\n",
    "                old_strain_name = line.split(\"\\t\")[2]\n",
    "                ID = line.split(\"\\t\")[1]\n",
    "                age = line.split(\"\\t\")[3]\n",
    "                vaccination_status = line.split(\"\\t\")[4].lower()                    \n",
    "                community_status = line.split(\"\\t\")[5]\n",
    "                \n",
    "                metadata[new_strain_name] = {}\n",
    "                metadata[new_strain_name]['age'] = age\n",
    "                metadata[new_strain_name]['vaccination_status'] = vaccination_status\n",
    "                metadata[new_strain_name]['community_status'] = community_status\n",
    "                \n",
    "    return(metadata)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 188,
   "metadata": {},
   "outputs": [],
   "source": [
    "def return_marshallese_status(strain_name, metadata):\n",
    "        \n",
    "    if \"Washington\" in strain_name: \n",
    "        community_status = metadata[strain_name][\"community_status\"]\n",
    "    else: \n",
    "        community_status = \"Not_Washington\"\n",
    "        \n",
    "    return(community_status)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 189,
   "metadata": {},
   "outputs": [],
   "source": [
    "def return_vaccination_status(strain_name, metadata):\n",
    "        \n",
    "    if \"Washington\" in strain_name: \n",
    "        vaccination_status = metadata[strain_name][\"vaccination_status\"]\n",
    "    else: \n",
    "        vaccination_status = \"Not_Washington\"\n",
    "        \n",
    "    return(vaccination_status)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 190,
   "metadata": {},
   "outputs": [],
   "source": [
    "def return_age(strain_name, metadata):\n",
    "        \n",
    "    if \"Washington\" in strain_name: \n",
    "        age = metadata[strain_name][\"age\"]\n",
    "    else: \n",
    "        age = \"Not_Washington\"\n",
    "        \n",
    "    return(age)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 191,
   "metadata": {},
   "outputs": [],
   "source": [
    "def return_proper_parent_node(node, cutoff):\n",
    "    \"\"\"given an internal node, traverse back up the tree to find a parental node with a \n",
    "    real branch length (basically, collapse the polytomy). This is necessary for most \n",
    "    tree software, including iqtree and treetime, which both normally atempt to resolve\n",
    "    polytomies, resulting in a fully bifurcating tree with lots of very tiny branches\"\"\"\n",
    "    \n",
    "    #print(node, node.length)\n",
    "    if abs(node.divergence - node.up.divergence) < cutoff: \n",
    "        \n",
    "        #print(\"going up 1 node\")\n",
    "        if node.up !=None:\n",
    "            parent_node = return_proper_parent_node(node.up, cutoff)\n",
    "        \n",
    "        else:\n",
    "            #print(\"root is proper parent\")\n",
    "            parent_node = node\n",
    "    \n",
    "    else: \n",
    "        #print(\"current node has proper length\")\n",
    "        parent_node = node\n",
    "    \n",
    "    return(parent_node)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 192,
   "metadata": {},
   "outputs": [],
   "source": [
    "def add_nodes(node):\n",
    "    \"\"\"Take node and add up branch lengths for total subtending tree from that node\"\"\"\n",
    "    total_lengths = 0\n",
    "    \n",
    "    branch_length = node.divergence - node.up.divergence\n",
    "    \n",
    "    if node.is_terminal() == True: \n",
    "        total_lengths += branch_length\n",
    "    \n",
    "    else:\n",
    "        total_lengths += branch_length\n",
    "        for child in node.clades:\n",
    "            \n",
    "            total_lengths += add_nodes(child)\n",
    "                            \n",
    "    return(total_lengths)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 193,
   "metadata": {},
   "outputs": [],
   "source": [
    "def return_branch_length_in_days(tip):\n",
    "    \n",
    "    branch_length = tip.numdate - tip.up.numdate\n",
    "    branch_length_days = branch_length * 365\n",
    "    \n",
    "    return(branch_length_days)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 195,
   "metadata": {},
   "outputs": [],
   "source": [
    "def return_descendants_dict(tree, metadata, cutoff):\n",
    "    \n",
    "    WA = []\n",
    "    output_dict = {}\n",
    "    not_polytomies = []\n",
    "    polytomies = []\n",
    "    \n",
    "    for i in tree.find_clades(): ## iterate over objects in tree            \n",
    "        if i.is_terminal() == True:\n",
    "                        \n",
    "            community_status = return_marshallese_status(i.strain, metadata)\n",
    "            vaccination_status = return_vaccination_status(i.strain, metadata)\n",
    "            age = return_age(i.strain, metadata)\n",
    "            \n",
    "            # calculate branch length in time \n",
    "            time_branch_length = return_branch_length_in_days(i)\n",
    "            \n",
    "            # now do descendants/no descendants with divergences\n",
    "            if abs(i.divergence - i.up.divergence < cutoff):\n",
    "                if \"Washington\" in i.strain:\n",
    "                    polytomies.append(i.strain)\n",
    "\n",
    "                proper_parent = return_proper_parent_node(i.up, cutoff)\n",
    "                branch_length = add_nodes(proper_parent) - (i.divergence - i.up.divergence)\n",
    "                children = len(proper_parent.get_terminals()) - 1\n",
    "\n",
    "                output_dict[i.strain] = {}\n",
    "                output_dict[i.strain]['branch_lengths'] = branch_length\n",
    "                output_dict[i.strain]['number_children'] = children\n",
    "                output_dict[i.strain]['community_status'] = community_status\n",
    "                output_dict[i.strain]['vaccination_status'] = vaccination_status\n",
    "                output_dict[i.strain]['age'] = age\n",
    "                output_dict[i.strain]['time_branch_length'] = time_branch_length\n",
    "\n",
    "\n",
    "            else:\n",
    "                if \"Washington\" in i.strain:\n",
    "                    not_polytomies.append(i.strain)\n",
    "                    output_dict[i.strain] = {}\n",
    "                    output_dict[i.strain]['branch_lengths'] = 0\n",
    "                    output_dict[i.strain]['number_children'] = 0\n",
    "                    output_dict[i.strain]['community_status'] = community_status\n",
    "                    output_dict[i.strain]['vaccination_status'] = vaccination_status\n",
    "                    output_dict[i.strain]['age'] = age\n",
    "                    output_dict[i.strain]['time_branch_length'] = time_branch_length\n",
    "\n",
    "    \n",
    "    return(polytomies,not_polytomies,output_dict)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Set paths, run it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 199,
   "metadata": {},
   "outputs": [],
   "source": [
    "\"\"\"set paths. This North American tree is distinct from thee tree on nextstrain.org/mumps/na because I retained all \n",
    "Washington sequences in it. Several Washington sequences are highly divergent, so they are dropped from Nextstrain,\n",
    "but retained here\"\"\"\n",
    "tree_path = \"../auspice/mumps_north-america-full-genomes.json\"\n",
    "metadata_path = \"../sample-metadata-complete-2020-10-06.txt\"\n",
    "\n",
    "\"\"\"set cutoff value for defining a near 0 branch length; here I'm using 1e-16. The value \n",
    "that seems to be the minimum set for 0 branch lengths using treetime seems to be 6e-18,\n",
    "so I have just it slightly higher. Alternatively, you can get the same results by setting the \n",
    "cutoff to be 1 mutation, which is 1/alignment length = 1/15265 = 0.000066667. I've tried both \n",
    "and they produce the same results\"\"\"\n",
    "\n",
    "#cutoff =  1e-16\n",
    "cutoff = float(1/15265)\n",
    "\n",
    "# run\n",
    "metadata = read_metadata(metadata_path)\n",
    "combo_json = read_json(tree_path)\n",
    "tree = combo_json['tree']\n",
    "tree = json_to_tree(tree)\n",
    "polytomies,not_polytomies,output_dict = return_descendants_dict(tree, metadata, cutoff)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Convert to dataframe and subset down to Washington tips"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 200,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>index</th>\n",
       "      <th>branch_lengths</th>\n",
       "      <th>number_children</th>\n",
       "      <th>community_status</th>\n",
       "      <th>vaccination_status</th>\n",
       "      <th>age</th>\n",
       "      <th>time_branch_length</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>Alabama.USA/11.17/FH90/G</td>\n",
       "      <td>0.005505</td>\n",
       "      <td>45</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>23.170908</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>Alabama.USA/19.17/FH92/G</td>\n",
       "      <td>0.005505</td>\n",
       "      <td>45</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>20.985626</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>Alabama.USA/7.17/FH164/G</td>\n",
       "      <td>0.005505</td>\n",
       "      <td>45</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>0.549958</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>Arkansas.USA/12.17/G</td>\n",
       "      <td>0.002883</td>\n",
       "      <td>21</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>34.496901</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>Arkansas.USA/38.16/G</td>\n",
       "      <td>0.007011</td>\n",
       "      <td>83</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>Not_Washington</td>\n",
       "      <td>11.630015</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                      index  branch_lengths  number_children community_status  \\\n",
       "0  Alabama.USA/11.17/FH90/G        0.005505               45   Not_Washington   \n",
       "1  Alabama.USA/19.17/FH92/G        0.005505               45   Not_Washington   \n",
       "2  Alabama.USA/7.17/FH164/G        0.005505               45   Not_Washington   \n",
       "3      Arkansas.USA/12.17/G        0.002883               21   Not_Washington   \n",
       "4      Arkansas.USA/38.16/G        0.007011               83   Not_Washington   \n",
       "\n",
       "  vaccination_status             age  time_branch_length  \n",
       "0     Not_Washington  Not_Washington           23.170908  \n",
       "1     Not_Washington  Not_Washington           20.985626  \n",
       "2     Not_Washington  Not_Washington            0.549958  \n",
       "3     Not_Washington  Not_Washington           34.496901  \n",
       "4     Not_Washington  Not_Washington           11.630015  "
      ]
     },
     "execution_count": 200,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = pd.DataFrame.from_dict(output_dict, orient=\"index\")\n",
    "df.reset_index(inplace=True)\n",
    "df.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 201,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "46 109\n"
     ]
    }
   ],
   "source": [
    "# subset dataframe to include only Washington tips and print the number of tips with descendants\n",
    "WA_df = df[df['index'].str.contains(\"Washington\")]\n",
    "print(len(WA_df[WA_df['number_children'] > 0]), len(WA_df))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 202,
   "metadata": {},
   "outputs": [],
   "source": [
    "WA_df.to_csv(\"WA_regression_dataframe.csv\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Manually generate the dataframe \n",
    "\n",
    "this seems like a good thing to read: https://www.datacamp.com/community/tutorials/logistic-regression-R\n",
    "and this https://www.r-bloggers.com/how-to-perform-a-logistic-regression-in-r/\n",
    "\n",
    "\n",
    "I would like the following predictors, coded as the following: \n",
    "\n",
    "**1. community status**: Not Marshallese = 0, Marshallese = 1\n",
    "\n",
    "**2. age**: continuous, leave as is \n",
    "\n",
    "**vaccination status:** \n",
    "for this, we need to break down into 2 categories: \n",
    "3. **vaccination unknown**: known vaccination status = 0, unknown status = 1\n",
    "4. **up-to-date**: up-to-date vaccination status = 0, not up-to-date vaccination status = 1\n",
    "\n",
    "The intercept will therefore reflect the odds of having descendants for a not marshallese person, of age 0, with a known vaccination status, who is vaccinated. This seems like a reasonable baseline. \n",
    "\n",
    "The dependent variable is has descendants. I will code this as 0 for no descendants and 1 for has descendants. \n",
    "\n",
    "\n",
    "I will use the `glm` function in R. This will fit a regression model to the input data. By specifying `binomial`, this will fit a logistic regression model, rather than a linear one. The significance values output are the result of a Wald test of significance. These are computed by comparing: \n",
    "\n",
    "for a model: `glm(y~x1+x2+x3, family=\"binomial\")`:\n",
    "* For coefficient of x1: `glm(y~x2+x3, family=\"binomial\")` vs. `glm(y~x1+x2+x3, family=\"binomial\")`\n",
    "* For coefficient of x2: `glm(y~x1+x3, family=\"binomial\")` vs. `glm(y~x1+x2+x3, family=\"binomial\")`\n",
    "* For coefficient of x3: `glm(y~x1+x2, family=\"binomial\")` vs. `glm(y~x1+x2+x3, family=\"binomial\")`\n",
    "\n",
    "The order of tests does not matter here. \n",
    "\n",
    "This is a nice description of the difference between the results output by the `summary` command vs. the `anova`. For the anova comparison, you evaluate the model during sequential addition of predictors. Here, the order DOES matter. https://stats.stackexchange.com/questions/59879/logistic-regression-anova-chi-square-test-vs-significance-of-coefficients-ano\n",
    "\n",
    "stepwise, this is what's occurring when you run `anova(model.final, test=\"Chisq\")` \n",
    "1. `glm(y~1, family=\"binomial\")` vs. `glm(y~x1, family=\"binomial\")`\n",
    "2. `glm(y~x1, family=\"binomial\")` vs. `glm(y~x1+x2, family=\"binomial\")`\n",
    "3. `glm(y~x1+x2, family=\"binomial\")` vs. `glm(y~x1+x2+x3, family=\"binomial\")`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 203,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/lmoncla/anaconda/envs/LHM-basics/lib/python3.7/site-packages/ipykernel_launcher.py:2: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy\n",
      "  \n"
     ]
    }
   ],
   "source": [
    "# convert age to numeric\n",
    "WA_df['age'] = pd.to_numeric(WA_df['age'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 204,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "21.321100917431192\n",
      "17.0\n"
     ]
    }
   ],
   "source": [
    "print(np.mean(WA_df['age']))\n",
    "print(np.median(WA_df['age']))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 206,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/lmoncla/anaconda/envs/LHM-basics/lib/python3.7/site-packages/ipykernel_launcher.py:2: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: http://pandas.pydata.org/pandas-docs/stable/indexing.html#indexing-view-versus-copy\n",
      "  \n"
     ]
    }
   ],
   "source": [
    "# generate age bins, where a 1 means you are at least 20 years of age, and a 0 means you are < 20 years of age\n",
    "WA_df['age_bin'] = np.where(WA_df['age'] >= 20, 1, 0)\n",
    "# WA_df.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 207,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "vaccination status unknown:\n",
      "total 32\n",
      "with descendants 12\n",
      "no descendants 20\n",
      "\n",
      "vaccination status not up to date:\n",
      "total 13\n",
      "with descendants 4\n",
      "no descendants 9\n",
      "\n",
      "vaccination status up to date:\n",
      "total 64\n",
      "with descendants 30\n",
      "no descendants 34\n",
      "\n",
      "age at least 20:\n",
      "total 50\n",
      "with descendants 17\n",
      "no descendants 33\n",
      "\n",
      "age under 20:\n",
      "total 59\n",
      "with descendants 29\n",
      "no descendants 30\n",
      "\n",
      "is Marshallese:\n",
      "total 57\n",
      "with descendants 32\n",
      "no descendants 25\n",
      "\n",
      "is not Marshallese:\n",
      "total 52\n",
      "with descendants 14\n",
      "no descendants 38\n"
     ]
    }
   ],
   "source": [
    "# print out the number of different categories that have descendants and don't have descendants\n",
    "print(\"\\nvaccination status unknown:\")\n",
    "print(\"total\",len(WA_df[WA_df['vaccination_status'] == \"unknown\"]))\n",
    "print(\"with descendants\", len(WA_df[(WA_df['vaccination_status'] == \"unknown\") & (WA_df['number_children'] > 0)]))\n",
    "print(\"no descendants\", len(WA_df[(WA_df['vaccination_status'] == \"unknown\") & (WA_df['number_children'] == 0)]))\n",
    "\n",
    "print(\"\\nvaccination status not up to date:\")\n",
    "print(\"total\",len(WA_df[WA_df['vaccination_status'] == \"not_up_to_date\"]))\n",
    "print(\"with descendants\", len(WA_df[(WA_df['vaccination_status'] == \"not_up_to_date\") & (WA_df['number_children'] > 0)]))\n",
    "print(\"no descendants\", len(WA_df[(WA_df['vaccination_status'] == \"not_up_to_date\") & (WA_df['number_children'] == 0)]))\n",
    "\n",
    "print(\"\\nvaccination status up to date:\")\n",
    "print(\"total\",len(WA_df[WA_df['vaccination_status'] == \"up_to_date\"]))\n",
    "print(\"with descendants\", len(WA_df[(WA_df['vaccination_status'] == \"up_to_date\") & (WA_df['number_children'] > 0)]))\n",
    "print(\"no descendants\", len(WA_df[(WA_df['vaccination_status'] == \"up_to_date\") & (WA_df['number_children'] == 0)]))\n",
    "\n",
    "print(\"\\nage at least 20:\")\n",
    "print(\"total\",len(WA_df[WA_df['age'] >= 20]))\n",
    "print(\"with descendants\", len(WA_df[(WA_df['age'] >= 20) & (WA_df['number_children'] > 0)]))\n",
    "print(\"no descendants\", len(WA_df[(WA_df['age'] >= 20) & (WA_df['number_children'] == 0)]))\n",
    "\n",
    "print(\"\\nage under 20:\")\n",
    "print(\"total\",len(WA_df[WA_df['age'] < 20]))\n",
    "print(\"with descendants\", len(WA_df[(WA_df['age'] < 20) & (WA_df['number_children'] > 0)]))\n",
    "print(\"no descendants\", len(WA_df[(WA_df['age'] < 20) & (WA_df['number_children'] == 0)]))\n",
    "\n",
    "print(\"\\nis Marshallese:\")\n",
    "print(\"total\",len(WA_df[WA_df['community_status'] == \"Marshallese\"]))\n",
    "print(\"with descendants\", len(WA_df[(WA_df['community_status'] == \"Marshallese\") & (WA_df['number_children'] > 0)]))\n",
    "print(\"no descendants\", len(WA_df[(WA_df['community_status'] == \"Marshallese\") & (WA_df['number_children'] == 0)]))\n",
    "\n",
    "print(\"\\nis not Marshallese:\")\n",
    "print(\"total\",len(WA_df[WA_df['community_status'] == \"Not_Marshallese\"]))\n",
    "print(\"with descendants\", len(WA_df[(WA_df['community_status'] == \"Not_Marshallese\") & (WA_df['number_children'] > 0)]))\n",
    "print(\"no descendants\", len(WA_df[(WA_df['community_status'] == \"Not_Marshallese\") & (WA_df['number_children'] == 0)]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 208,
   "metadata": {},
   "outputs": [],
   "source": [
    "reg_df2 = WA_df.copy()\n",
    "reg_df2['marshallese'] = np.where(reg_df2[\"community_status\"] == \"Not_Marshallese\", 0, 1)\n",
    "reg_df2['vaccination_status_unknown'] = np.where(reg_df2[\"vaccination_status\"] == \"unknown\", 1, 0)\n",
    "reg_df2['up_to_date'] = np.where(reg_df2[\"vaccination_status\"] == \"up_to_date\", 0, 1)\n",
    "reg_df2['has_descendants'] = np.where(reg_df2[\"number_children\"] == 0, 0, 1)\n",
    "#reg_df2.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 210,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "Call:\n",
       "glm(formula = has_descendants ~ vaccination_status_unknown + \n",
       "    up_to_date + age_bin + marshallese, family = binomial(link = \"logit\"), \n",
       "    data = reg_df2, weights = na.action(na.omit))\n",
       "\n",
       "Deviance Residuals: \n",
       "    Min       1Q   Median       3Q      Max  \n",
       "-1.3953  -0.8928  -0.7472   0.9904   1.6807  \n",
       "\n",
       "Coefficients:\n",
       "                           Estimate Std. Error z value Pr(>|z|)   \n",
       "(Intercept)                 -0.7140     0.3918  -1.823  0.06838 . \n",
       "vaccination_status_unknown   0.7151     0.7709   0.928  0.35362   \n",
       "up_to_date                  -0.7569     0.6886  -1.099  0.27165   \n",
       "age_bin                     -0.3774     0.5097  -0.740  0.45906   \n",
       "marshallese                  1.2129     0.4247   2.856  0.00429 **\n",
       "---\n",
       "Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
       "\n",
       "(Dispersion parameter for binomial family taken to be 1)\n",
       "\n",
       "    Null deviance: 148.44  on 108  degrees of freedom\n",
       "Residual deviance: 136.48  on 104  degrees of freedom\n",
       "AIC: 146.48\n",
       "\n",
       "Number of Fisher Scoring iterations: 4\n",
       "\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "               (Intercept) vaccination_status_unknown \n",
       "                 0.4896912                  2.0444126 \n",
       "                up_to_date                    age_bin \n",
       "                 0.4691056                  0.6856429 \n",
       "               marshallese \n",
       "                 3.3633113 \n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/lmoncla/anaconda/envs/LHM-basics/lib/python3.7/site-packages/rpy2/rinterface/__init__.py:146: RRuntimeWarning: Waiting for profiling to be done...\n",
      "\n",
      "  warnings.warn(x, RRuntimeWarning)\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "                               2.5 %    97.5 %\n",
       "(Intercept)                0.2202503  1.036226\n",
       "vaccination_status_unknown 0.4715406 10.152389\n",
       "up_to_date                 0.1101775  1.733805\n",
       "age_bin                    0.2472645  1.859645\n",
       "marshallese                1.4855957  7.911204\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "Analysis of Deviance Table\n",
       "\n",
       "Model: binomial, link: logit\n",
       "\n",
       "Response: has_descendants\n",
       "\n",
       "Terms added sequentially (first to last)\n",
       "\n",
       "\n",
       "                           Df Deviance Resid. Df Resid. Dev Pr(>Chi)   \n",
       "NULL                                         108     148.44            \n",
       "vaccination_status_unknown  1   0.4135       107     148.03 0.520200   \n",
       "up_to_date                  1   1.1693       106     146.86 0.279535   \n",
       "age_bin                     1   1.8120       105     145.05 0.178268   \n",
       "marshallese                 1   8.5648       104     136.48 0.003427 **\n",
       "---\n",
       "Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# now run with age as a bin instead of a continuous variable\n",
    "%R -i reg_df2\n",
    "%R model.has_descendants = glm(has_descendants~vaccination_status_unknown+up_to_date+age_bin+marshallese,data=reg_df2,family = binomial(link=\"logit\"),na.action(na.omit))\n",
    "%R print(summary(model.has_descendants))  # print the summary\n",
    "%R print(exp(coef(model.has_descendants)))  # exponentiate the coefficients\n",
    "%R print(exp(confint(model.has_descendants)))   # exponentiate the confidence intervals\n",
    "%R print(anova(model.has_descendants, test=\"Chisq\"))  # run a chi square?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Interpretation \n",
    "\n",
    "An odds ratio of 3 for Marshallese means that the odds of having descendants and being Marshallese is 3. So Marshallese people are 300% more likely to have descendants than non-Marshallese, holding all other variables constant (assuming you are under 20 years of age and have a known, vaccinated vaccination status). \n",
    "\n",
    "An odds ratio of 0.46 for up-to-date means that the odds of having descendants and being not up-to-date is 0.46, holding all other variables constant (being not Marshallese, having age 0, and having a known vaccination status). This is not statistically significant. \n",
    "\n",
    "An odds ratio of 2 for vaccination unknown means that the odds of having descendants and having an unknown vaccination status is 2, holding all other variables constant (being vaccinated, age <20, not Marshallese). However, this is not significant. \n",
    "\n",
    "Finally, we get an odds ratio of 0.69 for age bin, meaning that individuals who are older than 20 are less likely to have descendants in the tree. \n",
    "\n",
    "Significance testing by the Wald statistic is: b/se(b), where b is the coefficient and se is the standard error. \n",
    "\n",
    "For the ANOVA, I tried adding Marshallese as a predictor at the very beginning or very end, and it continues to be significant. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "LHM-basics (python3)",
   "language": "python",
   "name": "lhm-basics"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}