1_optimum_diet_init.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Matrix Porridge\n",
    "by Kyle Archie, M.Eng\n",
    "\n",
    "## The Problem\n",
    "\n",
    "Many COVID-19 survivors report losing their sense of smell and / or taste. While for most this returns over time, others continue to struggle with this issue months later. This problem made me remember an old project I worked on over a decade ago, way back in 2010. I was an entrepreneur at the time with very little time to eat during the day. I thought how, since I was basically just shoveling food into my face to satiate hunger rather than enjoy nicely prepared meals, it would be nice if those meals were at least healthy, even if they weren't all that tasty. \n",
    "\n",
    "I immediately thought back to the movie, \"The Matrix\" and the famous [\"tasty wheat scene\"](https://www.youtube.com/watch?v=v1EcrD5IyxM). Basically, in this dystopian future where the atmosphere had been ruined and nothing grew above ground, food had to be synthesized. The result was a gelatinous goop... \"a single cell protein combined with synthetic aminos, vitamins, and minerals. Everything the body needs.\" (The Matrix, 1999)\n",
    "\n",
    "Now there are a million products out there that claim to be healthy meals on the go for you, but looking at the nutrition labels on these, I never really found anything that quite delivered what I wanted. All focused too much on making something \"edible.\" But those days, I was looking for something you'd down like a shot of bad whiskey. Something that got the job done, but didn't necessarily taste good. Not a bar filled with chocolate and peanut butter.\n",
    "\n",
    "So back to these COVID-19 survivors. Losing their sense of taste is clearly not something anyone would ask for. But perhaps there's some lemonade they can make from these lemons. \n",
    "\n",
    "This project aims to find the optimal list of ingredients necessary to create the perfect meal for people who can't taste. It doesn't matter if the ingredients make sense to a chef, just as long as they solve the equation and meet the dietary constraints of a typical adult.\n",
    "\n",
    "## The Data\n",
    "For now, I'm going to use my original dataset, which is over a decade old. The reason for this is that the calories in a carrot shouldn't change over time, and the data, which came from the USDA, comes in several dozen separate files, which originally took me hours to stitch together via various joins and SQL queries to create what I wanted. Unfortunately, I no long have the original database I used; only the final export. I intend to revisit this eventually, but for now, it should let us take a first pass at this. \n",
    "\n",
    "The latest version of that dataset can be found [here](https://fdc.nal.usda.gov/download-datasets.html). At the time, I put them together using Microsoft Access, which is ideal for these one off sorts of joins, even today. Anyway, I have in the \"Matrix Porridge (2010).xlsx\" file a couple tabs of interest. First is the core ingredients list, which contains everything from protein content to the vitamins and minerals in 2331 different core ingredients. Second, is the dietary requirements tab. This was based on the needs of a male adult aged 18-30, because that was my demographic at the time. But that can be easily modified to any group, and the variation is mostly around quantity, so it should scale.\n",
    "\n",
    "The requirements tab is a bit more recent, since updating this was not very difficult. The source is [NIH Dietary Reference Intakes site](https://ods.od.nih.gov/HealthInformation/Dietary_Reference_Intakes.aspx). The reason for using this and not the FDA Daily Recommended Values was that the FDA's numbers are fixed numbers instead of ranges. Obviously some things we aim to maximize and others we hope to minimize, but while we may hope to reduce our sodium intake, for instance, we know that we need at least SOME salt in our diet, or we die. Interestingly enough, the NIH range data shows that the minimum recommended intake for even polyunsaturated fatty acids is not 0. It is 0.6% of our daily intake of calories. Anyway, I wanted to use ranges where applicable here, so I went with the NIH data. \n",
    "\n",
    "With that said, there are a lot of vitamins and minerals that can be found in daily supplements, so for the sake of not making things overly complex, I will only be considering core nutritional requirements as constraints here with the exception of sodium (because my early attempts at this had far too much salt in them), and I will be doing so based on an assumed 666.7 calories per meal (i.e. a 2000 calorie daily diet). Based on the NIH data and the assumption that there are [9 calories per gram of fat, and 4 calories for each gram of protein or carbohydrates](https://www.nal.usda.gov/fnic/how-many-calories-are-one-gram-fat-carbohydrate-or-protein), we can calculate the acceptable range in grams for Fat, Protein, Carbohydrates, which we'll need for an apples to apples comparison, since the ingredients data provides all this in grams.\n",
    "\n",
    "Also, the carbohydrate numbers from the [NIH ranges for macronutrients (fat, carbs, etc.)](https://www.ncbi.nlm.nih.gov/books/NBK56068/table/summarytables.t5/?report=objectonly) don't line up with the values in the [recommended daily allowances for macronutrients](https://www.ncbi.nlm.nih.gov/books/NBK56068/table/summarytables.t4/?report=objectonly) once we convert from calories to grams. After a bit of investigation, the reason appeas to be that these ranges don't break fiber out from the more general carbohydrate category. However, this is where things get a bit complicated. [The calorie content of fiber depends on whether it is dietary or insoluble](https://www.fiberfacts.org/fibers-count-calories-carbohydrates/). Dietary fiber are essentially complex carbohydrates, and so provide the same 4 calories per gram, but insoluble fiber goes through the system without being digested, and so provides no caloric value. A number of sources therefore estimate average fiber value to be about 2 calories per gram, which may be a bit off, but for now, that's what we'll use.\n",
    "\n",
    "One final issue is that there's a lot of disagreement amongst nutritionists regarding carbohydrates and whether these should be minimized or not. Also, there's been a fair amount of recent research surrounding the microbiome that has found that increasing fiber intake is the easiest way [to increase biodiversity in the gut biome](https://www.medicalnewstoday.com/articles/short-term-increase-in-fiber-alters-gut-microbiome#:~:text=Fiber%20promotes%20a%20healthy%20gut,significantly%20altered%20the%20gut%20microbiome.). Therefore, I will treat the recommended fiber intake as a minimum threshold with each gram accounting for 2 calories and have my fuction use the overall carb max threshold as the overal max for all carbs via a separate conditional statement."
   ]
  },
  {
   "attachments": {
    "bec3c86a-4696-47e8-912b-4a920ed1d068.png": {
     "image/png": "iVBORw0KGgoAAAANSUhEUgAABBUAAAESCAYAAACivO0YAAAgAElEQVR4Aey9TYskWZamqf8oVsNsh2FWQTskNthP6KxFTNOQzKIIfFEYSm9ikmKgAoaM2KoxuQmoLOiFZ9KbqjKne5Fk5SIgusnF4AvLJDf+D2R4Ve0xe+3YFRXRLzNR0dfAXUSu3HvuOc85cuXeo6Kiiy5/IRACIRACIRACIRACIRACIRACIRACIbAHgcUebdIkBEIgBEIgBEIgBEIgBEIgBEIgBEIgBLokFRIEIRACIRACIRACIRACIRACIRACIRACexFIUmEvbGkUAiEQAiEQAiEQAiEQAiEQAiEQAiHQTCr86U9/6vIvDBIDiYHEQGIgMZAYSAxMKwYydQ2BEAiBEAiBqRHoTSp8/vy5y78wSAwkBhIDmxjQwioscj0kBhIDbxkDGofyFwIhEAIhEAJTI5CkQpInWSglBhIDI2IgSYUsJt9yMZm+E3+KgSQVpjaNjj4hEAIhEAIikKTCiMVEJnOZzCUGEgNJKiQGMg4kBt46BpJUyOQ9BEIgBEJgigSSVEhSIZ9SJwYSAyNiIEmFLCjfekGZ/hODSSpMcSodnUIgBEIgBJJUGLGYyEQuE7nEQGIgSYXEQMaBxMBbx0CSCpm4h0AIhEAITJFAkgpJKuRT6sRAYmBEDCSpkAXlWy8o039iMEmFKU6lo1MIhEAIhECSCiMWE5nIZSKXGEgMJKmQGMg4kBh46xhIUiET9xAIgRAIgSkSSFIhSYV8Sp0YSAyMiIEkFbKgfOsFZfpPDCapMMWpdHQKgRAIgRBIUmHEYiITuUzkEgOJgZpUWHz9oVssFr3/Pny9OJtkzbvvf+rOSd9TXo8/ff7cLX723aNf58Dlw+fP3bvF++7zj9+dTUye0sfnLDtJhUzcQyAEQiAEpkggSYUkFTLJTAwkBkbEQE0qaGHy3Y+fu8XiXfd+oQTChzVHLeAWi/frRel3P1P5T5Pmq+TIOeg5ZiEof3z+7WELZ/F4t1h0H35UYuH9rNgoVj//9v2k43GMny+5TpIKU5xKR6cQCIEQCIEkFUYsJi55AhPb8wl9YmATA2OTCuKlT/71FIMWp1P+dFh6ekLknH3NEwaHLJqRMXW/7eun97/VEwvTjsl9bbuUdkkqZOIeAiEQAiEwRQJJKiSpkE+tEgOJgRExsE9SYcoLdp6y+On7d7PwP19HOSSpwFMmc114Y99cnky5lESC25mkwhSn0tEpBEIgBEIgSYURiwm/oWc/n9wnBi4zBsYmFVi46UmFusDlCQbexcDijk/Iebrh/fft7/R7PWSwpS99Gk0Z8lsxiy60q7Lff735Codk1eSI2+gL8Nr3+99untigvbeTXNfPz/XK/Hrx+L4Dr0NCAbvpb5vdqusyvH/kwKbKqXbqKy7OT3L7fChZQ+3V/0+/ff9oK4kffIZ+bP/2P+jdHu8e/a73QLg+/l4Ivt4x5SdoKu8cP425SSpk4h4CIRACITBFAkkqJKkwi08pM+l8mnSGxWlYbEsqsLhj21qM8mSAFtI/ff7pxYKR85KhRSDHvkBm8axFJotGXxizWPVFpS8oiY1WW53zhTXtkEkCgDroVb9C4QtfyeD8h896seX7xwQFtojVkEzqYivHcPjwo5IXmyRIiz12oxt1JAeZbr+X0ZYtPFqM8dk2H25rDwe1lw6810F20g4/YAsJB/rGL9JXfv7w/fN3TMAOBtiV7WnGjWNzTVJhilPp6BQCIRACIZCkQpIKSSokBhIDI2JgW1JBC7mftGh++NUAX9j1LSpY3NVFIQtaFokc+4JTC0ISA1qAIgOZLBh13NKFtshGx9qHyimjLovbmnSgTxa7LH6RXdtRrm09xzEysYs+6zE6igVtXL72qYMdKkNX5FJnmxz6ph9nXH1Wj9Xntvbb+qcduqK7c6aMeHj//U8v3ulR61ROOZ52ciFJhUzcQyAEQiAEpkggSYURi4lMsqY9yYp/4p/XiIGhpIJ+/YFFpBalvthz/Xzh6AkB2rLorceSwcJSi0YSAyyA/Vhl/EOe60Ddes51Y9FMGXVZlCKfLQtZzlf7KWdR7PpwDllskYndtK3H6Kh26O3ytd/iSb/oOiQHbujHFja1j3o81H5b/yRa0BXdYSIbaU+Sqz6loDq0g23llONpj6dJKkxxKh2dQiAEQiAEklRIUiGfUicGEgMjYmBMUkELMhZ/WnD6wo0FJU8OsDCmTl2A1mMWe7STfBazOof82i/t6lZyvL3OsyiVDBbnlFGXRSmL2yq37zzl2O/tONcnE5tZQNdjdHS9Xb7bhh0qo1/kDskZYlx9Vo+H2g/177HV52dsev+zdy+eUpDNsMO/lVOOk1TI1DgEQiAEQiAEdiWQpMKIxUQmWdOeZMU/8c9rxMC2pIIvVKULCzcvZ0GohfMHvVPh4R0AuyQVJIMFcMvmxwWlkgKfP6yTRe++/tBcXFLXF5e+qEUv11svJGSh7At41aE+cmuCoE+22g3JhCe21+Nni/Xfvu+++1rvEfjpRbIM3bBZctxHriN1KmdkeHIExtiBzHosWdvab+tf594v3j/6terFMX1W/pyvNlOe7XmMo3lSYddpbuqHQAiEQAi8BoEkFZJUeDHxzuTyPCaX8dPr+qkmFVjYanHNPxbWvjjUOS2GWezpWAvSdw+/8KAF6E8/fvf4Pob1+Z/pVw6efn2BBWKrT+QTDyxaVa5/fYtj9GGhrvau97ufPf2aAP3TB23pAxkkICiHx1A7nR8r03+Vwu3Dbhb09Fm3ztDruu3o38eOvryeJzZU/r7Hh9JnVHtLDKlNZUvfntzAVtnY0h0bqz9pl+3rjin78E5S4TWmxukjBEIgBEJgVwJJKiSpkKRCYiAxMCIGalJhnwXBoW3qwpuFpS+Od+lDi1tflLLolNzWonQX2al73AVqTVrge21fJG/0gsaHJ1XcD0pM7BsrLif7x/XtLjyTVNh1mpv6IRACIRACr0EgSYURi4ldbvip+3aTrbAP+1PGwFsnFVhU1sW+EgHv9Fj8j89/OnAsC32qzSfXSSpM9xpSAognQty38t+/fv9u/aSLEkTf/fih++l7fU3iuS0kpGr81Ho5fs5tajySVHiNqXH6CIEQCIEQ2JVAkgpl4jW1CUT0mfYEL/65HP+8dVJBsdZ6BN6fNNg3HrVg/X+/1lcuvnv8KkeeVphWbHvCx59SIEnA1zpIEHksHJp4clnZf9u4SFJh12lu6odACIRACLwGgSQVklR48YlWJo1vO2kM/2nyn0JSIbExzdiIX+KX14qBJBVeY2qcPkIgBEIgBHYlkKRCkgpJKiQGEgMjYiBJhSwcX2vhmH4Sa30xkKTCrtPc1A+BEAiBEHgNAkkqjFhM9N3cU56JX2LgcmIgSYXL8XWu6/h6qjGQpMJrTI3TRwiEQAiEwK4EklRIUiGfUicGEgMjYiBJhSw0p7rQjF6XE5tJKuw6zU39EAiBEAiB1yCQpMKIxUQmbJczYYuv4+u+GEhSIbHRFxspT2y8VgwkqfAaU+P0EQIhEAIhsCuBJBWSVMin1ImBxMCIGEhSIQvH11o4pp/EWl8MJKmw6zQ39UMgBEIgBF6DQJIKIxYTfTf3lGfilxi4nBjQZD7/wiAxkBh46xh4jclh+giBEAiBEAiBXQgkqZCkQj6lTgwkBkbEgBYSSSJdThIpvo6vpxgDGofyFwIhEAIhEAJTI5CkwojFxBQnFtEpE97EwOvGQJIKr8s78R3eiYGXMZCkwtSm0dEnBEIgBEJABJJUSFIhn74mBhIDI2IgSYWXC5ws+sIkMfC6MZCkQibvIRACIRACUySQpMKIxUQmTa87aQrv8J5iDCSpkLicYlxGp8uKyyQVpjiVjk4hEAIhEAJJKiSpkE+pEwOJgRExkKTCZS3esliPv6cYA0kqZOIeAiEQAiEwRQJJKoxYTExxYhGdMuFNDLxuDCSp8Lq8E9/hnRh4GQNJKkxxKh2dQiAEQiAEklRIUiGfUicGEgMjYkCT+fwLg8RAYuCtYyBT1xAIgRAIgRCYGoEkFUYsJvJpyctPS8IkTC4tBrSQuDSb52ZvfDj/cWvuPpZ9+QuBEAiBEAiBqRFIUiFJhSyUEgOJgRExMPfFytwSCC174sMkFVpxcU5lSSpMbRodfUIgBEIgBEQgSYURi4lzmnBE1/lPmuPjt/FxFqRvw/2Y8R4fnr8Ph+Jh7j5OUiGT9xAIgRAIgSkSSFIhSYV8Sp0YSAyMiIG5L1aGFmtzOB8fJqlw7nGcpMIUp9LRKQRCIARCIEmFEYuJc5+ERP/5T6Tj49P7OAvS0zM+dRzHh+fvw6EYmbuPk1TIxD0EQiAEQmCKBJJUSFIhn1InBhIDI2Jg7ouVocXaHM7Hh0kqnHscJ6kwxal0dAqBEAiBEEhSYcRi4twnIdF//hPp+Pj0Ps6C9PSMTx3H8eH5+3AoRubu4yQVMnEPgRAIgRCYIoEkFZJUyKfUiYHEwIgYmPtiZWixNofz8WGSCucex0kqTHEqHZ1CIARCIARGJRXeff9Tt1gsHv9997NF9/nzTy8WIouvPzzWUf0PX6ve/Ccx1cbvfvzcLRbvLtb+yiPHl3cNzNHn2xakHz7rmn//bPz76ft33U+fP3ffff3derz0cVTnWoyQwxg71IaxRuPtu8Wi+/yj+nq7eJO9i59998hhX538XnLM+0h8uHts4NO+mK3xRn3mDK0YIM631alyxx5v87FkoN8Ye6iLnmPaqA+3772uy88fjnZdJqmQiXsIhEAIhMAUCYxKKtSb5NAkoXV+7IRgqvU0ef/823ETdib6x5wMi8suOkyVY/TafVIfZtNg1rdYYeFfxz3KSRD4AqVvccJiZJc2jDe1/7eMG+zYRydPKLCY+/zb90dZlMWHu19L+KMvZvvirC8GajnXxbEW330+Rk+ul7H9oZ9icRcG73+rROOiq/0ceh9PUmGKU+noFAIhEAIhsFNS4bvv9QnU5tO4enP97vufunfr89P4xIwJxDG2TCqONbHdR6cp6LCP3mmz+yQ+zKbJrLVYYeHQt3jW+V0SBNX3XPfbFjR1kVZlvMXxvjppwfXh681THMiQ7cdK0MaHu11bxPe2+OuLL/xXr43Wol4JuFqvT+5QecvH3oYkiWwae0+nTZ33uNy6DztPKnA9j+23ytRxkgqZuIdACIRACEyRwE5JhZ+UNHj4ikO9UX74/ruOT+aONTlo3VDfogybD5kIHKr3FHQ41Ia0321CH17T4lUXKyyati16tYjg6w8sKHZZoI1pgx5TGnePpRPj3i6LuW3XTXw4/pqSD9/pQ4SfvVt/4r6rD/pigKSCXzfys88ptvlw6Fz1sdeXTu+/1gcjmw9Hxiar9onDVlIBOYfMJZJUmOJUOjqFQAiEQAjsnFTgRulZfh7n60sqMLlQG/7Vmypt/by304T5w49PT0qovU+41Y5PBDWJcD2/+3rx+D1fn3h7n5oweRsmG0wC0EsTnw+f9e6IzaSk6vXdrco3kzDXp9rCd59dh+++V7uNXJ9gtXQ45nc0fdKV/fGT7rC6LFZ1seKLozELLh+vPuhaf3j3AOOEjxFc/97G+2iNG4xtLqeOT61xk7HY26l/njzT2Od9j4l7ZKHTmDa1DuMxY3E9/y+fP3dffLnqvvjym/W/1X/8pvv8+U/dF3//X58de7v4cPw1q/vOT79VUmHzjoxjxsDYe9oxfOz+/+63n9fvHeH64TrzOtr3a1vX5+LrJwbEtq4Ltf9J84EHRh7vxC999NmMLswxhjgnqZCJewiEQAiEwBQJ7JxU8Bsqkz199UGLXG6OfmPVDZqbqSav9Uar87TTzZSbORNtjrmBf/ejFvTvHic73LDRi3beL/qgh9+0axnHsk0T+g8/6iWVm4U+k+866ZAO6KU62IMu6IauOs++LxpU5hMUydqmg/TIvzBIDLxODNQFKWPZ2EW3X+sak35SkvRhMcLYgszW+MC4RR3GX8Ybxrlt4xOLROpKFvtqh2xs4pixbGysMea57LFtVQ+bWGjBp8pg0fk3X37Tff4f//lxPFz9+r+uEwy1fnz48gXLlZGO5XfugcQo8deq3yrbFgN+LZzax+imPnlqyOcV1S705hps1aWs1vF459rhvt6aS9AXOiju2Ufvuk1SYYpT6egUAiEQAiGwc1JBNzgmfI83y+83LzCk3G+s9YZYb8ZMLvraUL+e54bNZBc5mqAwAfUEgfSox17GjbzW4abvctWuTy+dgwO6oSuLAI5JGvikze1Apz4d1Ff+hUFi4HVi4JgLUl3brWudsYHFSq1Tj+V7xgcfI/vGJ8oZmzhmrKn91+OxsdbSaWxb6iFDYy88OOfb//TPemLhm+5f/l5PKnzudPz5n1fNsTE+HP4VAsXEdz/Tey1+asaos9+2j/88Lr3+h99K/k+PibW+evh0Xx/Tp+z66fvNCz/9OuJaoF69f6uceQHXCdcNccmx21CvHXj4XML1QDZ69G2TVMjEPQRCIARCYIoE9koqcANdT/b0Uq2HX0XgZuw3Vr8xct4nidxo+9rQVz2PLCYErZszEwEW8/W4NVmoddDPJwJq16eXzlXdOJYM/9daWLTs6NPB2Wb/dRaW4Xy5nOuC1K9Lxpht8VGv7XqstnUhUuv4MclT9PAxsm98Qr6PQ9pHf86zWKrH2+zzc9t0qn1z7PojqyWHc2z/pK9B/Mf/3PG0Qt9TCqofHw4nFfruV/LT//WNvqK3+XoffvOt+7DPd5SziOZYcohpfMv2EB8jo88u11l1qcc1oTLmBejM9cV1wrHLqtdOn51+TXt79K7bJBWmOJWOTiEQAiEQAnslFfwmyE3Vb8b1xsjNVTdkbr608xstN2y/iVK/T2YrqcDEhIkAk4N6rH4oo2+OaeP6IVft+vRyDujGJIVjt89ZSod6rLp9Oric7F/uYje+fx3f1wWpuDNe1PEJn+ja1QtudVyv7XqsOoyVjI+1Tj1WG8YH16FvfKry0ZNtPV+PqTe0bek01KbvvMZPePTV4WmFv9F7Ff756WsQtX58eFhSgftk5do67osB7ofcTz2mKWvJO8TH0uXD189/lhT9lMzgfq9+0c/v11zn2M/1RVxy7NdgvXa8v5adyEBmi4HKklTIxD0EQiAEQmCKBHZKKvhNmRtm62bsN9Y6YeCG7TfOWqY273+mifiH3sU7cpHDDbs1EUBHJgYc6wbtZcj0SYaX6fvI9TuZbiuTAOxBFyYL/kmM+LWSCM/6+37z02rPykwH+sv2dRaV4XzZnFsLUr82GYuIE1337xe6hjcLOa/buvbVjnEVWbWN6tTxhTY+FjHmeJnaMk5qLPIFEj/jiCz6r8d9crGZ7dh61GfLeIxu2N9ahNFGW96t8J/0boXPep9CO1bjw/57aosZ/D1eWvVaZX0xQExxH6YeMdeSpbJDfKy4asUQ8ebXCfpQ1nfN6IkN6mATx9KXMux6xvLhPv7hsxJmSnZsxghd29Tv45CkwhSn0tEpBEIgBEJgVFKBSawmFtzwuNHqRs2+zvs/Jg3cuHXOf+GA87p5eh1uzFUui3RuttvOc0NHn83PSD3pxwSDCYTqqV9+PkvH1MH+MXq5PMnAxr5yZNP/u/XPXW30pD/ZSz0vg0O27QVEuITLMWOgtSBFfh1vdD3X8YpruO9a91+3of277//12ZjKmOTj5Xv9uo39aoz/Og1y9P14dO0bN2u5y5UcjWU+jjG2IVdbXzipDf8qC2/j+5XjLuPd3/z6T92ffv3/PNrpctmPD4d9CCtt3Z8kevx8a9/b4H9tPQaqn5lXtOR52a4+rrp4zPo1hJ5cX66fdFvYfZk63v6dfiHCfvVJyQKuSbedMYC41jXnP3Gpush3u30/SYVM3EMgBEIgBKZIYFRSwW9o2c9CLTGQGLjEGNi2IL00Hlp0jV1kvhablV7QaL8A0eo3Pnwau6bow5bPvCw+ztcfpjiRjk4hEAIhEAJdl6RCz2OyPpHJ/tNENCzC4lJjIAvSTezraYV3+lrHj5t3RbxlPHzx9/91/YLGf/kf/7n7l79/+sUH/VKA/v1/n/787MmF+HB6PhyKn/j4+T1HMZy/EAiBEAiBEJgagSQVklR4NukemuDl/PMJXnhcDo8sSDfvZHin74BPIKGga0+Pwyt5MPQuBa7T+HB6PsQ3fdv4+PkYm6TC1KbR0ScEQiAEQkAEklRIUiFJhcRAYmBEDGRB+nxx07cInHJ5fHj+PhyKr7n7OEmFTN5DIARCIASmSCBJhRGLiaFJTM7Pf6IaH8fHc1+sXEKMx4fzv47n7uMkFaY4lY5OIRACIRACSSokqZBPqRMDiYERMTD3xUqSCvNfcMfH5+/jJBUycQ+BEAiBEJgigSQVRiwmLmEiFhvPf7IZH57Wh0kqnJbva8RvfHj+PhyKk7n7OEmFKU6lo1MIhEAIhECSCkkq5FPqxEBiYEQMzH2xMrRYm8P5+DBJhXOP4yQVMnEPgRAIgRCYIoEkFUYsJs59EhL95z+Rjo9P7+MsSE/P+NRxHB+evw+HYmTuPk5SYYpT6egUAiEQAiHQm1TQjSv/wiAxkBhIDCQGEgOJgcTAdGIgU9cQCIEQCIEQmBqB3qTC1BSNPiEQAiHwlgS0qMrfeROID8/bf2O0n7uP527fGB+nTgiEQAiEwPQIJKkwPZ9EoxAIgQkSyGR+gk7ZUaX4cEdgZ1h97j6eu31nGHJROQRCIARCoOu6JBUSBiEQAiEwgkAm8yMgTbxKfDhxBx1Bvbn7eO72HSEEIiIEQiAEQuANCCSp8AbQ02UIhMD5Echk/vx8VjWODyuR+R3P3cdzt29+ERmLQiAEQuAyCCSpcBl+jpUhEAIHEshk/kCAE2geH07ACSdWYe4+nrt9Jw6PiA+BEAiBEDgRgSQVTgQ2YkMgBOZFIJP58/dnfHj+PhyyYO4+nrt9Q/7N+RAIgRAIgWkSSFJhmn6JViEQAhMjkMn8xByyhzrx4R7QzqzJ3H08d/vOLNyibgiEQAiEwAOBJBUSCiEQAiEwgkAm8yMgTbxKfDhxBx1Bvbn7eO72HSEEIiIEQiAEQuANCAwmFe71ExHXq26xWDT/XS0WXfdp9Qaqp8sQCIEQeD0Crcn8ncbHxfJxbOw+Ll8odHV7v/X8iwY9BcuPXbe61pCtUfnp71jynyTOd29fH64+yc9Xj368u2neOpvgxsTIIfKbnV5w4b4+PuQ6ek0ft+y7YHfH9BAIgRAIgYkQaM6MWjctTWiXSiB0un0+/S1v75NUeMKRvRAIgZkSqOOiRsLlQkmEzZioRUkdI33cZOHRSjwMIaNtTSocS/5Q/3M5X32o9MzyWknxJx8qgV59tLzR+U2dXVioxVCMSN6+8nfR5VLq7uPjQ66j1/Zxte9S/Bo7QyAEQiAEpk3g4KTCWPPWU7Kbl5/ijW2feiEQAiHwlgTqZP7+k7R5empAY9yVkgwPT27xlJcvULV4qYmBMTZp0Xl1u3rW9pjyx+gwhzrVh6uP3aO/8KaezPMnEfQUwf3tfveuoRhRn4fIn4NPjm3Drj4+9Dp6bR9X+47NL/JCIARCIARCYB8CeycVdk0SLG7unk3U9lE2bUIgBELgrQgMTea1OLy7uXpUT8dXi6tni9ZW2WODnh0lIpSY0JMQnpBoyWqV9Yi9yOIhHwpKvVfpmK//eYJoH4Dyj8cI/R1L/j46za3Nrj5uXTOtsrGcTu3jMfaN1TX1QiAEQiAEQuBYBPZOKrRunP6dRE2S7m+v1g+M+neOW4+WHsuYyAmBEAiBUxHYNpnXeLi6VkLh6ckFf6QanZSM1Xg4dnHqyduaVDiGfPS6lO02H4pB61Nr2Ii37l/+FAPnxmxbMeLtDpXvsi55f1cfH/M6eg0fD9l3yb6P7SEQAiEQAm9HYKekAp+msPXJlW6m/qmcf7+YiZrXfzuT03MIhEAI7E6gbzLvyVR/p8IxFiv+XfskFXb3WW3R50PqDS0Kd00KIbcvRjjPdl/5tM+263b18TGuU3F/LR8P2ZcYCIEQCIEQCIG3ILBTUsEnzJp81cc4MYDHRamfpAJksg2BEDhXAkOTeRYVekJLf63FisZN/YrA6h/1ixFPvyZAopatflVn9d9WXffx6Zd1xiQVkD/2SYhz9cW+eg/6UC9kHPg1I/lBCXJY47O6bf0yUo2Rlh3Ib51L2TCBXX186HVa4+XUPh6yb5hQaoRACIRACITA8QnsnVRQouDu9mnCK9W4mWpC6zfqJBWO77hIDIEQeF0CQ5N5xjmSClp0+tNb0rZV1mcFydm6WNWxxtiWrFZZn/xLLN/mQ92/xiRjdG/b96m7GiMtHxwivyXv0sp29XHrmmmVjeV4ah9vs2+sjqkXAiEQAiEQAscmsHdSoSriSQSd82NusvtOxGpfOQ6BEAiB1yYwZjKvRAALU8Y9jqWvFq7+ssVdbKhtjy1/F13OtW6fD3W/Ihk0ZNvY5EOfHI+RVp1D5bdkXlLZrj4+xXV0Sh/32XdJPo6tIRACIRAC0yOwU1Kh9TgnJmlS5ud1U+XrD6qjYx4Z9Ud6aZ9tCIRACEyZwNBkXp9utl7WSBLh0O/L16SCWGncPZb8KbM/lm4tH4phTXgv9dRC42sQLR/vottQ+6Hzu/R1qXX38fExr6MhHw6dH/Jby76hNjkfAiEQAiEQAqcmMJhUIIvvj+C2PtFhwky95Y2+M7x4TDRoQqxjJsCnNizyQyAEQuCYBOpkXgsRxrttYxtjn+r4Uwu76tZKKkjGseTvqs851q8+dHYtX2oB6O++2PX+NRQjh8o/Rx+cWuddfYw+Hgu7XKev7eNqH/pnGwIhEAIhEAJvSWAwqfCWyqXvEAiBEJgKgUzmp+KJ/fWID/dndy4t5+7judt3LnEWPUMgBEIgBJ4TSFLhOY8chUAIhECTQCbzTSxnVRgfnpW79lJ27j6eu317OT2NQiAEQiAE3pxAkgpv7oIoEAIhcJgWP2EAACAASURBVA4EMpk/By9t1zE+3M5nDmfn7uO52zeHGIwNIRACIXCJBJJUuESvx+YQCIGdCWQyvzOyyTWIDyfnkqMrNHcfz92+owdEBIZACIRACLwKgSQVXgVzOgmBEDh3ApnMn7sHuy4+PH8fDlkwdx/P3b4h/+Z8CIRACITANAkkqTBNv0SrEAiBiRHIZH5iDtlDnfhwD2hn1mTuPp67fWcWblE3BEIgBELggUCSCgmFEAiBEBhBIJP5EZAmXiU+nLiDjqDe3H08d/uOEAIREQIhEAIh8AYEklR4A+jpMgRC4PwIZDJ/fj6rGseHlcj8jufu47nbN7+IjEUhEAIhcBkEklS4DD/HyhAIgQMJZDJ/IMAJNI8PJ+CEE6swdx/P3b4Th0fEh0AIhEAInIhAb1JBN678C4PEQGIgMZAYSAwkBhID04mBE80HIzYEQiAEQiAE9ibQm1TYW2IahkAIhMAMCWhRlb/zJhAfnrf/xmg/dx/P3b4xPk6dEAiBEAiB6RFIUmF6PolGIRACEySQyfwEnbKjSvHhjsDOsPrcfTx3+84w5KJyCIRACIRA13VJKiQMQiAEQmAEgUzmR0CaeJX4cOIOOoJ6c/fx3O07QghERAiEQAiEwBsQSFLhDaCnyxAIgfMjkMn8+fmsahwfViLzO567j+du3/wiMhaFQAiEwGUQSFLhMvwcK0MgBA4kkMn8gQAn0Dw+nIATTqzC3H08d/tOHB4RHwIhEAIhcCICSSqcCGzEhkAIzItAJvPn78/48Px9OGTB3H08d/uG/JvzIRACIRAC0ySQpMI0/RKtzoDAX7uu++IXv+u++PKb9b8/fPvNGWgdFfclkMn8vuSm0y4+nI4vTqXJ3H08d/tOFReRGwIhEAIhcFoCSSqclm+kz5jAL3/fdT/8QomEv3Ta/+WX2v/jjC2+bNMymT9//8eH5+/DIQvm7uO52zfk35wPgRAIgRCYJoFXSyrc6acmFstusVh097dXJ6Fxrz6uVyeTfxKlJygUX3Ufl1u1W32ST6/WPl0uFEpqOf0/6X21uOq6T6umsmPt98Y//LnrvvryV13359958Yt99b26Vvzfd8uPXVe5LW7uurub5mX5QtaYglYf29odu/9tfXFu6Lp1f4jf3c1pxg/06dv2Teb9OtD4tqv/sH91Lb/r6Pnf1e39+hqT7KFr8nnLHFUC8WElMr/jufu4z775eTIWhUAIhEAInBOB5url2DctJsWnSiYIOBPzUyYtzsmx++rKAm5oAaN6y4WSDnfrVIISRqf077727NpurP0u96vf/GX99Yev9KTClqSCL/C1L8Y1qeByD91/jT4O1XGofcsfznGo/THP942LyxtpqX+7/7l9raSC20rdJBZ250yL+BAS893O3cd99s3Xo7EsBEIgBELgHAicPKmghMKVPmHr+VT4mJBILJzD4vaQhcgYZutlzs32Jw365IxZvGix01oEbZN5V/RZ3t6/Slz06dRXPsZ+tdVXHvQ+hb/+5lfd0JMK+jR7qacjbPHpC8Y+XQ4tf40+DtVxqH3LH7vG31AfY863JvPy6/3tfteZ96lxsl5PjGeeRHgLu13Pc9+PD8/dg8P6z93HLfuGqaRGCIRACIRACJyWwEmTCjwW/FqLfCbhr9Xfvq55jYXeIY+xtxZx1dbWIqjW8eOqj2Jj21cQvO1r74+xXzp98e0f1wkF7Q8lFWS/Lw7V5jXi4DX6OLV/Wv7gWq9MT6lLazIvv+qJk6Ene4b0al1PrWukVTYkO+efCMSHTyzmujd3H7fsm6svY1cIhEAIhMD5EBiVVCA5wMSZCb2Oty3gNVHW491Xt6vHife2+kyY7z6p/ub9C/XTOy2SmMTX/tFLfbjO/qQEi4Bf3Wgh8PwdD5zjO9EsZiRPtqg/lyU300Z1/OkD6qMrdlOf8p//+82CBLn0qfMsmGSzzi8fOMKk1Ye3pw+XQ1l97N653n0Ul+Vj/x7OMEaOtuhTbVO/LX2ubv72mQ+/+ndi8OULvsjDH66Hs0cX7NS5qqfLaC20t9mPDfKBYpMnLnhSQU8rfPUL/QLEr5ovaySu69M6VQ90JlY4lu6u3zY7W36tZR43+M6Z0b+Xwdg56rz/uVzVHysHvm6X29sXj+rP9XddTrG/bTKPvtv4bNOpZUuND7VvsdomN+eeE4gPn/OY49HcfbzNvjn6MzaFQAiEQAicB4FRSQWZosmsvkOvb9BfacH5abVeaPvCwU1mQXR1rcXF5vFgJt6tNpzTYoRFEIkBJup1kk0b5NFn37H0k0z0YYJOfZ3XQlb9cY7FlOog3/VhX4sC19s/hfdz6qMeyw6SCjpP31pkYaP0oC/s6Ouj6qn6kuPtfRFTF71Xt/r09aqZVJAs/Xl7HUs+9utYHFnw9enjNlPHfeG+Wndq/43pD1kwdJ6u65D9skVt0dE5mkq9u1VXKno5slmM+7H7vsUVO2njtnkf6ld+o76OiXfa0r/OUUZ954gNbF8yfLoetskh1tUvSYWXstrxWG1Dl1Nthybz2IIdu+hRrye1bdl3SB+76DPXuvHhXD37ZNfcfTxk3xOJ7IVACIRACITA6xEYnVTQBFeLKT6N10Lh6nrV+534upCQSZSx2Kxmqg9faOq8JttaJN13d81fdti2yFL7OjH37/EzQWfRpPossrQ/dL6lr9r5n+RpweQLPWza9LDR0e2mXxYnQ/3UPuDM4pdj6eH/6FPtnYEWdbsmFbDZ+8LPlKGP6rZs2uYr5NctsmUX/UkO+7U+fbs/huzX+Srvv9z99/X7FPSUgv/7D3/3w4suWwvGlh7Ygi84dm4uq/KSzOo7r0NceQxoHxb0R/9DHF8Y+lBQ43FIDnoR70P+oN+afKD8VNsxk3n5R/7CD5U1x1x76Op+pcx9RxlyYUV5tuMIxIfjOJ1zrbn7eIx95+y/6B4CIRACIXCeBEYnFTTRv9JP4T38zKAmt/7TeEyWtdWkui5QwNOaPHNOk+g62WZirSckWr8wIHl9iyLJZcGihZL2eXS9nkMH2ckiztu2zmMjtrMYU13ppXLxwgaSCK6z6la76RfW9Ty69PWBXtjBMfJory19ue6UterTVn3XhfY6Rh5eyunn6R99JKNlE/22fEW/vu3rz/v2+uy7P7xPzlOG/Szk5M8an7TZtu3Tx/VQezjhC46dm8vSPrFP/+iODO9D53jKiPq+pT/ael9er29f9VsxPyQHncWbfXRQX5ThD/qXX/xpHcpPtR0zmRdv99dYXVqMWva1ysb2kXpdFx/OPwrm7uMx9s3fy7EwBEIgBEJgagRGJRWY1PtEX5NgP24Z1poo+6K9ttGEvC7akHHf3a+fVKgLWc5rSVYXRchXHS2+lh+fvyyvZZfrN3Qe+dqqnRZUYuILOZ2rx+hDkqHaTb8sour5lkzvAw4sbjhu+Ytz1JXs2r/K6p9z17ltx60+WjYhp+WrXfp3FrWdjv18S7c++ymvMdrqw8u8v23l6IKfOHbfOGfJZRGPXHT02CHxUM/Rhi390X+f3tT3ba3rx77vbdh3vdDBbfbztNF2SK7XPcb+mMm8/AP7Xfp0v9IOFi6vVY/62Q4TiA+HGZ17jbn7eIx95+7D6B8CIRACIXB+BEYlFfTpGE8lyERN8vWJ5+rj6vH9BC3TWQywSBn6lE2LBBbm9OMvDBw6zyTcFyQupyYkan3pp8f+pYNkVP0ly5MO0gfbOKdjlfvCU21Y2KmezutYX+m4u109PirN4kELB+nAgrHKQ8a2PtBTNnUfV+s+kaf2+ltpAaR3ZNgnzCrnuNbftNr8rzpukx/D1XlXfTwWVrcirVZPCQ1v6/2yv60//OZxoP54l0blib19/K/0ToWHn0SV7G2f9qOfb91WL696wA2967HaSlfYcL76gfOqX/tQe48b1YE/8uh/iOM2Wzzmh+RwnmtpyB/06ywoO+V2aDIvP/s4uYsu4uV+o638RzmciFPqZDueQHw4ntW51py7j4fsO1e/Re8QCIEQCIHzJjAqqaDJOxN+mavJsxbfTHa3IaCuFqh1MVPbsQC6utn8KoPaeL+qrzoq5x8TbBZElLMwoo9qA+UuT/YsbvR2/6eEAvKkhyb+HKvu8uP9+ishXqbFMZN/ypcP9mA/TDiWLs9kq6+HX19w/bz+UB8szNxHlKGXs+3rH05sK2fJkhxs0rH0xIcseOkbfVx/fEgffb7ivLZD/bl86YQefTy32a+kwuJ6k2zCXtdlzL7kO++qx71+8eT66VdS/km/TmLHikkYuj3VH/CVTrUPEiNuq2SJf5XD9dPHsdr8ol6J+Rfn1e/6ta+6nl9e767jyq4H+kXfGjucP8W2TuY9BsXR2Y/tv3LxaxwZ7vfXtJf+57SND+fkzbYtc/dxta9NIaUhEAIhEAIh8LoERiUVXkslLYJak+pj9O8vaDyGvMg4HYE5+koL0KV+UWOdcjodu0uRrLFin0X8IXwymT+E3jTaxofT8MMptZi7j+du3yljI7JDIARCIAROR+Aikgr6NNBf0Hg6nJF8KIE5+0oLYZ6YOJTTJbd/K46ZzJ9/1MWH5+/DIQvm7uO52zfk35wPgRAIgRCYJoHJJBW0UNAjxPp3jKcVeDxa8rKQm2bwodUl+UpJk9X18vEdEjDIdhwBPfGx73sLxvXQXyuT+X4253ImPjwXT+2v59x9PHf79vd8WoZACIRACLwlgckkFd4SQvoOgRAIgSECmcwPEZr++fhw+j46VMO5+3ju9h3q/7QPgRAIgRB4GwJJKrwN9/QaAiFwZgQymT8zhzXUjQ8bUGZWNHcfz92+mYVjzAmBEAiBiyGQpMLFuDqGhkAIHEIgk/lD6E2jbXw4DT+cUou5+3ju9p0yNiI7BEIgBELgdASSVDgd20gOgRCYEYFM5s/fmfHh+ftwyIK5+3ju9g35N+dDIARCIASmSSBJhWn6JVqFQAhMjEAm8xNzyB7qxId7QDuzJnP38dztO7Nwi7ohEAIhEAIPBJJUSCiEQAiEwAgCmcyPgDTxKvHhxB10BPXm7uO523eEEIiIEAiBEAiBNyCQpMIbQE+XIRAC50cgk/nz81nVOD6sROZ3PHcfz92++UVkLAqBEAiByyDQm1TQjSv/wiAxkBhIDCQGEgOJgcTAdGLgMqansTIEQiAEQuCcCPQmFc7JiOgaAiEQAqcmoEVV/s6bQHx43v4bo/3cfTx3+8b4OHVCIARCIASmRyBJhen5JBqFQAhMkEAm8xN0yo4qxYc7AjvD6nP38dztO8OQi8ohEAIhEAJd1yWpkDAIgRAIgREEMpkfAWniVeLDiTvoCOrN3cdzt+8IIRARIRACIRACb0AgSYU3gJ4uQyAEzo9AJvPn57OqcXxYiczveO4+nrt984vIWBQCIRACl0EgSYXL8HOsDIEQOJBAJvMHApxA8/hwAk44sQpz9/Hc7TtxeER8CIRACITAiQgkqXAisBEbAiEwLwKZzJ+/P+PD8/fhkAVz9/Hc7Rvyb86HQAiEQAhMk0CSCtP0S7Q6AwJ/7brui1/8rvviy2/W//7w7TdnoHVU3JdAJvP7kptOu/hwOr44lSZz9/Hc7TtVXERuCIRACITAaQkkqXBavpE+YwK//H3X/fALJRL+0mn/l19q/48ztviyTctk/vz9Hx+evw+HLJi7j+du35B/cz4EQiAEQmCaBE6eVLi6ve8Wi8X63+pa3d2fjIQkL65X3f3t1cn6uATBd+K4WHbdx+VWc1efVO9q7dvlQr5Vy+n/Se+rxVXXfVo1lR1rvzf+4c9d99WXv+q6P//Oi1/sq+/VteJz/+tA15RfS8eO+6nLewG1p2D5seu2xaXHgeoOjRt9k3kf44auGVeVdtt09Prs4x+PAc7VLbHMGLyLflXWHI7jwzl4cbsNc/dxn33bqeRsCIRACIRACJyWwEmTCpo0M1Fncsvxsc1ioq3J86n6OLbOU5SHn8Rx2wJE9ZZKPHR361SCkhBz4D7WfvfdV7/5y/rrD1/pSYUtSYWhRa7L7NtnITpmQdknY1v5pVxHJMSuFOcPySWx3ca1NZl3nxI7264bsaferskEb6vrc5uuqitfLq+VOFOPXUfsDOm3rjzT/+LDmTrWzJq7j1v2mfnZDYEQCIEQCIE3IXDSpML9J01rn/40qb27aXb5VOmAPRZE57C4Xd5oor+Z7B9gcm/TtfSb7U8a9DVW26EnFbSYGlrUuPyWPsvb+8cFndd96/0x9ktHfeVB71P4629+1Q09qaBF7FJPRxzB50OL30P5ndN1dIit/qQCcraNUXUyDydfpA9dF8TWoePgmBhYfeyeXV/oe2jfsDrHbXx4jl7bTee5+7jatxud1A6BEAiBEAiB0xBorvBPddO60kK655HzY5jHpHnqSQX/dPMYdrdkLG7u9k7gsPDxxVLtY8yixttUfVoLOq//lvtj7Jd+X3z7x3VCQftDSQXZv43nLvbuyn4X2ap7LtfRrnbV+q0YlO+v9AROY5yq42KrfavM+1Uc7POEgsvQ/r4xUK/DKnfux/Hh3D3cdXP3cbVv/h6NhSEQAiEQAudAYFRSQRNlvjuvhRGLjl2+ajBmUcWE/O7Tav1JeesRXy3IVc4/TyCgl8pcZ3/EWXqo7a9uJGO53kcG5/gkj8WlzmsSr3YuSw6mjer40wfUr3pSn/Kf//uNLcilT51nESqbdX55Ky5Pjz23+vD29OFyKKsLG+d691Fc2u9UgDFyXJ9qm/pt6XN187drO5Dx1b8Tgy9f8EUe/qgXFOeRg52qV/V0Ga2kzjb7sUE+UGzePTwBwpMKelrhq1/oFyB+1XxZI3FdF6rIlf7Sj2O3o9qs47qgrO2Il3tdR9fPYwZ5bq/65xpwdq0yWPs56aN40hdhuKacd9VPfXj/xD661e02PyO7T6/6NA3XjPpcfVw1362h/lx/9KmTedlQryP0afmQc6tbvYdkM/bU9vQ1tK0xMFRf57kmWrqNaT+HOvHhHLy43Ya5+7jat51GzoZACIRACITA6xAYlVSQKpoQ6zv0WjjwSZ4mtr64aKnMRLo16ff6vshgok1igAl+ncTTBh2YNPcdqz/JvL/dfC0A3aiv8ywoOOd6I9/1YZ+FlUjVRaSfUx/1WHb4woq+NfnHRulBX9jhLxt0mVVP1Zccb6/6LLhe6qvF4dVjUkPt65+31znJx286Fkfk9+njNlPHfeG+qv2P6Q9ZMHSeruuQ/bJFbdHROVa9WsdVV9VBluu4ul02F7lVprMnVhQfbp+O4Sv73J9VH/igS9VNx1f6bv7Dp/f0qf6kC9cITCWPvqmLfrLF+6cvYqXa6nV1zuMK2dJDf9KFftHL5bqstU037aSC83V96mTe5VGv6kS5tqovDpsXdd49xgDcvO7Qfp+O29opDg59Seg2+edwLj48By8dpuPcfVztO4xWWodACIRACITAcQiMTipoQqzFFJ/G14XGkDpMtrdNoNUHiwLkafKsNvf6FLTxyw6+yGCBwuJIMurE37/Hj05eX/JYNA6db+mL3mwlTwsJtxub+H59lUO/LJbqeWSzrX3AATs4ZvHHFtZq7wzqIpR+fNu3qPG+WNBRhj6S07JJZc7JfeV9+z6yN4s1hfP9WjZ9e132az9D9ut8lfdf7v77+n0KekrB//2Hv/uBbh63LVY1kYEdtZ9HIbZT5Q3Fi5+nH/e3RLuNtY54ETO+xZ81nqttrf6JbTNr6y46uZ/VLzq0Gjsn2nu/rRiUnBofyK6T+Va9bddOq77K1jb9o55c2PySijNmn2sVXdw2yoa2p/762VD/UzgfH07BC6fVYe4+rvadlmakh0AIhEAIhMA4AqOTClp0XOmn8B4+FdTkmU+9mBgzAe6b6Lcm1a6mztfJM214tLouhnxBw8LB67CgUZn2eXRd/fo59JCd6D90nv6w2/uVXioXL2wgieA6q99qN/3Cup5H174+0As7OEYe7bWlL9edslZ92qrvugBex8jD2/T9PP2jj2S0bKLflq/o17d9/XnfXp9994f3yXnKsJ+FovxZ45M227YtfVr2e+zBjNjyfqu8qm+V7efZd39Ld8kkoUPf1JG86mu319uqfFtSgf5h63L69vv87LxabZ0TPvR+q57I8PigTNs6mW+1b5UhoyV3Hx6S57Yhf9tW9d32bXXnfC4+nLN3N7bN3cfVvvl7NBaGQAiEQAicA4FRSQUmviwyZJgmqX48xlhNuElEtOrXxRD9aEFz392vn1SoixufXNfFEH2ojhZMy4/PX5bXsssXKkPnka+t2mkBKCZ18VCP0YckQ7WbflkE1PPqr8r0YziwiOe45S/OUVeya/9uJ/vOXWXbjlt9tGxCTstX9Mt2W3/Ogvq+9fMt3frsp9wX+C63b9/7o05dfCIbn1Ovta2217aVrZ/H3l2uo5b+rpf0ISGh8m220X8rFl0m+9VWP/Z96vvWz8PA+6160tbbUaZtncxji/usr63ao4PX79PB+23tb+un1pf/3O56/pKO48P5e3vuPq72zd+jsTAEQiAEQuAcCIxKKmji68kATY71XgW96Iz3E4wxVgtvn1DXNpr8sjDXuToJHzrPJN8XyC6nLqRqfdnJI8iSQf8+IfekQ52s65zqqtwXnir3RZfO61hf6bi7Xa0XYf6ddy0Y+IRavKo82VTLah/oKZu6j6t1fcl0/it9eql3ZNhTFZLNca2vc/ypjtvkx3B13lUfX0ytbkVarZ587m3p07fb+sNvHgfqj1it7LAXNhxjvz82Ltm8U8T12bbvtlIPHYkt+iTu7m9XVH2xddt1ssqq9nEe+3RettF3PY//4Mex+0Rt0LHqU+1FPv2pvl8fsoGvVVVjXXbVQ/3oukHuWo5+RtHeseAxqhj0fnUsDl4mGcRq1aU1mRdLuGAn/df2OnZ7tvXVautl0pF+vbzuSz/8yLkxXyui7ty28eHcPPrSnrn7uGXfSwopCYEQCIEQCIHXJTAqqaCJsE/cmcxvm9SyANCknX/bJtsym8XQ1c3mzei+8AELC6Iqs/ZXJ9LVhpY82bO40dv9nxIK9CP7WYSoTHWXH+/XXwmhDjxYXFC+fLCHxQv8OJYuz2Q/vB1evNxerz/Uh+xFTxbslKGX+7SvfzixrZwlS3KwScfSEx+yqKPvFqMaF32+Qgdth/p7wUcJle6ul+c2+5VUWFw/fd/dublO2/Ylv7arNuhFiLLd/VxlwlGcxZavBeHTv7150lNy/JdU8JVkelypHB9U/3Id9fGs+qw+bV7yiT78kgjHMPB2rle1tzKqceXnJYf4cvnEYLWt9WJM2dmXNOqbzHtfcJQd9IdO2OaxBl/ODW2rHzxWan+uF/yd0VBfczwfH87Rq89tmruP++x7TiFHIRACIRACIfC6BEYlFV5LJZIKvGX+mP1e8qdzx+T4GrLm6Cstfpf6RY31cwWvQTF97ENAC/G+hX4m8/sQnVab+HBa/jiFNnP38dztO0VMRGYIhEAIhMDpCVxEUkGf7vkLGk+PNT3sS2DOvlLSjE/N9+WTdqcjoIRCfarAe8tk3mmc5358eJ5+20Xruft47vbt4uvUDYEQCIEQmA6BySQV/HFsf6R3X1Q8CqzHfbOQ25fi67S7JF8pabK6Xj6+Q+J1CKeXIQIaf/qeUKBtJvOQON9tfHi+vhur+dx9PHf7xvo59UIgBEIgBKZFYDJJhWlhiTYhEAIh8JxAJvPPeZzjUXx4jl7bTee5+3ju9u3m7dQOgRAIgRCYCoEkFabiiegRAiEwaQKZzE/aPaOUiw9HYTrrSnP38dztO+vgi/IhEAIhcMEEklS4YOfH9BAIgfEEMpkfz2qqNePDqXrmeHrN3cdzt+94kRBJIRACIRACr0kgSYXXpJ2+QiAEzpZAJvNn67pHxePDRxSz3Zm7j+du32wDM4aFQAiEwMwJJKkwcwfHvBAIgeMQyGT+OBzfUkp8+Jb0X6fvuft47va9TpSklxAIgRAIgWMTSFLh2EQjLwRCYJYEMpk/f7fGh+fvwyEL5u7juds35N+cD4EQCIEQmCaBJBWm6ZdoFQIhMDECmcxPzCF7qBMf7gHtzJrM3cdzt+/Mwi3qhkAIhEAIPBDoTSroxpV/YZAYSAwkBhIDiYHEQGJgOjGQGWwIhEAIhEAITI1Ab1JhaopGnxAIgRB4SwJaVOXvvAnEh+ftvzHaz93Hc7dvjI9TJwRCIARCYHoEklSYnk+iUQiEwAQJZDI/QafsqFJ8uCOwM6w+dx/P3b4zDLmoHAIhEAIh0HVdkgoJgxAIgRAYQSCT+RGQJl4lPpy4g46g3tx9PHf7jhACERECIRACIfAGBJJUeAPo6TIEQuD8CGQyf34+qxrHh5XI/I7n7uO52ze/iIxFIRACIXAZBJJUuAw/x8oQCIEDCWQyfyDACTSPDyfghBOrMHcfz92+E4dHxIdACIRACJyIQJIKJwIbsSEQAvMikMn8+fszPjx/Hw5ZMHcfz92+If/mfAiEQAiEwDQJJKkwTb9EqxAIgYkRyGR+Yg7ZQ534cA9oZ9Zk7j6eu31nFm5RNwRCIARC4IFAkgoJhRAIgRAYQSCT+RGQJl4lPpy4g46g3tx9PHf7jhACERECIRACIfAGBHZOKiw/dt1isXj8d3971d11XXd3s9xb/avb+265kCqStPvfoe1rj6tPXXe1uOq6T6t6KschEAIXSmBoMq+xcXWtcez+GSGNJ4vF1eOYeXfTHHaftWkdaJzT2FvHymPJb/U5t7L4cG4efWnP3H08ZN9LIikJgRAIgRAIgdMTaM5uWzctTZMX16vuSot/W2wz0T32RHms6fRfJ9pj29d6TNCrnbVejkMgBC6LQGtchIDSoYvFsplUWN6s065U3XmL7L4x7lD5Oyt0xg3iwzN23kjV5+7jbfaNRJRqIRACIRACIXB0AqOTClq89y209QndvkkFWSTZfRPmMRYf2r72eUbkpAAAIABJREFUcYwnFZa398+SL7WPHIdACJwXgW2TeS3sr25XL5IKGkvub/d/iouEQt/4eqj88/LA4drGh4cznLqEuft4m31T9030C4EQCIEQmC+BUUmFoYmtnmJY3eirAs8f+x2L7dCkwKHtq56HJhUObV/1yXEIhMDbE+ibzCup2n1crpOj9esPi5u7x689qM6uf2q/LeF6qPxd9Tn3+vHhuXtwWP+5+7jPvmEyqRECIRACIRACpyMwKqnAexT0/oQxfz7R1XeAmUyTnFjd6vvFy8cnH0gKLD8+TcDrJ3Oq4+9yQKb0of1dp/bLdT3auy6UadGv7zj7kxfIV9nqo77m8fROBeTrk0jpIDnYgk5M/JFDOYsMvj5C+ViWY3inTgiEwOkJtCbzGgd4n4yufa73qg1jKGNQPd86ZoxhvNTYwThT6+8jv8q4hOP4cP5enruPW/bN36uxMARCIARCYOoEjp5U0OTWJ75a1GuifdcpKbBZ8Pt5AWIhzoSbCTIL7zpZpz6JBY6Rq/aeMKg6qU//HrKf1+L/6uYpqYBsTejRR+1lF/oy+ed87X8t83r1+HUI6qP/1IMk+oVACHRdazLv40gdpyqzXa97jSMad1bXSuZqBN2814Zx7lD5tf0lHMeH8/fy3H3csm/+Xo2FIRACIRACUydw9KQCBjMB3kyK1c3946f7LL6pq8l4nSjXZIQvwJHNor62r18/oD79anLPp4ucc/k1KVDlo7e2PPXgSYfansWB6vg/9Hd52Q+BEJgmgTqZ17XffXz6hZihpIKsUh1d9z5u+JjAvpKiy39cvRgXGUt8vHJayPey7D8RiA+fWMx1b+4+rvbN1Y+xKwRCIARC4LwIjEoqMAEeuwhWQoAnBXyizSd1LO5BpTo1qUC71Sc94XD1+BUK2pB0ULKitpe+/vUFtfE6/hJFbPNJem3vbekfW/QUxr2ewrhePT7J0Eoq9D0WjbxsQyAEpk2gTuY1BpEEqFsfT9wqjQ1jx1HVreMi484x5Ltel7IfH87f03P3cbVv/h6NhSEQAiEQAudAYFRSQYZ4oqAapk/7l3q8v9Mb0J9/r9iPmRCPSSqoP02+aVMn4pyXLnXRX5MCqvMo5/bq2UslKXedavsqn6cbaFOPW0mFujioDHMcAiEwbQJDk3kf6/osUZ2+hEBtw9jk9evYVNvsIr+2vYTj+HD+Xp67j4fsm7+HY2EIhEAIhMAUCYxOKjDB5QkEjNGC+sreF+ALcBbbfEqPDBbjyFAbl1sX5TqvTwKZXNeJtfcpmfW89+NyKK8JEx2rHjpV+dhFokP96WkK7PL+V7f6LvTmSQY4qF+xuL99enQaXbINgRCYJoGhybzGCb/GqxUaFzbvR9AIMu6vjj2eTK0S9pFfZcz9OD6cu4fb7z5xq8/9Oh2KYbc1+yEQAiEQAiHwWgRGJxVQiAU3j/uy8OY8C2wW5Vc3m5cz/nz9PoGfPz4uzAKcdrrRI7P1qb6f9z69XO1Wn5QQ0K9LbN5fQCJC/Wghv1zoZ9209/RHkoA2q1v9MsXm1x+qfNo+K79edIvrp1+dIHkiefTvZSpv2fikUfZCIASmRmBoMq8xwZMKPhbqmvdzu9jmYy6JTLU/lvxddDn3uvHhuXtwWP+5+3jIvmFCqRECIRACIRACxyewc1Lh+Cq8nkQt7HlB4+v1mp5CIATmQCCT+fP3Ynx4/j4csmDuPp67fUP+zfkQCIEQCIFpEriopIK/oHGa7ohWIRACUyWQyfxUPTNer/hwPKtzrTl3H8/dvnONu+gdAiEQApdOYPZJBb2fga811K9cXLrzY38IhMB4ApnMj2c11Zrx4VQ9czy95u7judt3vEiIpBAIgRAIgdckMPukwmvCTF8hEALzJZDJ/Pn7Nj48fx8OWTB3H8/dviH/5nwIhEAIhMA0CSSpME2/RKsQCIGJEchkfmIO2UOd+HAPaGfWZO4+nrt9ZxZuUTcEQiAEQuCBQJIKCYUQCIEQGEEgk/kRkCZeJT6cuIOOoN7cfTx3+44QAhERAiEQAiHwBgSSVHgD6OkyBELg/AhkMn9+Pqsax4eVyPyO5+7juds3v4iMRSEQAiFwGQSSVLgMP8fKEAiBAwlkMn8gwAk0jw8n4IQTqzB3H8/dvhOHR8SHQAiEQAiciECSCicCG7EhEALzIpDJ/Pn7Mz48fx8OWTB3H8/dviH/5nwIhEAIhMA0CfQmFXTjyr8wSAwkBhIDiYHEQGIgMTCdGJjmdDJahUAIhEAIXDKB3qTCJUOJ7SEQAiFQCWhRlb/zJhAfnrf/xmg/dx/P3b4xPk6dEAiBEAiB6RFIUmF6PolGIRACEySQyfwEnbKjSvHhjsDOsPrcfTx3+84w5KJyCIRACIRA13VJKiQMQiAEQmAEgUzmR0CaeJX4cOIOOoJ6c/fx3O07QghERAiEQAiEwBsQSFLhDaCnyxAIgfMjkMn8+fmsahwfViLzO567j+du3/wiMhaFQAiEwGUQSFLhMvwcK0MgBA4kkMn8gQAn0Dw+nIATTqzC3H08d/tOHB4RHwIhEAIhcCICSSqcCGzEhkAIzItAJvPn78/48Px9OGTB3H08d/uG/JvzIRACIRAC0ySQpMI0/RKtQiAEJkYgk/mJOWQPdeLDPaCdWZO5+3ju9p1ZuEXdEAiBEAiBBwJJKiQUQiAEQmAEgUzmR0CaeJX4cOIOOoJ6c/fx3O07QghERAiEQAiEwBsQGEwq3OsnIq5X3WKxePHv7qbZ/KRmLD923XKhfu9O2k+Eh0AIhIAT6JvMrz513WJx9Tg+7jIuahRbLJaPbbuPS+9yvX91e7/1/IsGKegl0OdD7nOra91bdPT87xAfxMfPWZ76aKo+PmSccGZ99nmd7IdACIRACITAaxNoZgVaN626mNexEg27TKBf27j0FwIhEALHItAaFyV7eaNl4+5JTi1dl9erx7YsXD2x4OMui1M/fyzbLkVOy4dw1f2slVQ4xAeSvVwoUbSJD/m4JsUPkX8pftvFzin6WPrvO05U21v21To5DoEQCIEQCIHXJrB3UkGKLm7uuis9NfBJE+P9/451s91fg7QMgRAIge0EWpN5ffp4f/vy6YLtkjZnVx+7Z2Mnn5aTqOXYkwhagLYWvmP6S52ua/kQLlrwV7aH+uD+k6RLyuZPqYUrJRke7pmHykdutk8EpuZjaXbIOPFk2WZvm321bo5DIARCIARC4LUIHJRUaH2ytqvi/inNrm1TPwRCIARei0BrMq/EKl8N88X/vjpJHkkFLUSuFlePC1DJbJXt29cltmv5EA6tpEKLd6sMGUNbtb27uXqs1pLVKntskJ1BAlPzsRQ+5jixzb5BOKkQAiEQAiEQAicicFBSwZ9U4BOX1c3i8R0MTLJJPtTJt99o/asU/jiqypEjBnXip2M9TnrXaXK/+W4yk/ITMYvYEAiBCySwbTKv5KiPYfvgYQxlvGslXBkbqbNPP5fcZpsP671FnI7pAyULVtdKKDw9uXBM+ZfsV7d9aj523Y4xTmyzz/vKfgiEQAiEQAi8JoG9kwrcHPW46H13/5hIqF+HqBM1HXuigKQA3znVpNm/g+qf2tCWR1Q5ljy+pyq9qg6vCTR9hUAIzJPA0GT+0AV/XXRmwXn8ONrmQ91PuLfQ87F80LpXqY9jyUffbHf/isuxfNDn4+qTQ8eJbTFc+8pxCIRACIRACLwWgZ2SClq8+z+fgPEpmz8l0Lp51nq6EZMQkNF+Y/a+7m83j4zqvPdb23sS4rUgpp8QCIH5Exgzmdd4pDFQ45D/IoSPZdpvJT6v9MJHez9Na7GD3DypsF+8bfNhvbeoh20+WP2jnox7+tWPUT5+SKpzP9smPz6eh49bVjBOtM4NlW2L4aG2OR8CIRACIRACpyKwU1LBF/9VoZos0Pm+CbC+9kBioCYFhm62deJX2yepUD2T4xAIgWMQGDOZ1yLRE6tj+9U4VheRrbGsVTa2j9Tb/VPsFu9W2Vi23CdJKrRktcrGyk+96fm45ZN9xwnJGjMOtfpMWQiEQAiEQAicksBJkwo8qVAn2UoqUFaTAjVpUI2v52v7TMgqsRyHQAgcg8CYybzGo5ocGOpbCwwWmV6XBajLq+Of18/+MIFtPmyxPYUPdP/Dp6eQP0xh3jWm5uMW7X3GCeRss4862YZACIRACITAaxPYKanQemQXhZkckSygXDdPPRbKJKou+nn8816vWrxdrX/NWy9cdDlq033c/Gyl5PkTE/W4ykePbEMgBELgEAJDk3mNPfVFfEP9tT6xXCox8fA1CJ3nqS6StIylQ7Jz/iWBbT70J+i85TF90IqRY8p3vS91f4o+dl+0YsDPD+1vs2+obc6HQAiEQAiEwKkIDCYVSBb490Xrp2q1jicEpDiJhdZ3iXWD1fdSPWFBGX0yqXY5Sixc3f7T4zsedLz6pF+AePqOaybfpwqbyA2ByyNQJ/N949RYMj6eMdZpy3iHHK+XMQ0q+22rDyWFZA0+8HsRvezrAyUMkNvy7aHyaZ/tE4Gp+fjQceLJss1ey75aJ8chEAIhEAIh8NoEBpMKr61Q+guBEAiBKRLIZH6KXtlNp/hwN17nWHvuPp67fecYc9E5BEIgBEKg65JUSBSEQAiEwAgCmcyPgDTxKvHhxB10BPXm7uO523eEEIiIEAiBEAiBNyCQpMIbQE+XIRAC50cgk/nz81nVOD6sROZ3PHcfz92++UVkLAqBEAiByyCQpMJl+DlWhkAIHEggk/kDAU6geXw4ASecWIW5+3ju9p04PCI+BEIgBELgRASSVDgR2IgNgRCYF4FM5s/fn/Hh+ftwyIK5+3ju9g35N+dDIARCIASmSSBJhWn6JVqFQAhMjEAm8xNzyB7qxId7QDuzJnP38dztO7Nwi7ohEAIhEAIPBJJUSCiEQAiEwAgCmcyPgDTxKvHhxB10BPXm7uO523eEEIiIEAiBEAiBNyCQpMIbQE+XIRAC50cgk/nz81nVOD6sROZ3PHcfz92++UVkLAqBEAiByyCQpMJl+DlWhkAIHEggk/kDAU6geXw4ASecWIW5+3ju9p04PCI+BEIgBELgRAR6kwq6ceVfGCQGEgOJgcRAYiAxkBiYTgycaD4YsSEQAiEQAiGwN4FmUmFvaWkYAiEQAiEQAiEQAiEQAiEQAiEQAiFwMQSSVLgYV8fQEAiBEAiBEAiBEAiBEAiBEAiBEDgugSQVjssz0kIgBEIgBEIgBEIgBEIgBEIgBELgYggkqXAxro6hIRACIRACIRACIRACIRACIRACIXBcAkkqHJdnpIVACIRACIRACIRACIRACIRACITAxRBIUuFiXB1DQyAEQiAEQiAEQiAEQiAEQiAEQuC4BJJUOC7PSAuBEAiBEAiBEAiBEAiBEAiBEAiBiyGQpMLFuDqGhkAIhEAIhEAIhEAIhEAIhEAIhMBxCSSpcFyekRYCIRACIRACIRACIRACIRACIRACF0MgSYWLcXUMDYEQCIEQCIEQCIEQCIEQCIEQCIHjEkhS4bg8Iy0EQiAEQiAEQiAEQiAEQiAEQiAELoZAkgoX4+oYGgIhEAIhEAIhEAIhEAIhEAIhEALHJZCkwnF5RloIhEAIhEAIhEAIhEAIhEAIhEAIXAyBJBUuxtUxNARCIARCIARCIARCIARCIARCIASOSyBJhePyjLQQCIEQCIEQCIEQCIEQCIEQCIEQuBgCr5ZUWNzcdXc3r9bdaAcuP3bd6lp63Y9uc4yKV7f3e/NYfeq6q8VV131aPVPlruu6xWLZLRaLvWU/E3gGB2KxWFwN2izeb+HnbQgPiYFtcs/h3Fi/HdMWjUH3t1fHFHlSWfvErMaAq8XyxdhwTEU1Ui6uV+tr7mqxOGlfu+otHw+NfxrzVUf/phYP28aEc4vfXX2X+iEQAiEQAiEQAudL4FVW+SwglpqAdpr2zvdP1t3dLLcayKR2nyQLbetknol+93HT9zlNQGUTem8FV06K9VILqO6uayValrf36wUPzA5NKiCvqPHi8NQx8KLDMysY8tspzNFibYqLyD5bjxWzffIPKffrYJ/ExyF997VVTCmhOjSmapxYXSuxdN+Jcd89qV7D+45RffruWk6yZGpJkF3tSP0QCIEQCIEQCIF5EniVpIImoauP+mTraq/F4zmh1+RvaGIre7Z9IjVkrya89dNITXp90eyT5yF5b3mexcA+SQXZuFRMNRJVNclw6OKnytvG7DViYFv/Uz+3zW+n1P2cEm3icGjMnorl/aenp7rkSxbpp+pvSC4J1TEL7jpO9sn2a/iQMapP/j7l5xa/+9iYNiEQAiEQAiEQAudJoDepoAkMj4j6YnVXMzXhW92suvvufv3IbF1wM3G+ut30558cMbFDF/Rgkre8ef6ov2ShM/14mSadHPt59lk43n1SAmQjW4tdTUQl158OoMzL0VNl2CGZPJ6P/mIoPdQvtqgNeqBj3yRZbWpSAXn4p1VH5+jv7qPs29jo/bhd6AOXqwfeqs9EXnrrH0kBL3de+BLbxOJOMfGgg9uPDWydK7oih/4pVxtnrvM6p/rqc/lxEyP4R/X7dKb/ljydcx1g5brSh5e1YoB+fOs6OV/3zxDf1ld6+uSqb9cTnsTLveLl4ZF3nXMm+L7K8HJsc2bum+X1oltcLx+vMfpVHbeTWBxzjdInW9mHXSpzXfCfylVPx5x3n7nd0k3/aOu+oQx9/dpBH9+22Kv/vpj1viofEm0wrL7zftmnf/mMdm439XyrenUc8vPah6FzUrn6W96I3/JxrPS2u+gjFsSP+unTG5n4Tba2OHq9n699/PPH+8rf/Z/SefN1K2KJ+jXe67VGTMhO79f1lSzkOj/VWeg++vD1HdrreHUrT6i3/IVACIRACIRACITA2xBoJhU0EWZiyqS4fn9/rLpqf3+7eSRfE0wWWmrPxEiTPE2OmIT5hF7nfKL7b7dXj4sbn6RJNsdMipmAeb86x9cT6F/taKP+0FHtrrT4vdl8D5sJH3pKL/bpS3Vcj9W1bFdKZfM9ZOq5vjq3vNb7EaTF5v/72+fvS3DeqlUn8+hGvVYddHAbxYBFiXyFrfDoSz54f7KFifHVzd36KwfSg/Kr2397nJTLfmS7X+uEHDu8H+nnT7tId3xFfbY1bqWL7JZv4PDoi4bOdZLekud+1uIIeR4Datd3LXkMoDdbtxuOJOZq3GGb7Kt8kce2JVe2ejmcfTFKjKzj5VoLt80j5G6D9vHjNt/4Oe1Lb+RLzxq7kis/33VKPG4SYfhd51rXKPb61m10vYlHsVO59GmxpB42Sh4xIGZjrx3XSfuui7NBlxqzxG+NA/QTm3uxKu89cNlVh7UesufTav2VIrHmq0Wtui2/1XpuF7o5Y3xY23GssWSMPvLDSmP0J2mteHq6FpHFVjoxVvVxVF33LTLxu44Z11VXvu8+vhyvxYi+nL2PCS4bfzOOePvHa/Jh/KR/twcbsw2BEAiBEAiBEAiB1ybQTCq4Ej7p8fKx+ystNB4Wy0yMmJxJRp0U+eTLJ3Zel0+4keOTQ/RyOSqTLH2atVwv8mXV5k/9+8LAX4DokzrVrvogw8t9X+1ZoLBlkun9So6OmUyu9I6B8hJG+tJW2u+TVKCtJt2VnU/20VVbsZHfnIv6by04iBVvz2LRubi/aIM+Q3Z6vFQfe9uqs7dTPekjm+m/pXOfPNefOq6L28p5bekLW2sMUFf1Wnw5z9b78f2WfvTfkqv+ajzBq15r4kqSRDKpR8LDOWofW9FZW2flMkjk6DzXic5jj/xV/VrrOgfvU/s6JxnIc91cJ5fhdavt3rc4VNtb107VSfKvlFBsXO+whQv6Vxmur+tY483PVRmqq4Xq8maT3HTbal0/Vt+t5IBzo34fY877dqw+rX4qN5e77Zxz9P3KUfJ0nvjxe5z35fset9JBcdj6k9y+GOVcS5+WrJSFQAiEQAiEQAiEwGsR2JpU0CSQiTITqF0UY/KDDLYs4iWrTvJ84usTO9VlolsXOvTjOroctWXy6X3TP2U+8fP+fFJPXZ2nX9lFuess2yhXff+r5ySLBd/Q46yq21oEel+tOuofnZ0Vk1W2rqf2K5fKlvp95TrvXPCFdGjp4/I8maFyYkA+8UUKbdhWnWucYes2nZGlrctr6exy3FZktK6lGgPUdVmU+Zb+++LO+Xq7PrluG/XhXK+1KgOuqldjEll1W/2GDK6zeqz27i+PCfSkbYs9/SMDfn4NuF0uw1n6vmS6/5BNX2wrW7VhHFTCTZ/E93GrHGof2OFx4HZwHjv9HPqxJT6pq761uK361uSH+mjpX/tWP96/M0YH347V515fVVNS5uEFtZJRY8LlVqY6h67O0fXjvPchW/S0jpJpd1ueKpN8GHqCdSip0OLq/oePdHa93Nbsh0AIhEAIhEAIhMBrEWgmFTRh0mPmWqS2JlRjldPEp054VOafbNVJ3raJJ5PFutBhsu8TNSZ9LDZ4lLY+Gqv+WYyrzZgFC/1hh09Afb/a5ty8X8rXbfWVioFJamvCCRvsrcf00fInk1XXnfraVi6t/lv1XIbLhp9io6UP7VrnZBf+0j4+oA3bqnP1BTbXerSvW6+H/n3x5raqXd+11IoB9dvHl36x2fvxferVa69P7jbO9ZzsaT2pUK/Jys+Pq9+qb+p5tZV9ssf9oHLV9acanIP3iQz5DD5j/EddWOJPLebG9Fv1rTpV+X6+cpFtrn8rDtw/Q77zvpCtMunU9/SEt9lWF7vGMK4ydbyLPl5XbRUTjBFVtjNFxxZHyURG5SiZMNILiFtffUAPxYkYeBz0Xfdqgy2tPjnnNslWkhVenv0QCIEQCIEQCIEQeE0CzaSCT3pak5sxCmrSpRc0ssilDfKYbKovJnWq4xMn7TNxZwLYtxCVHJ9cVTl1UcCx2+oTP+nSt2DxiTs6MwFVv+yrnhaU2Cob+FTL+4VN7Z/yuhXD+umg80EnbPT2lb/bInv9ky/VVYKjpZfbKfl6BJjvceMzyrVM9vquK/q0dMUO4sPb6Zz05ZyO/a/qLN6uF/GBTD/Xepy5JW9bvBED7udqq59z3bXvvHQsnVaf9Kj5069deB3fx6YWU6+HXPwDS29fdfZ4UXvnKtnIWMvu+RpP9ZvLUDv6hKH3Wf3Qd41KTv2TflyL6nOM/yoLvsteZe9y7dS21X6+/lTL0d95SJb71M/BkTjwc66DbPQkgurxzgyv19qX3fipnh/LuLbbVR/pgC/dX1Wujp1p5eEcfb9yRK769XinnC3+0rH6ImHt+zon+SSSvY32kY8OGp/1nhvisLKi72xDIARCIARCIARC4DUJNJMKmvRoMawJjCZri+urx0nbGOV8AuSLNSZ8kqt/OqeX+OlFa/TnE1RNqvycFgQug4kkOqk+spGjSaTKdOxtVXZ187eP9X/18CZylWsit/z4JGuj59OxvyxyLef66W3gLC7QrbLU48PohF7oL/1aiRjOa+vy1J4Fg845dxZP3tbrXD3ojJ7UQ384tF6Op7qVJf21yp/Zay/aRH/85rGCPi15OucykUObqt/P//0m3lRPcUB/Oq4xRZnLqvKILZdDmerCUGxX/23zU6qSq2Oupat/+PYx9rwt/bbsrmX48Mv/48snWQ2+yKx2uK1VduWytuXhZ2HVrl4j2NDHBB2q35bfbV68KJnIUF2Pc2LU47v2X69R4pF+XS/OeRl9u353haW/uFL68o9rEL/Dp+/aQSffVl1clnTz833jz7f/oHFoMxb97TrW/m6tY8t3m9Fio4FYexKBvlvXo/ug+sztYd/1bjHGF9TXdhd9aOd+a8lUPddFTPjqBH7kelJ7GBB7tHUmYsFLf9HDt8iQfO5jJAn8HGX0ofrSwa9J1VncbF5AqX6vrnXdbHztXJHlemQ/BEIgBEIgBEIgBE5NoJlUOHWnLl+TQZ+o+TlNspgweflc9zVZ5BOrU9nIooCF0Kn6idwQmBuB1vWphR9PH83N3tiznYCSH31ffdjeMmdDIARCIARCIARCYF4EklSYkD+X+qWMxlvgj6likgrHpBlZl0RASc76Kbg+cU6C7pKi4MnWpd7zsX4+7KkseyEQAiEQAiEQAiFwiQTeNKngj4DWJxLGPM46B4f5I651wXJs+7wvHus9dh+RFwJzJUBCjsflta3j1lxtj10bAh4DSSYlKkIgBEIgBEIgBEJgQ+BNkwpxQgiEQAiEQAiEQAiEQAiEQAiEQAiEwPkSSFLhfH0XzUMgBEIgBEIgBEIgBEIgBEIgBELgTQkkqfCm+NN5CIRACIRACIRACIRACIRACIRACJwvgSQVztd30TwEQiAEQiAEQiAEQiAEQiAEQiAE3pRAkgpvij+dh0AIhEAIhEAIhEAIhEAIhEAIhMD5EkhS4Xx9F81DIARCIARCIARCIARCIARCIARC4E0JJKnwpvjTeQiEQAiEQAiEQAiEQAiEQAiEQAicL4EkFc7Xd9E8BEIgBEIgBEIgBEIgBEIgBEIgBN6UQJIKb4o/nYdACIRACIRACIRACIRACIRACITA+RLoTSqsPnXdYnHVLRaLbnW96Faf7rv72+UzS69u79fnVUf/7m6exNVz3cdN27tOcpdP7f7X66f9BzmSdX979awvHdS2qiM9/+122S2uV005kkXfLwSWguVH6baxxbctXUrTNzusnKW3+2FIMWxeitOa8FCL1z1/r3gx3+5i2zZNiaUr2f1pta3qzueOwRT9PA7ZPxaDnQ0b2cDHDvgub2SR/u32t7i5W1+Tu1yDxAx979bjpjY+1NjXdZL4en/oj7/rto5nxMoh9m6zrjXGoNOxYhEbkOvbU9m1zWb54O728HFhn/itelU28v/y9n49bhGnu/iB+JrSmE+M1dgeYkGc7DI+VJnndIzvZPcQq9XH7uj3tmOwestrS3Embs6xNb4Qj/Wct6v3BtrIN6141LWb8/zeAAAb0ElEQVSqvv16rtfgNvnerl7vfo7YUH8tPVo+rO25rthKju7rV4urScZUy6Yplek+UH02Jf2OrYtiZalY2XHOp2vokjgdi/sX3/6x++tvfnUscQfJecoCmBgGNm5aLBJ8gNJF4oMqbbxMg1odNOnGFxk1kBic6V9tWjpQD72qHNrVZAg6tLbS2YO61a/KxgyuGqjvbp4nYlp9Hlomu2GAvs6uTz7+ue/uuqvr1aRvFoo3bOyzh3KPLcpaW7FaXWvgU/Qe5++YTBU/y4XiR3ubP/c1ZVPZMgb4NU+Zjwu76Kv2+8TmMXx7CGsWfbvY6nVbdns81DHoGPZ6/3W/xULX5KhxZmRCye2j/7F2qa2PtWq3y7hPf2xlr98HKG9tt00YW35syegrk+y6sFGZX2Nj/eB9jOXqbdivrCk/ZLsLp1pXtugDkDGx2NJxzLWqPsbc81vyj13Wuk5qH6pzpXvHHglzxf6+43XVo3X8VteWrhPmEFcakx7Y1OtJ93Dsl999cfRMho0RaoNs7bOwd/tXt3fdfXffLTXX6rS/+cCEvlS3Tz51iXGvp3Z9Max6Y8YxxcvqWnMNabjRi77E4J/0wd1i+WIscvumuH+KsWpXOxXviocxfthV9hTrMx77PWqKes5FJyUUvvjym2knFepAKvgqY5Kmi6QVMLqANfAwuLbk4Ehf+NWbTJXDINe6KNXH3c3mqYYqh7522UoGgyntdJOoEzvObdu2ZG2rv+85v0HCrtpQZavevpOOKuvUx/I/N+KhvjwehurKP8TqUN0x54/N1G3Rfvdxtb4OtZ3inyYwrXFB10/r2h1jg+zeJ/FzqG8Vc/skM2ST4oAJ2hgbW3WqDCakqxv5Xto9/zvU3ufSXh75GMMnoH0TWW/tMezlrX2vKwv1pIA4eLKg1U5ldax1ffva9JVLj6vrp/tYX70x5fvGr2TLptZ9x2WKT008jtHrkHiprMf0N1THbRqqW68NxYqeZttnLK+yhvqewvkx46nif7lOmMvC6fy91bXlzNZEPj5x0Z7mrZozEUs+f9J9TcfS3RMMase1d2fyRLsu5lX3XuOZ7uOfnvpe83j49H+bfOnvyQfV5b7IWNnnZfTvO69y9NN+ZYB813WbrCmdq354K92mosdr2V9j+bX6vdR+Jv+kggYYDbI+iOEsBpy+RYIuHhYW2wJrW1JBMnwypQHVj9GlbjXZcb3Uzm8OfZM05Mi21kJCdvBJCGywUW05r2ykmCkHLX46dr2pp7bLj5snPRY3//T4eLf0/fn66xc/f5bp3sZR+vgCxvljl8qki/6Jhx/jY/yKzvefVo9fO1AbsaMudqiu2wsT2aFzPsmjDB3Qzcvdd5yHLzdQ6YFsdKZft4u+8Zf8sPq4ejYp0MRLZf41H91Sq30qcz3h4Dp635xHP+mLjsjWsZ4Q2TYZlq3Ery+SsAmZ6IbN0osy9S2uHLu/JJty1ePf//Q/b/axQ/Jkn8t327WPXejr53WOhCSy1nrdXq3tX/2jrpWrNaPKRPqtVO/henL5Lb70Cy/s+6eGDOlVJ4nqh1hRzNXrE+Zus8q4nq9u/+2RI/Hsevp4gK6trfSGt/Qkcaq61f8qWy8ibjdfAXMddQ4GQ7HQ0oP+GGNkSyvBx3Upm4kVvyawpa8PlXu8u/3eBpl3D7Fzp/HMxtrVf3u6nv+3h5h2HuoDv7hc9hlPr/SpYvkKnnMkDrHb61JP1w9y8Jl8MCYGqO9yXUcSLcQGMea2UqY+a3m9z6EzMaK+sI0ydKo21LYc17FGMvFfZaC+WrZis28l332IXkPjiPTWP/qp9qkPZ4aOyBfDO43Zi+UjT2xFpuvZ2pe/GOvqmCEZ9O/+khx05foi/lp9qK6uCxa8qgN3tUNn7WMbPDmnY3TlPoAO3jdyYdXShzL1pXGEa4JybelX/kE+/Tlb6u1ybYlpjXfvG+bwqB+2rHk+3ENhofbezuVpnzaUk4jlmK2Y0J9s65Nf5dV2xDbskK+t/Mj8ycv79vvskhw9rXOnuWHjfkwsuA+9D9pf3Wzmxh7vauMxJFnEus45F2K2tqnyF9f/y+O9uF5P0ktMJXf5cTOnpA52eJ+UoYvP8bl2JE/nPV6xf5s90rvOQWnHlrivNuMrlesf/h9T7ryRL931AYbPg9djgeai69Hi+bWK7dLT/bK61fjz8usPcIQ1tsNa5+HHuWobTP7Qdd0XX/768VP6H/6s4189Hv9V53/xu+5///Kb7osv/+91eff7X3e//L3qqeybTsfI+fKh7Jcq7/7Y8QQAx7T7w7c6//Snfr/68lfdH/78u0d9JFd/rhP9qfyr3/yl++EX33Rf/WbzlAF9PEl92kMPb48MlUnOF9/+bv2kguSqTDpif233JPk0e82vP6grd2grcHB8VUsXFsEqGa3AUpuaVCBwtK2yXWbtz4+5qF0WF5nX69uXvq3Bl4sFWboAYaJ9LhDV40ZeZdV6GpSxUxeSLsKr66vH/v3ikr7Oy/V3P8lu9FIdBhb6ER/2pQ/71T71LVvVXjfj5e1qrZtkVzuwVzI0YdAATT10Ub8MGt6vl1debqPrrfL1sRId+nTqk3TbfHUDfXmsUTJ5ZBUW6IS+ouTy++yDietV9922ylTn6Btbr26UuOl/bFd9EsvEGH2i/xBv+pKd0mF1s1jzUqwt/3H16BfJrTr6tSs5nhhAD7bSFR9TVrf4gLgTd+xCz8pEcqnj+m3jq3OKS03aayxKHv1LP64ryaafqufatnJ9en10Qa7qM1ZwjmM/V/n4sdjge22RTZ0WC+LL+5Ac/ALjVixwzSDft2rHpEe6II866o++pRfjf70eqd+3lRy3GWaqX33Sip1NredxLP5M2iVHxyzIW3oQD85Q9ZyjbMQf1UZnIW6MP+6vVr+1zO2r5/xY9ZRgqWNhjTvvX+dIEkmW20aM6CFoxlXvj/OwbrVtxRdJari9kGOPo3t/rX33DXFBTCLXxxHXseUv4qwy835a/Fpjb0tfyiSDMYa+4CEd73vuaZ54dhnI9a3s1zgNB/ykOiTRlaBcx+nt1TpBuPqksfvpa4But+So/lo/xdnDHALu6N/qz/XS/ltdW9LN72VVL/mCeUyrrmzW+MaW9jAgfijX1hnqmIWa11mX25xzm3zJg7XaSWcf15ClcZqxmL5qzFPet+2zS2xIiEkD6ct1p304tGKUeK/3MpL/0gV5npiXTNrK/mq32uiaItlX5a/j3N71hs3SkXuNeGGz5kX0CV+PCe0zpquNJ9ilm55Gaf1JD/WHbOb/asNYjA7uZ8miXG3Zp47HBfzUoq+8xVsJEnRAhn9Nx/0p+/mAA7+4XtrXn/pnrKs8dI56OkfCTX27/9wG7de4VlstnL/6xe+67s+/W3ejBMFXX/768fiH33frfZX/UuXdH9f1tNj/5Zd6B8HD8W/+st5X4sCTBlqkc6y+/vCbTT9rITq2xAaJARIGXfeXdWKCdx0oOSBZJCe02Nc5Fv/0g2xtJYsEhdrRh/aVTFAfJC5qP2qPbGS47FPt9yYV6BBHM1EkkAhq6rFVfer6xch5ttxgdKw2Cpg+2TrfF6DIczmU6WLw4KW8b6t+Wnb5QKK2XBTSl4FH5arHzdll1XoMDNKNc6v1Vzh0ZvMHEx2ptO+lYaonOd53S4bK0Lu1j93qi8nk2u5r3aQ2Xy9BVxYg3qdY6xNlJuwMCKrTuqHXcrVvDRquDyw0kCvRIT1cjuujdpqMuf/RSXLoT1tsH7KvpZ9k8Vf5IlfnXTf5jJsXbetWujDR0T6DPvXG8sZOtZN+60diH75L6jEKL3R+wWLLS79oO8THY1r6rG14eBFii4kzUH0xJBbdl5zj2pPcpcaTh/eZeF3XQXrziL1fx9jTd31WNrU+foM5TFWPa0vn+v5qPenscax2HmvuY52jvl8bKvd6NRb6dEEeNrgPkMnEUsfep/bxyTb5Old97ZMfnXe/6bjGDvrpnLNx3+icdOr7+pBkwtn7d5skw//cRtdJdbxdtc9ltPZlg3Nt1SFOlorzMhaqvTOpenKtuo7qw2NEzF2GznvZtrY1vqr/nIfvt+z0MtXlU1JfFFBH/fjYWnVUPeKj9uvMYMu9jjZqv/Zz416HDq2t+uobY9RX3z2tcqsx5n1JDl+Pkt2MlaqztvV62a1utDC5X98buReoD2KtMlmfu9UHIE8LGl0jVa/azvXSvvR+q2trGzPpJlvQrRUv2MoW28S7zjF0rsW+b7Hp48w2+R6b6qPGEzq1yvv0pE3d9tWXXSyo1QauxBPXI1uYIr+2l671WmbtIHsZo9RebBSjSuQRqypHV41TVb7OVzno4jI3UjZ1Ge/QrdpQy6UXbVgcex/sVz1gJ3vgxdbto71vkSVdfK5BnW3lLd48eeK2Olfkait70ZOt/IQ9sJQvWvN+yVBd961/dUi2wbP228fFH/1ngc0i+oeHZEFdXGtR/pWeSHhIRvzwGyUXNgt0TzZINot3JRD+WpIK0lF96kkFZPmC323wZIPvq44nDGiDzko++L+//v7X6ycwsFH1nQHJC5Ujw+si/1TbwaSCOibA5FQGEA8KV05BRwJga2ApW7UeGjeBijwFHAMLcltlnPOt+kaOl4/d9wHe2yjQCWi/YP3icLvhxUWqc36hVBlLLazse3bq22X3DVbqh0yp9n3B4sfqj8yoZHv/vi9/KRPtemN31cnt1Tkxwu/ed7VddfXn5c734fTjRjr5wgQd4emcJBP/e7mEuU46Vp8aFNFZZd6m2kcb+tWx/zlH36cOukkP/+SJ83XruqhNTSqN4e1cpZMmtfo0jD+V+U3GbXNerf6RoS02wd7Psd/HRH3S3vtXOzHwMmKmTxb997FBJvW4rpCLrj5uSYd6fdb6ro/qc814OXHLtUVfre2LmC8JHZer9j5uud9cz6FYaOmhMpenY/XNBLme03mPW/UP6z75lHs7la000X+4P1R7dR7bVMsn9qrrizfVle1w7/vEUO38WvB7EH2hq28512JR7ffYcBmtffdX67zKiCmuERj28cIXzgP9Ja/2WdnW4762+ICxpsXGWaB3n51evq0u+sFD7VxHHaObkgWugzODKzHj5yRDnLhvtGxzfdlv6cGCgP7QGxtrv/RNPWSzVR9M9LXF3zovmcR07U/2IJO+1QaedTxrvdjZWaIPW9nxlteWuGADOrEVC38CT7rW8UPtxbLaCB/iRDJVRmLnsY8yflMued52m3x0oG1LT8759e1l+Jiyvm3LLtUVq1ZSQZ90V2Yt2bW9jvsWnrKhxq98qKcYqi+J3ypfOlQ5rleNC+SoDjFb/cP1RbnqPV7HD18b8D7Yr3rga9njdlK/tUUnru0+fruWqy/pU22THJ5QJHacketY47NPB7VxZtrv7H0kVT6xKN2q3+lfC3I+5VcS4Ze/3zwRoAW1P1mghTyf5iuJoK8e6NiTBWrDkw8qV3uOeeqBftkOJRVY2CsxQIKiJhUkw5MZkq3+/akL+muVTz6pIEfzKR+GeJAogLhBcV5bAoAAVBsfhKhbB94akAosDyAuptbF57pWOfQ3Zis5DA5en4uNQUTH6KH+uBCxWW1dlmz1xb7O+wAjGciu/YqBBuzqC+qJry+4JRdZMJN+TIBa7SSDgcOZq64vVtbHW+z1m4rr5TrRv7Yqhx06+3n2nbfKqo9dft8+7bDP/eNt+vyJLt6OMrZus+/rPL6QndqvTx0gw7euF+Vqq1jQdhvvFlfp5J9cSaZik6QUfbDlWpbOY17IJ33hiwz6UEKk+lH6MDb0MfH48+toiG8fG+kjPfy6crnoq0WqX+Men7W+2sh2r89YID37ri216/tTHCKjVcdZVnbORnqNjYVWPypzebVOPSddtHiAVyuGqwyOt9V1e9HJY8fH7aqT6ku2eEoOutEvW39yTj7266JPN4+F2q+OPUFLP32yOM9WMVDHbc6RGOkbC6su7hftOy/p04qRx77sHjO2rfr3saa2EzdPrI5lIp22XRvqp46tVbbaM065LPwlFpyHgcef+tg2vtDGtx4nlEsvHzPY13l09n5VruPWvEvnpL8vjlXXxxDJ5FjnsLHaQ9/I9EUf+gz5U239762vLbfX9WrFC7HJOOHHYuw86rFk+z2LvrheOdZWbf0aoaxPfrWhHrvsqoPb4PX69vvqS2fGXbVFB94vALM+ubV9PfZ2fn14X/pkv46LxGxLXpXjffhYoHLV5RpRbHAvk1yN5bpGvRxZ6z701YSerz4g269x9a3jqgMyfYs/sBubpIuPRbTpK2/xoQ1b+bSOMdgv/9I39dlWO9TGY5l6bCVH8vzDA51TOT7QPrpU+cjRVotsfbXhr90f10kAjn/48x+77vdPX1fQwv2HXyiJoOTC79ZPGOj4h9//5fEpA8kj+UAS4fH44WkG71v725IKWuzzRIQnEnwfGa2kgt4ZUZ8ykH21/CySClxQAOQi0HENci/jZtVXpnIPHB27bB1z4bosBbtu+AQc9XyCVOWojnSti3qV179WW+mpPn2wVBnHauP6YJd0ZQCp/bvMes51YsBert/2LiIv/2r/tJHVauFs/CL3dir3xAQydKPw9urd2+kYe6sMr0cd1fd6Kocj8cQxlrb4eDu30ffRDXnSx2PH9UBXLZyv/uHbF/78t4fvnSoq1UdrEFd/yMFOBlRsIx48NrCzbqstOq8yrknXX+e8b+dD3+LgdehPunCjoowtbfWuD5+scr5upRM3Xs45L7ebulzffs7b+gKhxlEfX8n2yVq1W32pX64rjzH2ec8FxzwCLd1qmfRijKjnqp/oeyMFS59vq4znZzdH7mPJZAxy3qrp9SSXT/Qrk1YflG2r6/YRn8STjuv4gcy6Harr8TEUOy19xcG/mlX7d046V/XReRi7zdv21+8seXhRlXRiLHJZVQ8/hifXiM65D3XsiwfXWXr1XR+wJFHotiNfL/nj2nNZtEXPVtvWWCPdfNxUO+xCb32vf2icoe7GQ2jxtK366Yzz1nkmqfU6c1+qneqio9tZ67Xi7UmjzV7tS/IYM9CR+HAb3Z617tdKeDy9+4B+qk4qd71a/dOf2+l9S4b373qqXp8/0Ymts1NZ7cP943Zs29/12pIsrgf0EhNPHOqYJwycieuhtq6vM9Y5v851rGtMj5bXxaYYcH2pno754KhPvvRj/EYPfKhj/qSv3/9UXv1P3b6t9GGu4XUku5VUUA/Sm7FfbXgKcFt7bOI6W7d7eEJN8iinnuxFN+ZT7tuqn+RJDnVdF+3LX/RBXcZ5+lGfqocML0deq1/OsXV7XIbaat5Ev7K1PpXqNqIn+lT7WKS3ynkJttus+prZEo/Eit6B49eH5ImFrg0fu2SL4rvaIWaq533BQtv1NaYkTUnEqB+xkFwfY6qvXJb214/8/0bvUvjj4yP/P3z79G4F1eEphM3TC398fBJgXe/h3Qqq95h8eEhIrL+a8Itvuj9IfuOvL6mg5IU/aeCJBO3zfgSJ9KSAd6Fyr0eiw8tJMuhJCCUgJJsnN6S7yklsuOxT7Te//kCgyMEKjL7gIHCoQ6C7sgwI1NHWLyANlJzzQYng9foEbq3f6oM62qKX9GVSgY7b2tYLwuvqAvNj9aMy/cEFezhWHT0Sis16W3VrkiAZ2IpM9NWWQQkb4enlPgCoHrZI55pkcf1gJf7s0/c2e9FBdYkblWGH61BtqD7Refc/52WfL1SQLc6alGpwVj/VdunhOsle9KUfHbfsc6aSTTuYaDvE1NtIj5ZPkee+UH/+z+OpykQ31xduLf3WjHVTtV+/qP6WrsSN5Pogj76+9b6lD/qqjp+TTL5LqH39wkKLibOo5/2cs3Dfql/3u46JGZdHDIgXL5TTvn5JoHV9Puvbrmf/JQJ08rrwlQ51kivd0EPs8J3K/c9jVDb4sfNWG2eOvBoLfX71tsSW66H92jfXYOt6FAfsdznOp+pPPdelxg4v9aIt/iZuJUN9cB6Z2hILsg/dWj7wei6XupR536tPmlTpxbub+FYf+kfdPu6un/NVW3yoOmrfNxbqrHMlFtWOcsqcLfLdXvXLtUJbWLba1vjCHlhJHqx1jr7goj78fKs9/XOOrfijK2Wuo7fzclhgn+voPpBs1aG++sDnKtO+n0MHbV22zwHqrzTAQ7qyCEAfGDof+sd3ri/tJJNxTHZ73KCXZHrfiiHn6XLFAV3ow23VPrL8vLdBX69HDKg9dSnDTh0PXVtVlxqT7nvpxz/3HVw8ZiTXObT8gCxt17p+vF+/7wSd3F6vS9z2yVd719t1hQ19brSkx40v8D/9099Trc0e59ENG71vMdFLkKlTfaRy2iG/tt9Y85yn2mGXbLrS2Pkwp6Nc8lxH4qhPPnFEPfShHF2dIb8spHNq578opGN96ODy5DMSUshvbb0P59+yx9t7TKx10q9SPbCq5+C0S7nY6SeUYS3fqf36ZbcP/PGx9HJ2fn08K7/Rr3m9/PUH7MJfzsH5yA4/ViIRXyHDt75gV3k9pu46+WC/3lCPVY9FOk8I1GPVUTt9lUGJA36BQot/ffWC9x/o/Bff/uvj8Ve/2PwyhRb80u8rnX/4tQqSAC5b73jgGJnU869UqN8vvv3101c77Fcwfvj2m2fvXxBT9+W6gyP/10wqHLmPiNuDQP094z1ENJvoYiYr3qyQwp0JnCNTDS7cgDC4Tpj9BTqq45+M0ub/b+9echqGgQCA3pEdB4Ath0LiKNwEiSUbkNihKbixokJpndid9C1QC4nt+E0+zsSl2V7XOq5OcYgLdv0U4JSya6zbI64xcJo/QVujL4fqfI6B/W6YdGjpuL/1cB/Xu/NaHrmfnLfFU6kYVM/PqdNS70YKRGzmN7ojt6dn23Gtr2/eerZ9bluxzVniFeO/+WyUc/ut3HUI/Jb02ELvJRUuKIplQBUXwCWzSbus48/XLhrILhPwzKax7TFbZj4Ajgv5zd39/ml9PRDJODApkV7ruCr1n/Ja7zenlFtr3R5xDf/ypGytfszrjQRZ7L/Rv/l+Pl93xO893Ef0q6XNEftJy/bWZa/5prV2uOT3ccxluVFdyjHOgxn7nClW5aNkS8VMPf8TeP/43M8CKE/ye7zePjw2tRu9q5MKPbY52nh5ffsfbONakgqNgEsWj4HJX1N8WtqKk3TUXd8ottSn7DT1NaNpmX4W+0T5KdNy66l9W4jzmsfVFny22Afnuy1GVZ8ItAmUZGNbLTlKx3XvEhOqx/QiRvsxydP315kfK9N7eXkwE9uZ0bi3l/YmgfJ/DuJGv3ycYVqa/52kQv4Y6gEBAgQIECBAgAABAgQIEBgiIKkwhF2jBAgQIECAAAECBAgQIEAgv4CkQv4Y6gEBAgQIECBAgAABAgQIEBgiIKkwhF2jBAgQIECAAAECBAgQIEAgv4CkQv4Y6gEBAgQIECBAgAABAgQIEBgiIKkwhF2jBAgQIECAAAECBAgQIEAgv4CkQv4Y6gEBAgQIECBAgAABAgQIEBgi8AVxZDa9oYxrjAAAAABJRU5ErkJggg=="
    }
   },
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![image.png](attachment:bec3c86a-4696-47e8-912b-4a920ed1d068.png)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "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>Shrt_Desc</th>\n",
       "      <th>BUTTER,WITH SALT</th>\n",
       "      <th>BUTTER,WHIPPED,WITH SALT</th>\n",
       "      <th>BUTTER OIL,ANHYDROUS</th>\n",
       "      <th>CHEESE,BLUE</th>\n",
       "      <th>CHEESE,BRICK</th>\n",
       "      <th>CHEESE,BRIE</th>\n",
       "      <th>CHEESE,CAMEMBERT</th>\n",
       "      <th>CHEESE,CARAWAY</th>\n",
       "      <th>CHEESE,CHEDDAR</th>\n",
       "      <th>CHEESE,CHESHIRE</th>\n",
       "      <th>...</th>\n",
       "      <th>WHEAT FLR,WHITE (INDUSTRIAL),10% PROT,BLEACHED,UNENR</th>\n",
       "      <th>WHEAT FLR,WHITE (INDUSTRIAL),10% PROT,UNBLEACHED,ENR</th>\n",
       "      <th>WHEAT FLR,WHITE (INDUSTRIAL),11.5% PROT,BLEACHED,ENR</th>\n",
       "      <th>WHEAT FLR,WHITE (INDUSTRIAL),11.5% PROT,BLEACHED,UNENR</th>\n",
       "      <th>WHEAT FLR,WHITE (INDUSTRIAL),11.5% PROT,UNBLEACHED,ENR</th>\n",
       "      <th>WHEAT FLR,WHITE (INDUSTRIAL),13% PROT,BLEACHED,ENR</th>\n",
       "      <th>WHEAT FLR,WHITE (INDUSTRIAL),13% PROT,BLEACHED,UNENR</th>\n",
       "      <th>WHEAT FLR,WHITE (INDUSTRIAL),15% PROT,BLEACHED,ENR</th>\n",
       "      <th>WHEAT FLR,WHITE (INDUSTRIAL),15% PROT,BLEACHED,UNENR</th>\n",
       "      <th>VITAL WHEAT GLUTEN</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>Energ_Kcal</th>\n",
       "      <td>717.000</td>\n",
       "      <td>717.000</td>\n",
       "      <td>876.000</td>\n",
       "      <td>353.000</td>\n",
       "      <td>371.000</td>\n",
       "      <td>334.000</td>\n",
       "      <td>300.000</td>\n",
       "      <td>376.000</td>\n",
       "      <td>403.000</td>\n",
       "      <td>387.000</td>\n",
       "      <td>...</td>\n",
       "      <td>366.000</td>\n",
       "      <td>366.000</td>\n",
       "      <td>363.000</td>\n",
       "      <td>363.000</td>\n",
       "      <td>363.000</td>\n",
       "      <td>362.000</td>\n",
       "      <td>362.000</td>\n",
       "      <td>362.000</td>\n",
       "      <td>362.000</td>\n",
       "      <td>370.000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Protein</th>\n",
       "      <td>0.850</td>\n",
       "      <td>0.850</td>\n",
       "      <td>0.280</td>\n",
       "      <td>21.400</td>\n",
       "      <td>23.240</td>\n",
       "      <td>20.750</td>\n",
       "      <td>19.800</td>\n",
       "      <td>25.180</td>\n",
       "      <td>24.900</td>\n",
       "      <td>23.370</td>\n",
       "      <td>...</td>\n",
       "      <td>9.710</td>\n",
       "      <td>9.710</td>\n",
       "      <td>11.500</td>\n",
       "      <td>11.500</td>\n",
       "      <td>11.500</td>\n",
       "      <td>13.070</td>\n",
       "      <td>13.070</td>\n",
       "      <td>15.330</td>\n",
       "      <td>15.330</td>\n",
       "      <td>75.160</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Carbohydrt</th>\n",
       "      <td>0.060</td>\n",
       "      <td>0.060</td>\n",
       "      <td>0.000</td>\n",
       "      <td>2.340</td>\n",
       "      <td>2.790</td>\n",
       "      <td>0.450</td>\n",
       "      <td>0.460</td>\n",
       "      <td>3.060</td>\n",
       "      <td>1.280</td>\n",
       "      <td>4.780</td>\n",
       "      <td>...</td>\n",
       "      <td>76.220</td>\n",
       "      <td>76.220</td>\n",
       "      <td>73.810</td>\n",
       "      <td>73.810</td>\n",
       "      <td>73.810</td>\n",
       "      <td>72.200</td>\n",
       "      <td>72.200</td>\n",
       "      <td>69.880</td>\n",
       "      <td>69.880</td>\n",
       "      <td>13.790</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Fiber_TD</th>\n",
       "      <td>0.000</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.000</td>\n",
       "      <td>...</td>\n",
       "      <td>2.400</td>\n",
       "      <td>2.400</td>\n",
       "      <td>2.400</td>\n",
       "      <td>2.400</td>\n",
       "      <td>2.400</td>\n",
       "      <td>2.400</td>\n",
       "      <td>2.400</td>\n",
       "      <td>2.400</td>\n",
       "      <td>2.400</td>\n",
       "      <td>0.600</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Sugar_Tot</th>\n",
       "      <td>0.060</td>\n",
       "      <td>0.060</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.500</td>\n",
       "      <td>0.510</td>\n",
       "      <td>0.450</td>\n",
       "      <td>0.460</td>\n",
       "      <td>0.000</td>\n",
       "      <td>0.520</td>\n",
       "      <td>0.000</td>\n",
       "      <td>...</td>\n",
       "      <td>0.490</td>\n",
       "      <td>0.490</td>\n",
       "      <td>1.120</td>\n",
       "      <td>1.120</td>\n",
       "      <td>1.120</td>\n",
       "      <td>1.100</td>\n",
       "      <td>1.100</td>\n",
       "      <td>0.920</td>\n",
       "      <td>0.920</td>\n",
       "      <td>0.000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Lipid_Tot</th>\n",
       "      <td>81.110</td>\n",
       "      <td>81.110</td>\n",
       "      <td>99.480</td>\n",
       "      <td>28.740</td>\n",
       "      <td>29.680</td>\n",
       "      <td>27.680</td>\n",
       "      <td>24.260</td>\n",
       "      <td>29.200</td>\n",
       "      <td>33.140</td>\n",
       "      <td>30.600</td>\n",
       "      <td>...</td>\n",
       "      <td>1.480</td>\n",
       "      <td>1.480</td>\n",
       "      <td>1.450</td>\n",
       "      <td>1.450</td>\n",
       "      <td>1.450</td>\n",
       "      <td>1.380</td>\n",
       "      <td>1.380</td>\n",
       "      <td>1.410</td>\n",
       "      <td>1.410</td>\n",
       "      <td>1.850</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>FA_Sat</th>\n",
       "      <td>51.368</td>\n",
       "      <td>50.489</td>\n",
       "      <td>61.924</td>\n",
       "      <td>18.669</td>\n",
       "      <td>18.764</td>\n",
       "      <td>17.410</td>\n",
       "      <td>15.259</td>\n",
       "      <td>18.584</td>\n",
       "      <td>21.092</td>\n",
       "      <td>19.475</td>\n",
       "      <td>...</td>\n",
       "      <td>0.302</td>\n",
       "      <td>0.302</td>\n",
       "      <td>0.268</td>\n",
       "      <td>0.268</td>\n",
       "      <td>0.268</td>\n",
       "      <td>0.189</td>\n",
       "      <td>0.189</td>\n",
       "      <td>0.272</td>\n",
       "      <td>0.272</td>\n",
       "      <td>0.272</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>FA_Poly</th>\n",
       "      <td>3.043</td>\n",
       "      <td>3.012</td>\n",
       "      <td>3.694</td>\n",
       "      <td>0.800</td>\n",
       "      <td>0.784</td>\n",
       "      <td>0.826</td>\n",
       "      <td>0.724</td>\n",
       "      <td>0.830</td>\n",
       "      <td>0.942</td>\n",
       "      <td>0.870</td>\n",
       "      <td>...</td>\n",
       "      <td>0.845</td>\n",
       "      <td>0.845</td>\n",
       "      <td>0.749</td>\n",
       "      <td>0.749</td>\n",
       "      <td>0.749</td>\n",
       "      <td>0.683</td>\n",
       "      <td>0.683</td>\n",
       "      <td>0.712</td>\n",
       "      <td>0.712</td>\n",
       "      <td>0.810</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Sodium</th>\n",
       "      <td>0.576</td>\n",
       "      <td>0.827</td>\n",
       "      <td>0.002</td>\n",
       "      <td>1.395</td>\n",
       "      <td>0.560</td>\n",
       "      <td>0.629</td>\n",
       "      <td>0.842</td>\n",
       "      <td>0.690</td>\n",
       "      <td>0.621</td>\n",
       "      <td>0.700</td>\n",
       "      <td>...</td>\n",
       "      <td>0.002</td>\n",
       "      <td>0.002</td>\n",
       "      <td>0.002</td>\n",
       "      <td>0.002</td>\n",
       "      <td>0.002</td>\n",
       "      <td>0.002</td>\n",
       "      <td>0.002</td>\n",
       "      <td>0.002</td>\n",
       "      <td>0.002</td>\n",
       "      <td>0.029</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>9 rows × 2269 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "Shrt_Desc   BUTTER,WITH SALT  BUTTER,WHIPPED,WITH SALT  BUTTER OIL,ANHYDROUS  \\\n",
       "Energ_Kcal           717.000                   717.000               876.000   \n",
       "Protein                0.850                     0.850                 0.280   \n",
       "Carbohydrt             0.060                     0.060                 0.000   \n",
       "Fiber_TD               0.000                     0.000                 0.000   \n",
       "Sugar_Tot              0.060                     0.060                 0.000   \n",
       "Lipid_Tot             81.110                    81.110                99.480   \n",
       "FA_Sat                51.368                    50.489                61.924   \n",
       "FA_Poly                3.043                     3.012                 3.694   \n",
       "Sodium                 0.576                     0.827                 0.002   \n",
       "\n",
       "Shrt_Desc   CHEESE,BLUE  CHEESE,BRICK  CHEESE,BRIE  CHEESE,CAMEMBERT  \\\n",
       "Energ_Kcal      353.000       371.000      334.000           300.000   \n",
       "Protein          21.400        23.240       20.750            19.800   \n",
       "Carbohydrt        2.340         2.790        0.450             0.460   \n",
       "Fiber_TD          0.000         0.000        0.000             0.000   \n",
       "Sugar_Tot         0.500         0.510        0.450             0.460   \n",
       "Lipid_Tot        28.740        29.680       27.680            24.260   \n",
       "FA_Sat           18.669        18.764       17.410            15.259   \n",
       "FA_Poly           0.800         0.784        0.826             0.724   \n",
       "Sodium            1.395         0.560        0.629             0.842   \n",
       "\n",
       "Shrt_Desc   CHEESE,CARAWAY  CHEESE,CHEDDAR  CHEESE,CHESHIRE  ...  \\\n",
       "Energ_Kcal         376.000         403.000          387.000  ...   \n",
       "Protein             25.180          24.900           23.370  ...   \n",
       "Carbohydrt           3.060           1.280            4.780  ...   \n",
       "Fiber_TD             0.000           0.000            0.000  ...   \n",
       "Sugar_Tot            0.000           0.520            0.000  ...   \n",
       "Lipid_Tot           29.200          33.140           30.600  ...   \n",
       "FA_Sat              18.584          21.092           19.475  ...   \n",
       "FA_Poly              0.830           0.942            0.870  ...   \n",
       "Sodium               0.690           0.621            0.700  ...   \n",
       "\n",
       "Shrt_Desc   WHEAT FLR,WHITE (INDUSTRIAL),10% PROT,BLEACHED,UNENR  \\\n",
       "Energ_Kcal                                            366.000      \n",
       "Protein                                                 9.710      \n",
       "Carbohydrt                                             76.220      \n",
       "Fiber_TD                                                2.400      \n",
       "Sugar_Tot                                               0.490      \n",
       "Lipid_Tot                                               1.480      \n",
       "FA_Sat                                                  0.302      \n",
       "FA_Poly                                                 0.845      \n",
       "Sodium                                                  0.002      \n",
       "\n",
       "Shrt_Desc   WHEAT FLR,WHITE (INDUSTRIAL),10% PROT,UNBLEACHED,ENR  \\\n",
       "Energ_Kcal                                            366.000      \n",
       "Protein                                                 9.710      \n",
       "Carbohydrt                                             76.220      \n",
       "Fiber_TD                                                2.400      \n",
       "Sugar_Tot                                               0.490      \n",
       "Lipid_Tot                                               1.480      \n",
       "FA_Sat                                                  0.302      \n",
       "FA_Poly                                                 0.845      \n",
       "Sodium                                                  0.002      \n",
       "\n",
       "Shrt_Desc   WHEAT FLR,WHITE (INDUSTRIAL),11.5% PROT,BLEACHED,ENR  \\\n",
       "Energ_Kcal                                            363.000      \n",
       "Protein                                                11.500      \n",
       "Carbohydrt                                             73.810      \n",
       "Fiber_TD                                                2.400      \n",
       "Sugar_Tot                                               1.120      \n",
       "Lipid_Tot                                               1.450      \n",
       "FA_Sat                                                  0.268      \n",
       "FA_Poly                                                 0.749      \n",
       "Sodium                                                  0.002      \n",
       "\n",
       "Shrt_Desc   WHEAT FLR,WHITE (INDUSTRIAL),11.5% PROT,BLEACHED,UNENR  \\\n",
       "Energ_Kcal                                            363.000        \n",
       "Protein                                                11.500        \n",
       "Carbohydrt                                             73.810        \n",
       "Fiber_TD                                                2.400        \n",
       "Sugar_Tot                                               1.120        \n",
       "Lipid_Tot                                               1.450        \n",
       "FA_Sat                                                  0.268        \n",
       "FA_Poly                                                 0.749        \n",
       "Sodium                                                  0.002        \n",
       "\n",
       "Shrt_Desc   WHEAT FLR,WHITE (INDUSTRIAL),11.5% PROT,UNBLEACHED,ENR  \\\n",
       "Energ_Kcal                                            363.000        \n",
       "Protein                                                11.500        \n",
       "Carbohydrt                                             73.810        \n",
       "Fiber_TD                                                2.400        \n",
       "Sugar_Tot                                               1.120        \n",
       "Lipid_Tot                                               1.450        \n",
       "FA_Sat                                                  0.268        \n",
       "FA_Poly                                                 0.749        \n",
       "Sodium                                                  0.002        \n",
       "\n",
       "Shrt_Desc   WHEAT FLR,WHITE (INDUSTRIAL),13% PROT,BLEACHED,ENR  \\\n",
       "Energ_Kcal                                            362.000    \n",
       "Protein                                                13.070    \n",
       "Carbohydrt                                             72.200    \n",
       "Fiber_TD                                                2.400    \n",
       "Sugar_Tot                                               1.100    \n",
       "Lipid_Tot                                               1.380    \n",
       "FA_Sat                                                  0.189    \n",
       "FA_Poly                                                 0.683    \n",
       "Sodium                                                  0.002    \n",
       "\n",
       "Shrt_Desc   WHEAT FLR,WHITE (INDUSTRIAL),13% PROT,BLEACHED,UNENR  \\\n",
       "Energ_Kcal                                            362.000      \n",
       "Protein                                                13.070      \n",
       "Carbohydrt                                             72.200      \n",
       "Fiber_TD                                                2.400      \n",
       "Sugar_Tot                                               1.100      \n",
       "Lipid_Tot                                               1.380      \n",
       "FA_Sat                                                  0.189      \n",
       "FA_Poly                                                 0.683      \n",
       "Sodium                                                  0.002      \n",
       "\n",
       "Shrt_Desc   WHEAT FLR,WHITE (INDUSTRIAL),15% PROT,BLEACHED,ENR  \\\n",
       "Energ_Kcal                                            362.000    \n",
       "Protein                                                15.330    \n",
       "Carbohydrt                                             69.880    \n",
       "Fiber_TD                                                2.400    \n",
       "Sugar_Tot                                               0.920    \n",
       "Lipid_Tot                                               1.410    \n",
       "FA_Sat                                                  0.272    \n",
       "FA_Poly                                                 0.712    \n",
       "Sodium                                                  0.002    \n",
       "\n",
       "Shrt_Desc   WHEAT FLR,WHITE (INDUSTRIAL),15% PROT,BLEACHED,UNENR  \\\n",
       "Energ_Kcal                                            362.000      \n",
       "Protein                                                15.330      \n",
       "Carbohydrt                                             69.880      \n",
       "Fiber_TD                                                2.400      \n",
       "Sugar_Tot                                               0.920      \n",
       "Lipid_Tot                                               1.410      \n",
       "FA_Sat                                                  0.272      \n",
       "FA_Poly                                                 0.712      \n",
       "Sodium                                                  0.002      \n",
       "\n",
       "Shrt_Desc   VITAL WHEAT GLUTEN  \n",
       "Energ_Kcal             370.000  \n",
       "Protein                 75.160  \n",
       "Carbohydrt              13.790  \n",
       "Fiber_TD                 0.600  \n",
       "Sugar_Tot                0.000  \n",
       "Lipid_Tot                1.850  \n",
       "FA_Sat                   0.272  \n",
       "FA_Poly                  0.810  \n",
       "Sodium                   0.029  \n",
       "\n",
       "[9 rows x 2269 columns]"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import pandas as pd\n",
    "\n",
    "ingredients = pd.read_excel(\"Matrix Porridge (filtered).xlsx\", \"core ingredients\")\n",
    "ingredients['Sodium']=ingredients['Sodium']/1000 #values are in mg... convert to g\n",
    "inputs=ingredients[['Energ_Kcal','Protein','Carbohydrt','Fiber_TD','Sugar_Tot','Lipid_Tot','FA_Sat','FA_Poly','Sodium']].T\n",
    "inputs.columns=ingredients['Shrt_Desc']\n",
    "# inputs.dropna(inplace=True,axis=1) #some of these ingredients have null values. Remove those. \n",
    "inputs.fillna(0,inplace=True)\n",
    "inputs.head(10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Index(['NDB_No', 'Shrt_Desc', 'Water', 'Energ_Kcal', 'Protein', 'Lipid_Tot',\n",
      "       'Ash', 'Carbohydrt', 'Fiber_TD', 'Sugar_Tot', 'Calcium', 'Iron',\n",
      "       'Magnesium', 'Phosphorus', 'Potassium', 'Sodium', 'Zinc', 'Copper',\n",
      "       'Manganese', 'Selenium', 'Vit_C', 'Thiamin', 'Riboflavin', 'Niacin',\n",
      "       'Panto_Acid', 'Vit_B6', 'Folate_Tot', 'Folic_Acid', 'Food_Folate',\n",
      "       'Folate_DFE', 'Choline_Tot', 'Vit_B12', 'Vit_A_IU', 'Vit_A_RAE',\n",
      "       'Retinol', 'Alpha_Carot', 'Beta_Carot', 'Beta_Crypt', 'Lycopene',\n",
      "       'Lut+Zea', 'Vit_E', 'Vit_D_mcg', 'ViVit_D_IU', 'Vit_K', 'FA_Sat',\n",
      "       'FA_Mono', 'FA_Poly', 'Cholestrl', 'GmWt_1', 'GmWt_Desc1', 'GmWt_2',\n",
      "       'GmWt_Desc2', 'Refuse_Pct', 'Expr1053', 'FdGrp_CD'],\n",
      "      dtype='object')\n"
     ]
    }
   ],
   "source": [
    "print(ingredients.columns)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                                          Unnamed: 0    min (g)     max (g)\n",
      "0                                                Fat  14.814815   25.925926\n",
      "1   n-6 polyunsaturated fatty acidsa (linoleic acid)   3.703704    7.407407\n",
      "2  n-3 polyunsaturated fatty acidsa (α-linolenic ...   0.444444    0.888889\n",
      "3                                       Carbohydrate  75.000000  108.333333\n",
      "4                                            Protein  16.666667   58.333333\n",
      "5                                             Sodium   0.500000    0.766667\n"
     ]
    }
   ],
   "source": [
    "requirements=pd.read_excel(\"Matrix Porridge (filtered).xlsx\", \"requirements (2021)\")\n",
    "requirements['min (g)']/=3\n",
    "requirements['max (g)']/=3\n",
    "\n",
    "print(requirements)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note: our data does not split the two kinds of polyunsaturated fat, so we'll have to sum these up. There's another dataset I've been looking at which might work better for this, but we'll get to that later."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "fat_min=requirements.iloc[0,1]\n",
    "fat_max=requirements.iloc[0,2]\n",
    "fat_half_range=(fat_max-fat_min)/2 #calculate this once so we don't need to do it repeatedly in our function later\n",
    "fat_opt=(fat_min+fat_max)/2\n",
    "\n",
    "pufat_min=requirements.iloc[1,1]+requirements.iloc[2,2]\n",
    "pufat_max=requirements.iloc[1,2]+requirements.iloc[2,2]\n",
    "pufat_half_range=(pufat_max-pufat_min)/2\n",
    "pufat_opt=(pufat_min+pufat_max)/2\n",
    "\n",
    "carb_min=requirements.iloc[3,1]\n",
    "carb_max=requirements.iloc[3,2]\n",
    "carb_half_range=(carb_max-carb_min)/2\n",
    "carb_opt=(carb_min+carb_max)/2\n",
    "\n",
    "protein_min=requirements.iloc[4,1]\n",
    "protein_max=requirements.iloc[4,2]\n",
    "protein_half_range=(protein_max-protein_min)/2\n",
    "protein_opt=(protein_min+protein_max)/2\n",
    "\n",
    "sodium_min=requirements.iloc[5,1]\n",
    "sodium_max=requirements.iloc[5,2]\n",
    "sodium_half_range=(sodium_max-sodium_min)/2\n",
    "sodium_opt=(sodium_min+sodium_max)/2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The Approach\n",
    "This is clearly an optimization problem. However, it is a bit more complicated than what you'd typically use Linear Programming to solve. We could potentially frame it that way... with an A matrix 2330 columns wide. But our objective function here isn't written easily as a function of pure X (our ingredients vector), no matter what we decide to optimize.\n",
    "\n",
    "There are, however, many modern machine learning approaches that can help us here. The trick is to use a solver / algorithm that works with a custom function instead of a vector. This way, we can to optimize for a custom value function. SciPy's Linear Programming functionality requires that our objective be a vector, but that doesn't really work here. But there are plenty of alternatives, so we'll try a few of those. \n",
    "\n",
    "First, however, we need to define what it is we seek to optimize. Eventually, we may wish to make this something the user could select from a list of options (which would also alter constraints), to accommodate different nutrition guides, such as Atkins, or a high fiber diet. Or, maybe we seek to maximize quantity of food while still meeting nutrition guidelines. For now, I'm going to take bit of a fuzzy logic approach, with our function outputting a value that's most optimal when all nutrition requirements are exactly in the center of the ranges and where we impose serious (but linear with a slope) penalties if any nutritional requirements fall outisde the acceptable ranges."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[8.3000e+02 5.1940e+01 2.6800e+01 6.0000e+00 5.6400e+00 5.8120e+01\n",
      " 3.6592e+01 2.0980e+00 7.6300e-01]\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "x=np.zeros(2269)\n",
    "x[4]=1\n",
    "x[39]=1\n",
    "x[554]=1\n",
    "x[1000]=1\n",
    "A=inputs.values #get A matrix\n",
    "y=A.dot(x)\n",
    "print(y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "def diet_function(ingredients_vector):\n",
    "    y=A.dot(ingredients_vector)\n",
    "    calorie_penalty=10*abs(666.7-y[0]) #weighting calories very highly here\n",
    "    \n",
    "    protein_penalty=abs(y[1]-protein_opt)\n",
    "    if abs(y[1]-protein_opt)>protein_half_range:\n",
    "        protein_penalty*=20\n",
    "    \n",
    "    carb_penalty=abs((y[2]+y[3])-carb_opt)\n",
    "    if carb_penalty>carb_half_range:\n",
    "        carb_penalty*=20    \n",
    "    fiber_bonus=y[3]\n",
    "    sugar_penalty=y[4]\n",
    "\n",
    "        \n",
    "    fat_penalty=abs(y[5]-fat_opt)\n",
    "    if fat_penalty>fat_half_range:\n",
    "        fat_penalty*=20\n",
    "       \n",
    "    sat_fat_penalty=y[6]\n",
    "    \n",
    "    pufat_penalty=abs(y[7]-pufat_opt)\n",
    "    if pufat_penalty>pufat_half_range:\n",
    "        pufat_penalty*=20\n",
    "           \n",
    "    sodium_penalty=abs(y[8]-sodium_opt)\n",
    "    if sodium_penalty>sodium_half_range:\n",
    "        sodium_penalty*=100            \n",
    "            \n",
    "    value=-calorie_penalty-protein_penalty-carb_penalty+fiber_bonus-sugar_penalty-fat_penalty-sat_fat_penalty-pufat_penalty-sodium_penalty\n",
    "    return -value #return negative since this is a minimization problem\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy.optimize import minimize,dual_annealing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "bounds=tuple([(0,10) for i in range(len(ingredients))])\n",
    "x0=[1]*len(ingredients)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "res=minimize(diet_function,x0,method='TNC',bounds=bounds,tol=1e-6,options={'maxiter':int(1e9),'minfev':-100})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Local vs Global Optimization\n",
    "SciPy's minimize function is a local optimization algorithm, with many different methods you can choose from to find a local minima based on various methods depending on whether you have bounded inputs or other constraints. \n",
    "\n",
    "For those that unfamiliar with local vs global optimization concepts, Mathworks (the makers of Matlab) explains it quite well [here](https://www.mathworks.com/help/gads/what-is-global-optimization.html).\n",
    "\n",
    "In a nutshell, because of the way we set up our value function with different slopes / contributions to the overall value for each nutritional category depending on whether they are inside our outside our acceptable ranges, we have made this into a nonlinear problem. What that means is that if we start at a random initial position on our value function and use something like gradient descent or Newton's method (or various other approaches) to follow the slope to the local minima, we can't be sure that this is the same as the overall or global minima, which is the true optimal solution. Our starting point and various other hyperparameters (like learning rate) matter. Check out the following graph for a visual explanation:\n",
    "\n",
    "<img src=\"https://www.mathworks.com/help/gads/local_vs_global.png\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So the way that we typically go about finding the true optimal solution for these sorts of problems is to use a local optimizer with a multitude of initial starting conditions. Those starting conditions can be purely random, or they can follow some sort of search logic. SciPy offers several options. Generally, I find the dual_annealing offers a good overall performance here. But, I'm going to hold off on running that for now because it's quite slow, and as you'll see later, we have a lot of work left to do before we're ready to go for the final run."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wall time: 0 ns\n"
     ]
    }
   ],
   "source": [
    "# %%time\n",
    "# res=dual_annealing(diet_function,x0=x0,bounds=bounds,maxfun=1e8,local_search_options={'method':'TNC','options':{'maxiter':int(1e9),'minfev':-100}})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "     fun: 2109.4068206642883\n",
      "     jac: array([5569.21813768, 5581.50914003, 6758.8240654 , ..., 5037.38206135,\n",
      "       5037.38206135, 4010.88204853])\n",
      " message: 'Converged (|x_n-x_(n-1)| ~= 0)'\n",
      "    nfev: 14307810\n",
      "     nit: 2117\n",
      "  status: 2\n",
      " success: True\n",
      "       x: array([0., 0., 0., ..., 0., 0., 0.])\n"
     ]
    }
   ],
   "source": [
    "print(res)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "solution=res.x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Energ_Kcal    666.700000\n",
      "Protein        48.623464\n",
      "Carbohydrt    132.148711\n",
      "Fiber_TD       43.081574\n",
      "Sugar_Tot      47.571787\n",
      "Lipid_Tot       7.716130\n",
      "FA_Sat          1.469776\n",
      "FA_Poly         2.997426\n",
      "Sodium          4.613850\n",
      "dtype: float64\n"
     ]
    }
   ],
   "source": [
    "output=pd.Series(A.dot(solution),index=inputs.index)\n",
    "print(output)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note: this solution looks pretty good, but it's a bit high in both carbs and sodium. Sugar content is high too, thought it's offset by a very nice fiber content. On the sugar side, they do appear to be natural sugars from vegetables, at least. We may want to revise our value function a bit, but this is not a bad initial result."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Shrt_Desc\n",
      "ACEROLA JUICE,RAW                                            0.418687\n",
      "BUTTERBUR,CANNED                                             0.859419\n",
      "CABBAGE,CHINESE (PAK-CHOI),RAW                               0.594643\n",
      "CABBAGE,CHINESE (PAK-CHOI),CKD,BLD,DRND,WO/SALT              0.507121\n",
      "CUCUMBER,WITH PEEL,RAW                                       0.439779\n",
      "CUCUMBER,PEELED,RAW                                          0.556533\n",
      "GOURD,WHITE-FLOWERED (CALABASH),CKD,BLD,DRND,WO/SALT         0.541778\n",
      "LETTUCE,BUTTERHEAD (INCL BOSTON&BIBB TYPES),RAW              0.640030\n",
      "PUMPKIN FLOWERS,CKD,BLD,DRND,WO/SALT                         0.447101\n",
      "PUMPKIN,CKD,BLD,DRND,WO/SALT                                 0.551579\n",
      "SQUASH,SMMR,CROOKNECK&STRAIGHTNECK,RAW                       0.430714\n",
      "SQUASH,SMMR,CROOKNECK&STRAIGHTNECK,CND,DRND,SOLID,WO/SALT    0.408181\n",
      "SQUASH,SUMMER,SCALLOP,RAW                                    0.434071\n",
      "SQUASH,SMMR,SCALLOP,CKD,BLD,DRND,WO/SALT                     0.449769\n",
      "SQUASH,SMMR,ZUCCHINI,INCL SKN,RAW                            0.644876\n",
      "SQUASH,SMMR,ZUCCHINI,INCL SKN,CKD,BLD,DRND,WO/SALT           0.595942\n",
      "SWAMP CABBAGE,CKD,BLD,DRND,WO/SALT                           0.471019\n",
      "TARO SHOOTS,RAW                                              0.724743\n",
      "TURNIP GRNS,CND,SOL&LIQUIDS                                  0.435474\n",
      "WATERCRESS,RAW                                               0.627250\n",
      "RADISHES,WHITE ICICLE,RAW                                    0.469647\n",
      "SQUASH,SMMR,ALL VAR,RAW                                      0.467083\n",
      "ASPARAGUS,CND,NO SALT,SOL&LIQUIDS                            0.614277\n",
      "BEANS,SNAP,GRN,CND,NO SALT,DRND SOL                          0.528400\n",
      "BUTTERBUR,CKD,BLD,DRND,W/SALT                                0.577284\n",
      "NEW ZEALAND SPINACH,CKD,BLD,DRND,W/SALT                      0.471180\n",
      "PUMPKIN,CKD,BLD,DRND,W/SALT                                  0.557505\n",
      "SQUASH,SMMR,SCALLOP,CKD,BLD,DRND,W/SALT                      0.513966\n",
      "WAXGOURD,(CHINESE PRESERVING MELON),CKD,BLD,DRND,W/SALT      0.653377\n",
      "PEPPERS,SWT,RED,FRZ,CHOPD,CKD,BLD,DRND,WO/SALT               0.453757\n",
      "PEPPERS,SWT,RED,FRZ,CHOPD,CKD,BLD,DRND,W/SALT                0.431042\n",
      "PICKLES,CUCUMBER,DILL OR KOSHER DILL                         0.477346\n",
      "PICKLES,CUCUMBER,SOUR                                        0.483702\n",
      "PICKLES,CUCUMBER,SOUR,LO NA                                  0.611388\n",
      "SQUASH,ZUCCHINI,BABY,RAW                                     0.525956\n",
      "CABBAGE,NAPA,COOKED                                          0.565690\n",
      "COCONUT H2O (LIQ FROM COCONUTS)                              0.411036\n",
      "dtype: float64\n"
     ]
    }
   ],
   "source": [
    "solution_ds=pd.Series(solution,index=ingredients.Shrt_Desc)\n",
    "solution_ds.to_csv('matrix_porridge_init_recipe.csv')\n",
    "solution_ds=solution_ds[solution_ds>0.4]\n",
    "print(solution_ds)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "An interesting note here is that we have clusters of similar ingredients. We could probably greatly reduce the work of the algorithm by first using an unsupervised clustering algorithm on the input data to trim our ingredients list down. We could then add some complexity on the value side by adding in more dietary considerations. We'll try that next."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "jokia",
   "language": "python",
   "name": "jokia"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}