{ "cells": [ { "cell_type": "code", "execution_count": 1, "id": "5bf19128-f39c-43ce-8503-9eff8df18434", "metadata": { "slideshow": { "slide_type": "" }, "tags": [ "hide-cell" ] }, "outputs": [], "source": [ "# Install the necessary dependencies\n", "\n", "import os\n", "import sys\n", "!{sys.executable} -m pip install --quiet pandas scikit-learn numpy matplotlib jupyterlab_myst ipython\n" ] }, { "cell_type": "markdown", "id": "cb578cb6-9845-4eb3-94ab-51c5e677031a", "metadata": { "slideshow": { "slide_type": "skip" }, "tags": [ "remove-cell" ] }, "source": [ "---\n", "license:\n", " code: MIT\n", " content: CC-BY-4.0\n", "github: https://github.com/ocademy-ai/machine-learning\n", "venue: By Ocademy\n", "open_access: true\n", "bibliography:\n", " - https://raw.githubusercontent.com/ocademy-ai/machine-learning/main/open-machine-learning-jupyter-book/references.bib\n", "---" ] }, { "cell_type": "markdown", "id": "663dea79", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "# Linear and polynomial regression\n", "\n", ":::{figure} https://static-1300131294.cos.ap-shanghai.myqcloud.com/images/ml-regression/linear-polynomial.png\n", "---\n", "name: 'Linear vs polynomial regression infographic'\n", "width: 100%\n", "---\n", "Infographic by [Dasani Madipalli](https://twitter.com/dasani_decoded)\n", ":::" ] }, { "cell_type": "markdown", "id": "ed080429-60c9-4f56-a3b5-7f3c533ec961", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "```{seealso}\n", "Run the notebook accompanying this lesson and look at the Month to Price scatterplot. Does the data associating Month to Price for pumpkin sales seem to have high or low correlation, according to our visual interpretation of the scatterplot? Does that change if we use more fine-grained measure instead of `Month`, eg. *day of the year* (i.e. number of days since the beginning of the year)?\n", "```" ] }, { "cell_type": "markdown", "id": "d0525b63-aae3-492a-b251-8131622648d0", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "## Build a regression model using Scikit-learn: regression four ways" ] }, { "cell_type": "code", "execution_count": 2, "id": "f1dfeda2", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "text/html": [ "\n", "

\n", "\n", "A demo of linear-regression. [source]\n", "

\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "from IPython.display import HTML\n", "\n", "display(\n", " HTML(\n", " \"\"\"\n", "

\n", "\n", "A demo of linear-regression. [source]\n", "

\n", "\"\"\"\n", " )\n", ")" ] }, { "cell_type": "code", "execution_count": 2, "id": "28b44e74", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "text/html": [ "\n", "
\n", "

\n", " \n", " A demo of gradient_react_3D. [source]\n", "

\n", "
\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "from IPython.display import HTML\n", "\n", "display(\n", " HTML(\n", " \"\"\"\n", "
\n", "

\n", " \n", " A demo of gradient_react_3D. [source]\n", "

\n", "
\n", "\"\"\"\n", " )\n", ")" ] }, { "cell_type": "code", "execution_count": 4, "id": "a14c82b5", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "text/html": [ "\n", "
\n", "

\n", " \n", " A demo of gradient-descent-visualiser. [source]\n", "

\n", "
\n" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "from IPython.display import HTML\n", "\n", "display(\n", " HTML(\n", " \"\"\"\n", "
\n", "

\n", " \n", " A demo of gradient-descent-visualiser. [source]\n", "

\n", "
\n", "\"\"\"\n", " )\n", ")" ] }, { "cell_type": "markdown", "id": "d08017d2", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "\n", "\n", "### Introduction\n", "\n", "So far we have explored what regression is with sample data gathered from the pumpkin pricing dataset that we will use throughout this lesson. We have also visualized it using Matplotlib.\n", "\n", "Now we are ready to dive deeper into regression for Machine Learning. While visualization allows we to make sense of data, the real power of Machine Learning comes from _training models_. Models are trained on historical data to automatically capture data dependencies, and they allow we to predict outcomes for new data, which the model has not seen before.\n", "\n", "In this lesson, we will learn more about two types of regression: _basic linear regression_ and _polynomial regression_, along with some of the math underlying these techniques. Those models will allow us to predict pumpkin prices depending on different input data.\n", "\n", ":::{note}\n", "Throughout this curriculum, we assume minimal knowledge of math, and seek to make it accessible for students coming from other fields, so watch for notes, callouts, diagrams, and other learning tools to aid in comprehension.\n", ":::\n", "\n", "### Prerequisite\n", "\n", "We should be familiar by now with the structure of the pumpkin data that we are examining. We can find it preloaded and pre-cleaned in this section's [linear and polynomial regression.ipynb](../../assignments/ml-fundamentals/linear-and-polynomial-regression.ipynb) file. In the file, the pumpkin price is displayed per bushel in a new data frame. Make sure we can run these notebooks in kernels in Visual Studio Code.\n", "\n", "### Preparation\n", "\n", "As a reminder, we are loading this data so as to ask questions about it.\n", "\n", "- When is the best time to buy pumpkins? \n", "- What price can I expect of a case of miniature pumpkins?\n", "- Should I buy them in half-bushel baskets or by the 1 1/9 bushel box?\n", "Let's keep digging into this data.\n", "\n", "In the previous lesson, we created a Pandas data frame and populated it with part of the original dataset, standardizing the pricing by bushel. By doing that, however, we were only able to gather about 400 data points and only for the fall months.\n", "\n", "Take a look at the data that we preloaded in this lesson's accompanying notebook. The data is preloaded and an initial scatterplot is charted to show monthly data. Maybe we can get a little more detail about the nature of the data by cleaning it more.\n", "\n", "## A linear regression line\n", "\n", "As we learned in section 1, the goal of a linear regression exercise is to be able to plot a line to:\n", "\n", "- **Show variable relationships**. Show the relationship between variables\n", "- **Make predictions**. Make accurate predictions on where a new data point would fall in relationship to that line.\n", "\n", "It is typical of **Least-Squares Regression** to draw this type of line. The term 'least-squares' means that all the data points surrounding the regression line are squared and then added up. Ideally, that final sum is as small as possible, because we want a low number of errors or `least-squares``.\n", "\n", "We do so since we want to model a line that has the least cumulative distance from all of our data points. We also square the terms before adding them since we are concerned with their magnitude rather than their direction.\n", "\n", ":::{seealso}\n", "**Show me the math** \n", " \n", "This line, called the _line of best fit_ can be expressed by [an equation](https://en.wikipedia.org/wiki/Simple_linear_regression): \n", "\n", "> ```\n", "> Y = a + bX\n", "> ```\n", "\n", "`X` is the 'explanatory variable'. `Y` is the 'dependent variable'. The slope of the line is `b` and `a` is the y-intercept, which refers to the value of `Y` when `X = 0`. \n", "\n", ">:::{figure} https://static-1300131294.cos.ap-shanghai.myqcloud.com/images/ml-regression/slope.png\n", ">---\n", ">name: 'calculate the slope'\n", ">width: 60%\n", ">---\n", ">Infographic by [Jen Looper](https://twitter.com/jenlooper)\n", ">:::\n", "\n", "First, calculate the slope `b`.\n", "\n", "In other words, and referring to our pumpkin data's original question: \"predict the price of a pumpkin per bushel by month\", `X` would refer to the price and `Y` would refer to the month of sale.\n", "\n", ">:::{figure} https://static-1300131294.cos.ap-shanghai.myqcloud.com/images/ml-regression/calculation.png\n", ">---\n", ">name: 'complete the equation'\n", ">width: 60%\n", ">---\n", ">Infographic by [Jen Looper](https://twitter.com/jenlooper)\n", ">:::\n", "\n", "Calculate the value of Y. If we're paying around $4, it must be April!\n", "\n", "The math that calculates the line must demonstrate the slope of the line, which is also dependent on the intercept, or where `Y` is situated when `X = 0`.\n", "\n", "We can observe the method of calculation for these values on the [Math is Fun](https://www.mathsisfun.com/data/least-squares-regression.html) website. Also, visit [this Least-squares calculator](https://www.mathsisfun.com/data/least-squares-calculator.html) to watch how the numbers' values impact the line.\n", ":::\n", "\n", "## Correlation\n", "\n", "One more term to understand is the **Correlation Coefficient** between given X and Y variables. Using a scatterplot, we can quickly visualize this coefficient. A plot with data points scattered in a neat line has a high correlation, but a plot with data points scattered everywhere between X and Y has a low correlation.\n", "\n", "A good linear regression model will be one that has a high (nearer to 1 than 0) Correlation Coefficient using the Least-Squares Regression method with a line of regression.\n", "\n", ":::{seealso}\n", "Run the notebook accompanying this lesson and look at the Month to Price scatterplot. Does the data associating Month to Price for pumpkin sales seem to have a high or low correlation, according to wer visual interpretation of the scatterplot? Does that change if we use a more fine-grained measure instead of `Month`, eg. *day of the year* (i.e. number of days since the beginning of the year)?\n", ":::\n", "\n", "In the code below, we will assume that we have cleaned up the data, and obtained a data frame called `new_pumpkins`, similar to the following:" ] }, { "cell_type": "code", "execution_count": 4, "id": "93a478b1", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [ "output_scroll" ] }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
City NameTypePackageVarietySub VarietyGradeDateLow PriceHigh PriceMostly Low...Unit of SaleQualityConditionAppearanceStorageCropRepackTrans ModeUnnamed: 24Unnamed: 25
0BALTIMORENaN24 inch binsNaNNaNNaN4/29/17270.0280.0270.0...NaNNaNNaNNaNNaNNaNENaNNaNNaN
1BALTIMORENaN24 inch binsNaNNaNNaN5/6/17270.0280.0270.0...NaNNaNNaNNaNNaNNaNENaNNaNNaN
2BALTIMORENaN24 inch binsHOWDEN TYPENaNNaN9/24/16160.0160.0160.0...NaNNaNNaNNaNNaNNaNNNaNNaNNaN
3BALTIMORENaN24 inch binsHOWDEN TYPENaNNaN9/24/16160.0160.0160.0...NaNNaNNaNNaNNaNNaNNNaNNaNNaN
4BALTIMORENaN24 inch binsHOWDEN TYPENaNNaN11/5/1690.0100.090.0...NaNNaNNaNNaNNaNNaNNNaNNaNNaN
\n", "

5 rows × 26 columns

\n", "
" ], "text/plain": [ " City Name Type Package Variety Sub Variety Grade Date \\\n", "0 BALTIMORE NaN 24 inch bins NaN NaN NaN 4/29/17 \n", "1 BALTIMORE NaN 24 inch bins NaN NaN NaN 5/6/17 \n", "2 BALTIMORE NaN 24 inch bins HOWDEN TYPE NaN NaN 9/24/16 \n", "3 BALTIMORE NaN 24 inch bins HOWDEN TYPE NaN NaN 9/24/16 \n", "4 BALTIMORE NaN 24 inch bins HOWDEN TYPE NaN NaN 11/5/16 \n", "\n", " Low Price High Price Mostly Low ... Unit of Sale Quality Condition \\\n", "0 270.0 280.0 270.0 ... NaN NaN NaN \n", "1 270.0 280.0 270.0 ... NaN NaN NaN \n", "2 160.0 160.0 160.0 ... NaN NaN NaN \n", "3 160.0 160.0 160.0 ... NaN NaN NaN \n", "4 90.0 100.0 90.0 ... NaN NaN NaN \n", "\n", " Appearance Storage Crop Repack Trans Mode Unnamed: 24 Unnamed: 25 \n", "0 NaN NaN NaN E NaN NaN NaN \n", "1 NaN NaN NaN E NaN NaN NaN \n", "2 NaN NaN NaN N NaN NaN NaN \n", "3 NaN NaN NaN N NaN NaN NaN \n", "4 NaN NaN NaN N NaN NaN NaN \n", "\n", "[5 rows x 26 columns]" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import pandas as pd\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "from datetime import datetime\n", "\n", "pumpkins = pd.read_csv(\n", " \"https://static-1300131294.cos.accelerate.myqcloud.com/data/us-pumpkins.csv\"\n", ")\n", "\n", "pumpkins.head()" ] }, { "cell_type": "code", "execution_count": 5, "id": "351f7c01", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
MonthDayOfYearVarietyCityPackageLow PriceHigh PricePrice
709267PIE TYPEBALTIMORE1 1/9 bushel cartons15.015.013.636364
719267PIE TYPEBALTIMORE1 1/9 bushel cartons18.018.016.363636
7210274PIE TYPEBALTIMORE1 1/9 bushel cartons18.018.016.363636
7310274PIE TYPEBALTIMORE1 1/9 bushel cartons17.017.015.454545
7410281PIE TYPEBALTIMORE1 1/9 bushel cartons15.015.013.636364
\n", "
" ], "text/plain": [ " Month DayOfYear Variety City Package Low Price \\\n", "70 9 267 PIE TYPE BALTIMORE 1 1/9 bushel cartons 15.0 \n", "71 9 267 PIE TYPE BALTIMORE 1 1/9 bushel cartons 18.0 \n", "72 10 274 PIE TYPE BALTIMORE 1 1/9 bushel cartons 18.0 \n", "73 10 274 PIE TYPE BALTIMORE 1 1/9 bushel cartons 17.0 \n", "74 10 281 PIE TYPE BALTIMORE 1 1/9 bushel cartons 15.0 \n", "\n", " High Price Price \n", "70 15.0 13.636364 \n", "71 18.0 16.363636 \n", "72 18.0 16.363636 \n", "73 17.0 15.454545 \n", "74 15.0 13.636364 " ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "pumpkins = pumpkins[pumpkins[\"Package\"].str.contains(\"bushel\", case=True, regex=True)]\n", "\n", "new_columns = [\n", " \"Package\",\n", " \"Variety\",\n", " \"City Name\",\n", " \"Month\",\n", " \"Low Price\",\n", " \"High Price\",\n", " \"Date\",\n", "]\n", "pumpkins = pumpkins.drop([c for c in pumpkins.columns if c not in new_columns], axis=1)\n", "\n", "price = (pumpkins[\"Low Price\"] + pumpkins[\"High Price\"]) / 2\n", "\n", "month = pd.DatetimeIndex(pumpkins[\"Date\"]).month\n", "day_of_year = pd.to_datetime(pumpkins[\"Date\"]).apply(\n", " lambda dt: (dt - datetime(dt.year, 1, 1)).days\n", ")\n", "\n", "new_pumpkins = pd.DataFrame(\n", " {\n", " \"Month\": month,\n", " \"DayOfYear\": day_of_year,\n", " \"Variety\": pumpkins[\"Variety\"],\n", " \"City\": pumpkins[\"City Name\"],\n", " \"Package\": pumpkins[\"Package\"],\n", " \"Low Price\": pumpkins[\"Low Price\"],\n", " \"High Price\": pumpkins[\"High Price\"],\n", " \"Price\": price,\n", " }\n", ")\n", "\n", "new_pumpkins.loc[new_pumpkins[\"Package\"].str.contains(\"1 1/9\"), \"Price\"] = price / 1.1\n", "new_pumpkins.loc[new_pumpkins[\"Package\"].str.contains(\"1/2\"), \"Price\"] = price * 2\n", "\n", "new_pumpkins.head()" ] }, { "cell_type": "markdown", "id": "75dff0fd", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "A basic scatterplot reminds us that we only have monthly data from August through December. We probably need more data to be able to draw conclusions in a linear fashion." ] }, { "cell_type": "code", "execution_count": 6, "id": "62162c46", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAh8AAAGdCAYAAACyzRGfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABENklEQVR4nO3df3hT9d0//mfa0qSFNrNgm9SWEpTiugofEJAiBYYWC6z+vO6hKMLc3TF+bFbuDQWmlE0psns4lUumvZxMmIPvfQsORqn0HqOUUaRKmYU6RE2hSmu1ZUlLf9HmfP/oEkmbpOekJyfnpM/HdeW6zDmvhNfx3TavvM/7h04QBAFERERECgkLdgJEREQ0uLD4ICIiIkWx+CAiIiJFsfggIiIiRbH4ICIiIkWx+CAiIiJFsfggIiIiRbH4ICIiIkVFBDuB3hwOBy5duoSYmBjodLpgp0NEREQiCIKA5uZmJCYmIizMd9+G6oqPS5cuITk5OdhpEBERkR9qa2uRlJTkM0Z1xUdMTAyAnuRjY2ODnA0RERGJYbfbkZyc7Poc90V1xYfzVktsbCyLDyIiIo0RM2SCA06JiIhIUSw+iIiISFEsPoiIiEhRLD6IiIhIUSw+iIiISFEsPoiIiEhRLD6IiIhIUSw+iIiISFGqW2SMyJtuh4CT1iY0NLcjPsaAKZY4hIdx/x8iIq1h8UGaUHymDhv2V6PO1u46ZjYasD4nDdnp5iBmRkREUvG2C6le8Zk6LNt5yq3wAIB6WzuW7TyF4jN1QcqMiIj8weKDVK3bIWDD/moIHs45j23YX41uh6cIIiJSIxYfpGonrU19ejyuJQCos7XjpLVJuaSIiGhAWHyQqjU0ey88/IkjIqLgY/FBqhYfY5A1joiIgo/FB6naFEsczEYDvE2o1aFn1ssUS5ySaRER0QCw+CBVCw/TYX1OGgD0KUCcz9fnpHG9DyIiDWHxQaqXnW7GtkcmIiFW73Y8IVaPbY9M5DofREQaw+KDNMRb3wcREWkJiw9SPeciY/V29xktX9q5yBgRkRax+CBV4yJjREShh8UHqRoXGSMiCj2Sio/8/HzodDq3h8lkcp0XBAH5+flITExEVFQUZs2ahbNnz8qeNA0eXGSMiCj0SO75+M53voO6ujrXo6qqynVu8+bN2LJlC7Zu3YqKigqYTCZkZWWhublZ1qRp8OAiY0REoUdy8REREQGTyeR6XH/99QB6ej1++9vfYt26dbj//vuRnp6OP/zhD2htbcVbb70le+I0OHCRMSKi0CO5+Dh//jwSExNhsVjw4IMP4rPPPgMAWK1W1NfXY86cOa5YvV6PmTNn4vjx417fr6OjA3a73e1B5MRFxoiIQo+k4uO2227Dm2++iXfffReFhYWor6/HtGnT0NjYiPr6egBAQkKC22sSEhJc5zwpKCiA0Wh0PZKTk/24DAplzkXGTEb3Wysmo4GLjBERaZBOEAS/5yheuXIFN954I1avXo2pU6fi9ttvx6VLl2A2f/NhkJubi9raWhQXF3t8j46ODnR0dLie2+12JCcnw2azITY21t/UKAR1OwSctDahobkd8TE9t1rY40FEpA52ux1Go1HU53fEQP6hoUOH4pZbbsH58+dx7733AgDq6+vdio+GhoY+vSHX0uv10Ov1Xs8TOYWH6ZBx4/Bgp0FERAM0oHU+Ojo68NFHH8FsNsNiscBkMqGkpMR1vrOzE6WlpZg2bdqAEyUiIqLQIKnn42c/+xlycnIwcuRINDQ04Nlnn4XdbsfixYuh0+mQl5eHjRs3YsyYMRgzZgw2btyI6OhoLFy4MFD5ExERkcZIKj4+//xzPPTQQ/j6669x/fXXY+rUqThx4gRSUlIAAKtXr0ZbWxuWL1+Oy5cv47bbbsOhQ4cQExMTkOSJiIhIewY04DQQpAxYISIiInWQ8vnNvV2IiIhIUSw+iIiISFEsPoiIiEhRLD6IiIhIUQNaZIxISW2d3dhYVI2axlaMGh6NtfPSEBUZHuy0iIhIIhYfpAm5b1agpLrB9bzsPLDjxEVkpcWj8NHJQcyMiIik4m0XUr3ehce1SqobkPtmhcIZERHRQLD4IFVr6+z2Wng4lVQ3oK2zW6GMiIhooFh8kKptLKqWNY6IiIKPxQepWk1jq6xxREQUfCw+SNVGDY+WNY6IiIKPxQep2tp5abLGERFR8LH4IFWLigxHVlq8z5istHiu90FEpCEsPkj1Ch+d7LUA4TofRETaw+KDNOGBiUkwxerdjpli9XhgYlKQMiIiIn+x+CDVKz5Th2U7T6He3uF2/Et7B5btPIXiM3VByoyIiPzB4oNUrdshYMP+aggezjmPbdhfjW6HpwgiIlIj7u1CqnbS2oQ6W7vX8wKAOls7TlqbkHHjcOUSI8la2rvwxO5KXLzchpHXReGFBRMwzMA/QUSDEX/zSdUamr0XHv7EUXDcvbUMH35udz0/V9+M9Px3MS4pFvtWZgYxMyIKBt52IVWLjzHIGkfK6114XOvDz+24e2uZwhkRUbCx+CBVm2KJg9logM7LeR0As9GAKZY4JdMikVrau7wWHk4ffm5HS3uXQhkRkRqw+CBVCw/TYX2O79VL1+ekITzMW3lCwfTE7kpZ44goNLD4INXLTjfjTi+LjN2ZFo/sdLPCGZFYFy+3yRpHRKGBxQepXkFRNUqqGzyeK6luQEFRtcIZkVgjr4uSNY6IQgOLD1K1zi4HCsusPmMKy6zo7HIolBFJ8cKCCbLGEVFoYPFBqrajvAb9rR/mEHriSH2GGSIwLinWZ8y4pFiu90E0yLD4IFW70NQqaxwpb9/KTK8FCNf5IBqc+HWDVC0lLlrWOAqOfSszucIpEbnoBEFQ1aYYdrsdRqMRNpsNsbG+u2sp9HV2OXDz0wd93noJ0wH//NVcREawI4+IKFikfH7zrzWpWmREGHIzLT5jcjMtLDyIFNDtEFD+aSP+fPoLlH/ayA0dyW8D+otdUFAAnU6HvLw817ElS5ZAp9O5PaZOnTrQPGkQWzMvDUtnWNB7HbEwHbB0hgVr5vlehIyIBq74TB2mP38YDxWewOO7TuOhwhOY/vxhFJ+pC3ZqpEF+33CtqKjAa6+9hnHjxvU5l52djTfeeMP1PDIy0t9/hghATwHyX3Nuxo7yGlxoakVKXDQWZYxijweRAorP1GHZzlPo3c9Rb2vHsp2nsO2RiVzsjyTxq/hoaWnBww8/jMLCQjz77LN9zuv1ephMpgEnR3StyIgw/DBzdLDTIBpUuh0CNuyv7lN4AICAnv2VNuyvRlaaidsckGh+fW1csWIF5s+fjzvvvNPj+SNHjiA+Ph6pqanIzc1FQ4Pn1SkBoKOjA3a73e1BRETqcNLahDpbu9fzAoA6WztOWpuUS4o0T3LPx65du3Dq1ClUVFR4PD937lz8x3/8B1JSUmC1WvH0009j9uzZ+OCDD6DX6/vEFxQUYMOGDdIzJyKigGto9l54+BNHBEgsPmpra/H444/j0KFDMBgMHmMWLFjg+u/09HRMmjQJKSkpOHDgAO6///4+8WvWrMGqVatcz+12O5KTk6WkRUREARIf4/lvvb9xRIDE4uODDz5AQ0MDbr31Vtex7u5uHD16FFu3bkVHRwfCw8PdXmM2m5GSkoLz5897fE+9Xu+xR4SIiIJviiUOZqMB9bZ2j+M+dABMRgOmWOKUTo00TNKYjzvuuANVVVU4ffq06zFp0iQ8/PDDOH36dJ/CAwAaGxtRW1sLs5kjoYmItCY8TIf1OT3T2XsPJ3U+X5+TxsGmJImkno+YmBikp6e7HRs6dCiGDx+O9PR0tLS0ID8/Hw888ADMZjNqamqwdu1ajBgxAvfdd5+siRMRkTKy083Y9shEbNhf7Tb41GQ0YH1OGqfZkmSybqwQHh6OqqoqvPnmm/jXv/4Fs9mM7373u9i9ezdiYmLk/KeIiEhB2elmZKWZcNLahIbmdsTH9NxqYY8H+YN7uxAREdGAcW8XIiIiUi0WH0RERKQoFh9ERESkKBYfREREpCgWH0RERKQoFh9ERESkKBYfREREpCgWH0RERKQoFh9ERESkKBYfREREpCgWH0RERKQoFh9ERESkKBYfREREpKiIYCdAJFa3Q+B23hp27lIz5r18FN0CEK4Din4yA2MTY4KdFknQ2eXAjvIaXGhqRUpcNBZljEJkBL/DknQ6QRCEYCdxLSlb8tLgUXymDhv2V6PO1u46ZjYasD4nDdnp5iBmRmKMeuqA13M1m+YrmAn5q6CoGoVlVjiu+cQI0wG5mRasmZcWvMRINaR8frNkJdUrPlOHZTtPuRUeAFBva8eynadQfKYuSJmRGL4KDzHnKfgKiqrx6lH3wgMAHALw6lErCoqqg5MYaRaLD1K1boeADfur4al7znlsw/5qdPf+q0iqcO5Ss6xxpLzOLgcKy6w+YwrLrOjsciiUEYUCFh+kaietTX16PK4lAKizteOktUm5pEi0eS8flTWOlLejvKZPj0dvDqEnjkgsFh+kag3N3gsPf+JIWd0iO6TExpHyLjS1yhpHBLD4IJWLjzHIGkfKChc5GUlsHCkvJS5a1jgigMUHqdwUSxzMRgO8fTbp0DPrZYolTsm0SKSin8yQNY6UtyhjFPqb0R6m64kjEovFB6laeJgO63N6pvH1/vvnfL4+J43rfaiU2HU8uN6HekVGhCE30+IzJjfTwvU+SBL+tJDqZaebse2RiTAZ3W+tmIwGbHtkItf5ULn+1vHgOh/qt2ZeGpbOsPTpAQnTAUtncJ0Pko6LjJFmcIVTbeMKp9rHFU7JFymf3yw+iIiIaMC4wikRERGpFosPIiIiUhSLDyIiIlIUiw8iIiJSFIsPIiIiUtSAio+CggLodDrk5eW5jgmCgPz8fCQmJiIqKgqzZs3C2bNnB5onEbodAso/bcSfT3+B8k8buZOtxnzR1IZx+cW4cc0BjMsvxhdNbcFOiSSytV7FA6/8HRkFf8UDr/wdttarwU6JNCrC3xdWVFTgtddew7hx49yOb968GVu2bMH27duRmpqKZ599FllZWTh37hxiYjinn/xTfKYOG/ZXu+1wazYasD4njYuMaUDquiJ0XrN7nL29G7dvPozIcB0+fm5eEDMjsWb++jAuNH5TMNbZ2jH+l4eQMjwKpT+fHcTMSIv86vloaWnBww8/jMLCQlx33XWu44Ig4Le//S3WrVuH+++/H+np6fjDH/6A1tZWvPXWW7IlTYNL8Zk6LNt5yq3wAIB6WzuW7TyF4jN1QcqMxOhdeFyrs1tA6roihTMiqXoXHte60NiGmb8+rHBGpHV+FR8rVqzA/Pnzceedd7odt1qtqK+vx5w5c1zH9Ho9Zs6ciePHjw8sUxqUuh0CNuyvhqePLuexDfureQtGpb5oavNaeDh1dgu8BaNittarXgsPpwuNbbwFQ5JILj527dqFU6dOoaCgoM+5+vp6AEBCQoLb8YSEBNe53jo6OmC3290eRE4nrU19ejyuJaCn+/ektUm5pEi0uS+VyhpHynts+0lZ44gAicVHbW0tHn/8cezcuRMGg8FrnE7nvt+GIAh9jjkVFBTAaDS6HsnJyVJSohDX0Oy98PAnjpR1paNb1jhS3iUfxb8/cUSAxOLjgw8+QENDA2699VZEREQgIiICpaWleOmllxAREeHq8ejdy9HQ0NCnN8RpzZo1sNlsrkdtba2fl0KhKD7Ge5HrTxwpa6g+XNY4Ul6iUdzvltg4IkBi8XHHHXegqqoKp0+fdj0mTZqEhx9+GKdPn8bo0aNhMplQUlLiek1nZydKS0sxbdo0j++p1+sRGxvr9iBymmKJg9logLe9a3XomfUyxRKnZFok0sGfzpQ1jpT3+yVTZI0jAiROtY2JiUF6errbsaFDh2L48OGu43l5edi4cSPGjBmDMWPGYOPGjYiOjsbChQvly5oGjfAwHdbnpGHZzlPQAW4DT50FyfqcNISHeStPKJhuiItCZLjO56DTyHAdboiLUjArksIYPQQpw6N8DjpNGR4FY/QQBbMirZN9hdPVq1cjLy8Py5cvx6RJk/DFF1/g0KFDXOOD/Jadbsa2RybC1Ktb12Q0YNsjE7nOh8p9/Nw8RIZ7Lg65zoc2lP58NlKGey4Quc4H+UMnCIKq5ija7XYYjUbYbDbegiE33Q4BJ61NaGhuR3xMz60W9nhoxxdNbZj7UimudHRjqD4cB386kz0eGmNrvYrHtp/EJVs7Eo0G/H7JFPZ4kIuUz28WH0RERDRgUj6/ubEcERERKYrFBxERESmKxQcREREpisUHERERKUrSOh9axpkS2veVvQP3vXIMTVeuIm7oEOxdPh3Xx+qDnRaJ1NnlwI7yGlxoakVKXDQWZYxCZAS//xANRoNitkvxmTps2F/ttkGZ2WjA+pw0rhGhEePy34W9vavP8VhDBD7MvysIGZEUBUXVKCyz4trNh8N0QG6mBWvmpQUvMSKSDWe7XKP4TB2W7TzVZ2fUels7lu08heIzdUHKjMTyVngAgL29C+Py31U4I5KioKgarx51LzwAwCEArx61oqCoOjiJEVHQhHTx0e0QsGF/NTx17TiPbdhfje7efxVJNb6yd3gtPJzs7V34yt6hUEYkRWeXA4VlVp8xhWVWdHY5FMqIiNQgpIuPk9amPj0e1xIA1NnacdLapFxSJMl9rxyTNY6UtaO8pk+PR28OoSeOiAaPkC4+Gpq9Fx7+xJHymq5clTWOlHWhqVXWOCIKDSFdfMTHGPoPkhBHyosbKm7fCLFxpKyUuGhZ44goNIR08THFEgez0QBvE2p16Jn1MsUSp2RaJMHe5dNljSNlLcoYhf5mtIfpeuKIaPAI6eIjPEyH9Tk90/h6//1zPl+fk8b1PlTs+lg9Yg2+l6OJNURwvQ+ViowIQ26mxWdMbqaF630QDTIh/xufnW7GtkcmwmR0v7ViMhqw7ZGJXOdDAz7Mv8trAcJ1PtRvzbw0LJ1h6dMDEqYDls7gOh9Eg9GgWGQM4AqnoYArnGobVzglCm1SPr8HTfFBREREgcMVTomIiEi1Bs3GckQUXLz1SUROLD6IKOC4uSMRXYu3XYgooLi5IxH1xuKDiAKGmzsSkSe87UKa8cFnl/HAa8ddz9/+0TTcOvq6IGZE/ZGyuWPGjcOVS4z8wnE7JBcWH6QJo5460OeYsxCp2TRf6XRIJG7uGDo4bofkxNsupHqeCg8p5yl4uLljaOC4HZIbiw9StQ8+uyxrHClLuCpuLIfYOFIex+1QILD4IFW7doyHHHGkrIXb35M1jpQnZdwOkVgsPoiIyCuO26FAYPFBRERecdwOBQKLD1K1t380TdY4Utaux6bKGkfKm2KJg9logLcJtTr0zHqZYolTMi3SOBYfpGpi1/Hgeh/qNDVV3NodYuNIeeFhOqzPSQOAPgWI8/n6nDSu90GSSCo+tm3bhnHjxiE2NhaxsbHIyMjAwYMHXeeXLFkCnU7n9pg6ld9oaGD6W8eD63yoG9tP+7LTzbgzLb7PjBcBwJ1p8VzngySTVHwkJSVh06ZNeP/99/H+++9j9uzZuOeee3D27FlXTHZ2Nurq6lyPoqIi2ZOmwSX3zYoBnafg4jot2ldQVI2S6gaP50qqG1BQVK1wRqR1koqPnJwczJs3D6mpqUhNTcVzzz2HYcOG4cSJE64YvV4Pk8nkesTF8T4g+a+ts9vrHz2nkuoGtHV2K5QRSXHuUrOscaS8zi4HCsusPmMKy6zo7HIolBGFAr/HfHR3d2PXrl24cuUKMjIyXMePHDmC+Ph4pKamIjc3Fw0Nvj84Ojo6YLfb3R5EThtFfqMSG0fKmvfyUVnjSHk7ymvQ3/phDqEnjkgsycVHVVUVhg0bBr1ejx//+MfYu3cv0tJ6BiPNnTsXf/zjH3H48GH85je/QUVFBWbPno2Ojg6v71dQUACj0eh6JCcn+381FHJqGltljSNldYtc9FJsHCnvQpO43y2xcUSAH8XH2LFjcfr0aZw4cQLLli3D4sWLUV3d861zwYIFmD9/PtLT05GTk4ODBw/i448/xoED3u/prlmzBjabzfWora31/2oo5IwaHi1rHCkrXOQECLFxpLyUOHG/W2LjiAA/io/IyEjcdNNNmDRpEgoKCjB+/Hi8+OKLHmPNZjNSUlJw/vx5r++n1+tds2ecDyKntfPSZI0jZRX9ZIascaS8RRmj0N8s2jBdTxyRWANe50MQBK+3VRobG1FbWwuzmdOwyD9RkeHISov3GZOVFo+oyHCFMiIpxibGyBpHyouMCENupsVnTG6mBZERXDaKxJP007J27VqUlZWhpqYGVVVVWLduHY4cOYKHH34YLS0t+NnPfoby8nLU1NTgyJEjyMnJwYgRI3DfffcFKn8aBAofnTyg8xRcXOdD+9b007PY33mi3iQVH19++SUWLVqEsWPH4o477sB7772H4uJiZGVlITw8HFVVVbjnnnuQmpqKxYsXIzU1FeXl5YiJ4bca8t/oNb7XgejvPAUX1/nQPv4OktwipAS//vrrXs9FRUXh3XffHXBCRNeyNlwRNc3P2nAFlvihyiRFolVdtImOu2WkMcDZkD/4O0iBwJt0pGrZL5bKGkfKuvuVY7LGkfL4O0iBwOKDVK1D5AIQYuNIWWJbha2nXvwdpEBg8UGqphe5AITYOFKW2FZh66kXfwcpEFh8kKoVPz5T1jhS1r7l02WNI+Xxd5ACgcUHqZolfqioBY440E2dxA4i5WBT9eLvIAUCiw9Svc8K5nv94xem6zlP6sV1PrSPv4MkNxYfpAnDh0ZKOk5E8kqOi5J0nMgXFh+kepOfLcFXLZ0ez33V0onJz5YonBFJwUXGtG/mrw/jQmObx3MXGtsw89eHFc6ItI7FB6laU0un18LD6auWTjT1E0PB8Zf3xO1SLTaOlGdrveq18HC60NgGW+tVhTKiUMDig1TtwdeOyxpHylq590NZ40h5j20/KWscEcDig1SuoVlcj4bYOCKS5pKtXdY4IoDFB6lcfIy4AaVi44hImkSjQdY4IoDFB6ncrh9NkzWOlLX1vnGyxpHyfr9kiqxxRACLD1K5uGGRuH6Y716N64dFIq6fGAqO792WLGscKc8YPQQpw31Pp00ZHgVj9BCFMqJQwOKDVK/iF1mIGuL5RzVqSBgqfpGlcEYkBRcZ077Sn8/2WoCkDI9C6c9nK5wRaR2LD1K9gqJqtF11eDzXdtWBgqJqhTMiKYrP1HndOE737/OkfqU/n41/PDMHt478FsxGA24d+S3845k5LDzILzpBEFS1D7LdbofRaITNZkNsbGyw06Eg6+xy4OanD8Lh46c0TAf881dzERnBWlptuh0Cpj9/GHVeZkLoAJiMBhx7cjbC+9tAhIhUTcrnN/9ak6rtKK/xWXgAgEPoiSP1OWlt8lp4AIAAoM7WjpPWJuWSIqKgY/FBqnahqVXWOFJWQ7O4tR/ExhFRaGDxQaqWEhctaxwpKz5G3NoPYuOIKDSw+CBVW3hbiqxxpKwpljiYjQafA07NRgOmWOKUTIuIgozFB6na6dp/yRpHygoP02F9TprPmPU5aRxsqhGdXQ68XvYZnvnzGbxe9hk6uzzPQiPqT0SwEyDyhWMGtC873YwYQwTs7V19zsUYIpCdbg5CViRVQVE1CsusbgPAnyv6CLmZFqyZ57vAJOqNPR+kahwzoH2Tny3xWHgAgL29C5OfLVE4I5KqoKgarx619pl55hCAV49audYOScbig1Tt/yV/S9Y4UlZTSye+avG94/BXLZ1o6ieGgqezy4HCMqvPmMIyK2/BkCQsPkjV3nrvgqxxpKwHXzsuaxwpj2vtUCCw+CBV+/TrFlnjSFkNzeJ6NMTGkfK41g4FAosPUrWv7B2yxpGy4mPE7TYsNo6Ux7V2KBBYfJCqJcSKG0gqNo6UtetH02SNI+UtyhiF/mZCh+l64ojEYvFBqmYZMVTWOFJW3LBIXD/Md6/G9cMiEddPDAVPZEQYcjMtPmNyMy3c2JEkkfTTsm3bNowbNw6xsbGIjY1FRkYGDh486DovCALy8/ORmJiIqKgozJo1C2fPnpU9aRo8+K1L+yp+kTWg8xR8a+alYekMS5/fxTAdsHQG1/kg6SQtMpaUlIRNmzbhpptuAgD84Q9/wD333IPKykp85zvfwebNm7FlyxZs374dqampePbZZ5GVlYVz584hJiYmIBdAoS0yIkzUSHt+61KvUU8d6Pd8zab5CmVD/lozLw3/Nedm7CivwYWmVqTERWNRxij+7pFfJBUfOTk5bs+fe+45bNu2DSdOnEBaWhp++9vfYt26dbj//vsB9BQnCQkJeOutt7B06VL5sqZB41j1V6LjpqddH+BsSKqqizbRcbeMNAY4GxqoyIgw/DBzdLDToBDgd8na3d2NXbt24cqVK8jIyIDVakV9fT3mzJnjitHr9Zg5cyaOH/c+h7+jowN2u93tQeT0yJsnZY0jZd39yjFZ44goNEguPqqqqjBs2DDo9Xr8+Mc/xt69e5GWlob6+noAQEJCglt8QkKC65wnBQUFMBqNrkdycrLUlIhIpfq5YyY5johCg+TiY+zYsTh9+jROnDiBZcuWYfHixaiu/mZdf53OfUSSIAh9jl1rzZo1sNlsrkdtba3UlIhIpcTuVcs9bYkGF8nFR2RkJG666SZMmjQJBQUFGD9+PF588UWYTCYA6NPL0dDQ0Kc35Fp6vd41e8b5IHLa+egUWeNIWfuWT5c1johCw4CHKQuCgI6ODlgsFphMJpSUfLNDZWdnJ0pLSzFtGhcQIv+IHUTKwabqJHYQKQebEg0ukoqPtWvXoqysDDU1NaiqqsK6detw5MgRPPzww9DpdMjLy8PGjRuxd+9enDlzBkuWLEF0dDQWLlwYqPxpEOhvGianaaob24+IepM01fbLL7/EokWLUFdXB6PRiHHjxqG4uBhZWT2LBK1evRptbW1Yvnw5Ll++jNtuuw2HDh3iGh80YDWb5uNY9Vdus1p2PjqFPR4aUbNpPqou2nD3K8cgoGeMx77l09njQTRI6QRBUNVAc7vdDqPRCJvNxvEf5Kazy8EFjoiIVErK57ekng+iYCkoqkZhmdVttdPnij5CbiaXdiYi0hoWH6R6BUXVePWotc9xhwDXcRYgRETawT5rUrXOLgcKy/oWHtcqLLOis8uhUEZERDRQLD5I1XaU14jaWG5HeY0i+RAR0cCx+CBVu9DUKmscEREFH8d8kKqlxEXLGkfBc/HrVmS/WIq2qw5EDQlD8eMzMXIE241oMOJUW1K1zi4Hbn76oM9bL2E64J+/mstptyp209oD8DQsJyIM+GQjFxkjCgVSPr/515pULTIiDOk3+P4hTr8hloWHinkrPACgy9FznogGF/7FJlXr7HLgzBd2nzFnvrBztotKXfy61Wvh4dTl6IkjosGDxQepGme7aFv2i6WyxhFRaGDxQarG2S7a1nZVXI+U2DgiCg0sPkjVONtF26KGiPsTIzaOiEIDf+NJ1RZMHilrHCmr+PGZssYRUWhg8UGqtrvioqxxpKyRI6LR30SkiDBwvQ+iQYbFB6kax3xo3ycb53stQLjOB9HgxBVOSdU45iM0fLJxPlc4JSIXrnBKqsYVTomItIErnFLIiIwIQ26mxWdMbqaFhQcRkYbwtgup3pp5aQCAwjKrWw9ImK6n8HCeJyIibeBtF9KMzi4HdpTX4EJTK1LiorEoYxR7PIiIVELK5zd7PkgzIiPC8MPM0cFOg4iIBohfG4mIiEhRLD6IiIhIUSw+iIiISFEc80GawQGnREShgcUHaUJBUXWfqbbPFX3EqbZERBrE4oNUr6CoGq8etfY57hDgOs4ChIhIO9hnTarW2eVAYVnfwuNahWVWdHY5FMqIiIgGisUHqdqO8hqf+7oAPT0gO8prFMmHiIgGjsUHqdqFplZZ44iIKPhYfJCqpcSJ23JdbBwREQWfpOKjoKAAkydPRkxMDOLj43Hvvffi3LlzbjFLliyBTqdze0ydOlXWpGnwWJQxCmE63zFhup44IiLSBknFR2lpKVasWIETJ06gpKQEXV1dmDNnDq5cueIWl52djbq6OtejqKhI1qRp8IiMCMMd3473GXPHt+O53gcRkYZImmpbXFzs9vyNN95AfHw8PvjgA8yYMcN1XK/Xw2QyyZMhDWrdDgFnvrD7jDnzhR3dDgHh/XWREBGRKgzo66LNZgMAxMXFuR0/cuQI4uPjkZqaitzcXDQ0NHh9j46ODtjtdrcHkdNJaxPqbO0+Y+ps7ThpbVIoIyIiGii/iw9BELBq1SpMnz4d6enpruNz587FH//4Rxw+fBi/+c1vUFFRgdmzZ6Ojo8Pj+xQUFMBoNLoeycnJ/qZEIaih2XfhITWOiIiCTycIQj+rKHi2YsUKHDhwAMeOHUNSUpLXuLq6OqSkpGDXrl24//77+5zv6OhwK0zsdjuSk5Nhs9kQGxvrT2oUQso/bcRDhSf6jftT7lRk3DhcgYyIiMgTu90Oo9Eo6vPbr+XVf/KTn2Dfvn04evSoz8IDAMxmM1JSUnD+/HmP5/V6PfR6vT9p0CAwxRIHs9GAels7PFXJOgAmowFTLHEezhIRkRpJuu0iCAJWrlyJPXv24PDhw7BYLP2+prGxEbW1tTCbzX4nSYNXeJgO63N879uyPieNg02JiDREUvGxYsUK7Ny5E2+99RZiYmJQX1+P+vp6tLW1AQBaWlrws5/9DOXl5aipqcGRI0eQk5ODESNG4L777gvIBVDoy043Y+TwKI/nRg6PQnY6C1siIi2RVHxs27YNNpsNs2bNgtlsdj12794NAAgPD0dVVRXuuecepKamYvHixUhNTUV5eTliYmICcgEU+u7eWoYLjW0ez11obMPdW8sUzoiIiAZC0piP/samRkVF4d133x1QQkTXamnvwoef+55+/eHndrS0d2GYwa8hTEREpDAuC0mq9sTuSlnjiIgo+Fh8kKpdvOz5dou/cUREFHwsPkjVRl7neaCpv3FERBR8LD5I1V5YMEHWOCIiCj4WH6RqwwwRGJfke6W8cUmxHGxKRKQhLD5I9fatzPRagIxLisW+lZkKZ0RERAPBr4ukCftWZqKlvQtP7K7ExcttGHldFF5YMIE9HkREGsS/3KQZwwwRKFw8OdhpEBHRAPG2CxERESmKxQcREREpisUHERERKYrFBxERESmKxQcREREpirNdSDM6uxzYUV6DC02tSImLxqKMUYiMYP1MRKQ1LD5IEwqKqlFYZoVD+ObYc0UfITfTgjXz0oKXGBERScbig1SvoKgarx619jnuEOA6zgKEiEg72GdNqtbZ5UBhWd/C41qFZVZ0djkUyoiIiAaKxQep2o7yGrdbLZ44hJ44IiLSBhYfpGoXmlpljSMiouBj8UGqlhIXLWscEREFH4sPUrVFGaMQpvMdE6briSMiIm1g8UGqFhkRhtxMi8+Y3EwL1/sgItIQTrUl1XNOo+29zkeYDlzng4hIg3SCIPQzl0BZdrsdRqMRNpsNsbGxwU6HVIQrnBIRqZeUz2/2fJBmREaE4YeZo4OdBhERDRCLDyJSREt7F57YXYmLl9sw8roovLBgAoYZ+CdIS9iG2tftEHDS2oSG5nbExxgwxRKH8P5G9QcAb7sQUcDdvbUMH35u73N8XFIs9q3MDEJGJBXbUPuKz9Rhw/5q1NnaXcfMRgPW56QhO9084PeX8vnNG+ZEFFDePrQA4MPP7bh7a5nCGZFUbEPtKz5Th2U7T7kVHgBQb2vHsp2nUHymTtF8WHwQUcC0tHd5/dBy+vBzO1rauxTKiKRiG2pft0PAhv3V8HSbw3lsw/5qdPe3l4WMWHyQZrR1duPpd6qw6PX38PQ7VWjr7A52StSPJ3ZXyhpHymMbat9Ja1OfHo9rCQDqbO04aW1SLCdJxUdBQQEmT56MmJgYxMfH495778W5c+fcYgRBQH5+PhITExEVFYVZs2bh7NmzsiZNg0/umxX49jPF2HHiIsrOf40dJy7i288UI/fNimCnRj5cvNwmaxwpj22ofQ3N3gsPf+LkIKn4KC0txYoVK3DixAmUlJSgq6sLc+bMwZUrV1wxmzdvxpYtW7B161ZUVFTAZDIhKysLzc3NsidPg0PumxUoqW7weK6kuoEFiIqNvC5K1jhSHttQ++JjDLLGyWFAs12++uorxMfHo7S0FDNmzIAgCEhMTEReXh6efPJJAEBHRwcSEhLw/PPPY+nSpf2+J2e70LXaOrvx7WeK+4376JfZiIoMVyAjkqKlvQvp+e/2G3cm/y5O2VQptqH2dTsETH/+MOpt7R7HfegAmIwGHHty9oCm3So228VmswEA4uLiAABWqxX19fWYM2eOK0av12PmzJk4fvy4x/fo6OiA3W53exA5bSyqljWOlDXMEIFxSb7/CI1LiuWHloqxDbUvPEyH9Tk921D0Li2cz9fnpCm63offxYcgCFi1ahWmT5+O9PR0AEB9fT0AICEhwS02ISHBda63goICGI1G1yM5OdnflCgE1TS2yhpHytu3MtPrhxfXiNAGtqH2Zaebse2RiTAZ3W+tmIwGbHtkoizrfEjhd6m6cuVKfPjhhzh27Fifczqde/UkCEKfY05r1qzBqlWrXM/tdjsLEHIZNTwaZefFxZF67VuZydUxNY5tqH3Z6WZkpZlUscKpXz81P/nJT7Bv3z4cPXoUSUlJruMmkwlATw+I2fxNFdXQ0NCnN8RJr9dDr9f7kwYNAmvnpWHHiYui4kjdhhkiULh4crDToAFgG2pfeJgOGTcOD3Ya0m67CIKAlStXYs+ePTh8+DAsFovbeYvFApPJhJKSEtexzs5OlJaWYtq0afJkTINKVGQ4stLifcZkpcVzsCkRkYZI6vlYsWIF3nrrLfz5z39GTEyMaxyH0WhEVFQUdDod8vLysHHjRowZMwZjxozBxo0bER0djYULFwbkAij0FT462et026y0eBQ+ym9iWmBrvYrHtp/EJVs7Eo0G/H7JFBijhwQ7LZJALZuSkfZJKj62bdsGAJg1a5bb8TfeeANLliwBAKxevRptbW1Yvnw5Ll++jNtuuw2HDh1CTEyMLAnT4PTxl57XifF2nNRl5q8P40LjN4tQ1dnaMf6Xh5AyPAqlP58dxMxIrEBvSkaDC3e1JdXr/cHVGz/A1I3tp33OTcl6f1g4+zyCMVuC1Ie72lLIsLVe9fnBBQAXGttga72qUEYkBdtP+9S4KRlpH4sPUrXHtp+UNY6UxfbTPjVuSkbax+KDVO2Sjz96/sSRsth+2qfGTclI+wZN8dHtEFD+aSP+fPoLlH/ayC5CjUg0itvoSGwcKYvtp31q3JSMtG9QLE3HUdra9fslUzD+l4dExZH6sP20b4olDmajod9NyaZY4pROjTQs5Hs+nKO0e9+zrLe1Y9nOUyg+UxekzEgMY/SQPhsh9ab7dxypjzF6CFKG+95qPWV4FNtPxdS4KRlpX0gXHxylrX1f2Ts8tt+1hH/HkTqV/ny21wKE02y1QW2bkpH2hfQ6H+WfNuKhwhP9xv0pd6oq1rqnvqZv+is+/1f/A9mSvmXAsafuUCAj8hdXONU+rnBKvkj5/A7pMR8cpa19TVfErf8gNo6Cxxg9BG8vvz3YadAAqGVTMtK+kL7twlHa2hc3VNw3Y7FxREQUfCFdfDhHaXvrFNShZ9YLR2mr197l02WNIyKi4Avp4uPaUdrecJS2ul0fq0eswffdwVhDBK6P1SuUERERDVRIFx9AzyjtH82woHd9EaYDfjTDwlHaGvBh/l2I8PKTGhHWc56IiLQj5IuP4jN1eO2oFb1n0woC8NpRK9f50IDcNyvQ5fB8rsvRc57Ur62zG0+/U4VFr7+Hp9+pQltnd7BTIok6uxx4vewzPPPnM3i97DN0evvFJOpHSE+17XYImP78Ya+bIjlX5jv25GzeelGpts5ufPuZ4n7jPvplNqIiwxXIiPyR+2YFSqob+hzPSotH4aOTg5ARSVVQVI3CMvcvcmE6IDfTgjXzfN/epsFByud3SPd8cDdG7dtYVC1rHCnPW+EBACXVDey50oCComq86qEH2SEArx61ooC/fyRRSBcfXOdD+2oaW2WNI2W1dXZ7LTycSqobeAtGxTq7HCgss/qMKSyz8hYMSRLSxQfX+dC+UcOjZY0jZbHnSvt2lNf06fHozSH0xBGJFdLFB9f50L61Iu8li40jZbHnSvsuNIlrG7FxRECIFx/cjVH7oiLDkZUW7zMmKy2eg01Vij1X2pcSJ65txMYRASFefADcjTEUlJ77akDnKXjYc6V9izJG9VknqbcwXU8ckVghvbGcU3a6GVlpJu7GqEFfNLWhs9v3DefObgFfNLXhhjjP27ZT8NhaxW34Z2u9yt4rlYqMCENupgWvHvU+6DQ304JIbysBEnkwaH5anLsx3vP/bkDGjcNZeGjE3JdKZY0jZX3v5aOyxlFwrJmXhqVeVopeOoPrfJB0g6Lng7TrSoe4KZhi40hZtvYuWeMoeNbMS8N/zbkZO8prcKGpFSlx0ViUMYo9HuQXFh+kakP14bC3919YDNWzy16NxH4s8eNLGyIjwvDDzNHBToNCAH/nSdX+vDxT1jhS1pjrxBWFYuOIKDSw54NU7fC5L0XH/TCe38jU5lyTuNspYuMouD6pb8Hcl0px1QEMCQMO/nQmbjINC3ZaJEFLexee2F2Ji5fbMPK6KLywYAKGGZQvBVh8kKpxkSpt629lTKlxFDyWpw7g2ma66gDu/G0pdACsm+YHKy2S4O6tZfjwc7vr+bn6ZqTnv4txSbHYt1LZ3mPediGVE/upxE8vNRI7FodjdtStd+FxLeHf50ndehce1/rwczvu3lqmaD4sPkjV0swxssaRsrbcPV7WOFLeJ/Ut/Zb2wr/jSJ1a2ru8Fh5OH35uR4uCs84kFx9Hjx5FTk4OEhMTodPp8M4777idX7JkCXQ6ndtj6tSpcuVLg8zhf4pbvVRsHCnrP/+/U7LGkfK41o72PbG7UtY4OUguPq5cuYLx48dj69atXmOys7NRV1fnehQVFQ0oSRq82q6K26ZbbBwRSSP2V4u/gup18XKbrHFykDzgdO7cuZg7d67PGL1eD5PJ5HdSRE6WEdE49om4OCKS35AwcYXFEN7EV62R10XhXH2zqDilBOTH5ciRI4iPj0dqaipyc3PR0NAQiH+GBgFuTKZtb/9omqxxpLyDP50paxwp74UFE2SNk4PsxcfcuXPxxz/+EYcPH8ZvfvMbVFRUYPbs2ejo6PAY39HRAbvd7vYgcoqKDEdWWrzPmKy0eG5KplKdgri+eLFxpLybTMPQ305Yun/HkToNM0RgXFKsz5hxSbGKrvche/GxYMECzJ8/H+np6cjJycHBgwfx8ccf48ABz1OxCgoKYDQaXY/k5GS5UyKN+9LePqDzFDwNzeLaRmwcBYd103yvBQjX+dCGfSszvRYgwVjnI+BljtlsRkpKCs6fP+/x/Jo1a7Bq1SrXc7vdzgKEXKRMEQvGKn3kW3yMQdY4Ch7rpvlc4VTj9q3MVM0KpwEfItTY2Ija2lqYzWaP5/V6PWJjY90eRE5qnCJG4k2xxMFsNPj81mw2GjDFEqdkWuSnYYYIGKOGIDJcB2PUEBb8GhQepoPJaEB8jB4mowHhYf3dVAsMyT85LS0t+OSTb6YfWK1WnD59GnFxcYiLi0N+fj4eeOABmM1m1NTUYO3atRgxYgTuu+8+WROnwUGNU8RIvPAwHdbnpOHHOz2v4yEAWJ+TFrQ/gCTet58+6Dal/esrVzF1018RNSQMH/3K9wxIUofcNytQUv3NBJCy88COExeRlRaPwkcnK5qL5J6P999/HxMmTMCECT2jYletWoUJEybgmWeeQXh4OKqqqnDPPfcgNTUVixcvRmpqKsrLyxETwxUoSTqxU7+UnCJG0qx8y/cCYv2dp+DrXXhcq+2qA99++qDCGZFUvQuPa5VUNyD3zQpF85Hc8zFr1iwIgvfFdt99990BJUR0recfGI+SZ0tExZH6XPy6FV39TGTpcvTEjeRaLapU/6/2fhfxa7vqQP2/2mH6FsfuqFFbZ7fXwsOppLoBbZ3dis0c5LIwpGp7Kz+XNY6Ulf2iuCW3xcaR8r738lFZ40h5G4uqZY2TA4sPUrULTa2yxpGyuDy+9tlFbjYmNo6UV9Mo7u+j2Dg5sPggVUuJE9cVLzaOlBUlcs1tsXGkvFiRM1rExpHyRg0X9/dRbJwc+BtPqrZg8khZ40hZ/7v0dlnjSHl/+ckMWeNIeWrcpoLFB6na7oqLssaRsl746zlZ40h5pm95X6fFSffvOFInNW5TweKDVI1jPrSN67RoX0t7F7zPb+wh/DuO1Kvw0cleC5BgrPPBm3SkahzzoW1q3MqbpJGyynDhYmU/wEiawkcno62zGxuLqlHT2IpRw6Oxdl5aUDbmZPFBqrYoYxSeK/oIDh9fvcJ0PXGkPi8smID0/P7X/lFyK2+Shr1XoSUqMhy/uveWYKfB2y6kbpERYcjNtPiMyc20IDKCP8pqpMatvEkarjJMgcC/2KR6a+alYekMC3pv/xGmA5bOsGCNgiO0STq1beVN0ojtlWLvFUnBrxukCWvmpeG/5tyMHeU1uNDUipS4aCzKGMUeD41Q01beJI2z9+rDz+1eY9h7RVLpBF8btQSB3W6H0WiEzWZDbKzv7loiIlLG3VvLPBYg7L0iJymf3yxViYioX+y9Ijnxp4aIiEQZZojgdFqSBW+YExERkaLY80FEiuh2CDhpbUJDczviYwyYYolDeO8pTEQ0KLD4IKKAKz5Thw37q1Fna3cdMxsNWJ+Thux0cxAzI6Jg4G0XIgqo4jN1WLbzlFvhAQD1tnYs23kKxWfqgpQZEQULiw8iCphuh4AN+6s9bkzmPLZhfzW6fa2fT0Qhh8UHEQXMSWtTnx6PawkA6mztOGltUi4pIgo6Fh9EFDANzd4LD3/iiCg0sPggooCJjzHIGkdEoYGzXYgoYKZY4mA2GlBva/c47kMHwGTsmXZL6sfp0iQXFh9EFDDhYTqsz0nDj3ee8nheALA+J40fYBrA6dIkJ952IaKAeuXIJwM6T8HH6dIkNxYfRBQwLe1dPrdiB4APP7ejpb1LoYxIKk6XpkBg8UFEAfPE7kpZ40h5nC5NgcDig4gC5uLlNlnjSHmcLk2BwOKDiAJm5HVRssaR8jhdmgKBxQcRBcwLCybIGkfKc06X9jYfSYeeWS+cLk1SsPggooAZZojAuKRYnzHjkmIxzMBZ/2rlnC4NoE8B4nzO6dIkleTi4+jRo8jJyUFiYiJ0Oh3eeecdt/OCICA/Px+JiYmIiorCrFmzcPbsWbnyJSKN2bcy02sBMi4pFvtWZiqcEUmVnW7GtkcmwmR0v7ViMhqw7ZGJXOeDJJP8dePKlSsYP348fvCDH+CBBx7oc37z5s3YsmULtm/fjtTUVDz77LPIysrCuXPnEBMTI0vSRKQt+1ZmoqW9C0/srsTFy20YeV0UXlgwgT0eGpKdbkZWmokrnJIsdIIg+D05W6fTYe/evbj33nsB9PR6JCYmIi8vD08++SQAoKOjAwkJCXj++eexdOnSft/TbrfDaDTCZrMhNtZ3dy0RERGpg5TPb1nHfFitVtTX12POnDmuY3q9HjNnzsTx48c9vqajowN2u93tQURERKFL1uKjvr4eAJCQkOB2PCEhwXWut4KCAhiNRtcjOTlZzpSIiIhIZQIy20Wnc78HKAhCn2NOa9asgc1mcz1qa2sDkRIRERGphKyjvUwmE4CeHhCz+ZvRzw0NDX16Q5z0ej30er2caRAREZGKydrzYbFYYDKZUFJS4jrW2dmJ0tJSTJs2Tc5/ioiIiDRKcs9HS0sLPvnkmy2wrVYrTp8+jbi4OIwcORJ5eXnYuHEjxowZgzFjxmDjxo2Ijo7GwoULZU2ciIiItEly8fH+++/ju9/9ruv5qlWrAACLFy/G9u3bsXr1arS1tWH58uW4fPkybrvtNhw6dIhrfBARERGAAa7zEQhc54OIiEh7grbOBxEREVF/VLe2sbMjhouNERERaYfzc1vMDRXVFR/Nzc0AwMXGiIiINKi5uRlGo9FnjOrGfDgcDly6dAkxMTFeFybzl91uR3JyMmpra0NyPEmoXx8Q+tfI69O+UL/GUL8+IPSvMVDXJwgCmpubkZiYiLAw36M6VNfzERYWhqSkpID+G7GxsSH5A+UU6tcHhP418vq0L9SvMdSvDwj9awzE9fXX4+HEAadERESkKBYfREREpKhBVXzo9XqsX78+ZPeSCfXrA0L/Gnl92hfq1xjq1weE/jWq4fpUN+CUiIiIQtug6vkgIiKi4GPxQURERIpi8UFERESKYvFBREREigqZ4qOrqwu/+MUvYLFYEBUVhdGjR+OXv/wlHA6Hz9eVlpbi1ltvhcFgwOjRo/G73/1OoYyl8ef6jhw5Ap1O1+fxz3/+U8HMpWlubkZeXh5SUlIQFRWFadOmoaKiwudrtNKGgPTrU3sbHj16FDk5OUhMTIROp8M777zjdl4QBOTn5yMxMRFRUVGYNWsWzp492+/7vv3220hLS4Ner0daWhr27t0boCvwLRDXt337do9t2t7eHsAr8ay/69uzZw/uuusujBgxAjqdDqdPnxb1vmppPyAw16iVNrx69SqefPJJ3HLLLRg6dCgSExPx6KOP4tKlS/2+b6DbMGSKj+effx6/+93vsHXrVnz00UfYvHkzfv3rX+Pll1/2+hqr1Yp58+YhMzMTlZWVWLt2LX7605/i7bffVjBzcfy5Pqdz586hrq7O9RgzZowCGfvnP//zP1FSUoIdO3agqqoKc+bMwZ133okvvvjCY7yW2hCQfn1Oam3DK1euYPz48di6davH85s3b8aWLVuwdetWVFRUwGQyISsry7WHkyfl5eVYsGABFi1ahH/84x9YtGgRvv/97+O9994L1GV4FYjrA3pWlry2Pevq6mAwGAJxCT71d31XrlzB7bffjk2bNol+TzW1HxCYawS00Yatra04deoUnn76aZw6dQp79uzBxx9/jLvvvtvneyrShkKImD9/vvDYY4+5Hbv//vuFRx55xOtrVq9eLdx8881ux5YuXSpMnTo1IDkOhD/X97e//U0AIFy+fDnA2cmjtbVVCA8PF/7yl7+4HR8/frywbt06j6/RUhv6c31aakMAwt69e13PHQ6HYDKZhE2bNrmOtbe3C0ajUfjd737n9X2+//3vC9nZ2W7H7rrrLuHBBx+UPWcp5Lq+N954QzAajQHM1D+9r+9aVqtVACBUVlb2+z5qbT9BkO8atdiGTidPnhQACBcuXPAao0QbhkzPx/Tp0/HXv/4VH3/8MQDgH//4B44dO4Z58+Z5fU15eTnmzJnjduyuu+7C+++/j6tXrwY0X6n8uT6nCRMmwGw244477sDf/va3QKfqt66uLnR3d/f59hAVFYVjx455fI2W2tCf63PSShtey2q1or6+3q199Ho9Zs6ciePHj3t9nbc29fWaYPD3+gCgpaUFKSkpSEpKwve+9z1UVlYGOl3FaKX9BkqrbWiz2aDT6fCtb33La4wSbRgyxceTTz6Jhx56CDfffDOGDBmCCRMmIC8vDw899JDX19TX1yMhIcHtWEJCArq6uvD1118HOmVJ/Lk+s9mM1157DW+//Tb27NmDsWPH4o477sDRo0cVzFy8mJgYZGRk4Fe/+hUuXbqE7u5u7Ny5E++99x7q6uo8vkZLbejP9WmtDa9VX18PAB7bx3nO2+ukviYY/L2+m2++Gdu3b8e+ffvwpz/9CQaDAbfffjvOnz8f0HyVopX2GwittmF7ezueeuopLFy40OeGckq0oep2tfXX7t27sXPnTrz11lv4zne+g9OnTyMvLw+JiYlYvHix19fpdDq358K/F3ztfTzY/Lm+sWPHYuzYsa7nGRkZqK2txX//939jxowZSqUuyY4dO/DYY4/hhhtuQHh4OCZOnIiFCxfi1KlTXl+jlTYEpF+fFtuwN0/t01/b+POaYJGa69SpUzF16lTX89tvvx0TJ07Eyy+/jJdeeilgeSpJS+3nDy224dWrV/Hggw/C4XDglVde6Tc+0G0YMj0fP//5z/HUU0/hwQcfxC233IJFixbhiSeeQEFBgdfXmEymPpVcQ0MDIiIiMHz48ECnLIk/1+fJ1KlTVV2d33jjjSgtLUVLSwtqa2tx8uRJXL16FRaLxWO8ltoQkH59nqi9DZ1MJhMAeGyf3t+qer9O6muCwd/r6y0sLAyTJ0/WRJuKoZX2k5Pa2/Dq1av4/ve/D6vVipKSEp+9HoAybRgyxUdrayvCwtwvJzw83OdU1IyMDJSUlLgdO3ToECZNmoQhQ4YEJE9/+XN9nlRWVsJsNsuZWkAMHToUZrMZly9fxrvvvot77rnHY5yW2vBaYq/PE620ocVigclkcmufzs5OlJaWYtq0aV5f561Nfb0mGPy9vt4EQcDp06c10aZiaKX95KTmNnQWHufPn8f//d//ifpSpkgbyjZ0NcgWL14s3HDDDcJf/vIXwWq1Cnv27BFGjBghrF692hXz1FNPCYsWLXI9/+yzz4To6GjhiSeeEKqrq4XXX39dGDJkiPC///u/wbgEn/y5vhdeeEHYu3ev8PHHHwtnzpwRnnrqKQGA8PbbbwfjEkQpLi4WDh48KHz22WfCoUOHhPHjxwtTpkwROjs7BUHQdhsKgvTrU3sbNjc3C5WVlUJlZaUAQNiyZYtQWVnpGkm/adMmwWg0Cnv27BGqqqqEhx56SDCbzYLdbne9x6JFi4SnnnrK9fzvf/+7EB4eLmzatEn46KOPhE2bNgkRERHCiRMnQuL68vPzheLiYuHTTz8VKisrhR/84AdCRESE8N5776nu+hobG4XKykrhwIEDAgBh165dQmVlpVBXV+f1+tTUfoIQmGvUShtevXpVuPvuu4WkpCTh9OnTQl1dnevR0dHh9fqUaMOQKT7sdrvw+OOPCyNHjhQMBoMwevRoYd26dW7/gxcvXizMnDnT7XVHjhwRJkyYIERGRgqjRo0Stm3bpnDm4vhzfc8//7xw4403CgaDQbjuuuuE6dOnCwcOHAhC9uLt3r1bGD16tBAZGSmYTCZhxYoVwr/+9S/XeS23oSBIvz61t6FzKnDvx+LFiwVB6JmOun79esFkMgl6vV6YMWOGUFVV5fYeM2fOdMU7/c///I8wduxYYciQIcLNN98ctGIrENeXl5cnjBw5UoiMjBSuv/56Yc6cOcLx48cVvKpv9Hd9b7zxhsfz69evd72HmttPEAJzjVppQ+f0YU+Pv/3tb673CEYb6gTh36PziIiIiBQQMmM+iIiISBtYfBAREZGiWHwQERGRolh8EBERkaJYfBAREZGiWHwQERGRolh8EBERkaJYfBAREZGiWHwQERGRolh8EBERkaJYfBAREZGiWHwQERGRov5/gvuBO8ytiM8AAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.scatter(\"Month\", \"Price\", data=new_pumpkins)" ] }, { "cell_type": "markdown", "id": "025f014e", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ ":::{note}\n", "We have performed the same cleaning steps as in the previous section, and have calculated `DayOfYear` column using the following expression: \n", ":::" ] }, { "cell_type": "code", "execution_count": 7, "id": "bc1d5493", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "day_of_year = pd.to_datetime(pumpkins[\"Date\"]).apply(\n", " lambda dt: (dt - datetime(dt.year, 1, 1)).days\n", ")" ] }, { "cell_type": "markdown", "id": "83b08051", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Now that we have an understanding of the math behind linear regression, let's create a Regression model to see if we can predict which package of pumpkins will have the best pumpkin prices. Someone buying pumpkins for a holiday pumpkin patch might want this information to be able to optimize their purchases of pumpkin packages for the patch.\n", "\n", "## Looking for Correlation\n", "\n", "From the previous section, we have probably seen that the average price for different months looks like this:\n", "\n", ":::{figure} https://static-1300131294.cos.ap-shanghai.myqcloud.com/images/ml-regression/barchart.png\n", "---\n", "name: 'Average price by month'\n", "width: 70%\n", "---\n", "Average price by month{cite}`Average_price_by_month`\n", ":::\n", "\n", "This suggests that there should be some correlation, and we can try training a linear regression model to predict the relationship between `Month` and `Price`, or between `DayOfYear` and `Price`. Here is the scatter plot that shows the latter relationship:\n", "\n", ":::{figure} https://static-1300131294.cos.ap-shanghai.myqcloud.com/images/ml-regression/scatter-dayofyear.png\n", "---\n", "name: 'Scatter plot of Price vs. Day of Year'\n", "width: 70%\n", "---\n", "Scatter plot of Price vs. Day of Year{cite}`Scatter_plot_of_Price_vs._Day_of_Year`\n", ":::\n", "\n", "It looks like there are different clusters of prices corresponding to different pumpkin varieties. To confirm this hypothesis, let's plot each pumpkin category using a different color. By passing an `ax` parameter to the `scatter` plotting function we can plot all points on the same graph:" ] }, { "cell_type": "code", "execution_count": 26, "id": "9864740b", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjMAAAGwCAYAAABcnuQpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABqG0lEQVR4nO3deVxU9f4/8NcIA26AisqSqCRLuVeW4a9cQVBvaXlLJQ3SulrpdS3TumXdFNu7bVpdBb3dMMvlmqkoKli5oUaZ+RVEXEq4FCi4sQif3x/nMjIMs8Gc+czyej4e89A558x5fz4fZpg353wWjRBCgIiIiMhJNZNdACIiIqKmYDJDRERETo3JDBERETk1JjNERETk1JjMEBERkVNjMkNEREROjckMEREROTVP2QVQW01NDc6fPw8fHx9oNBrZxSEiIiILCCFw6dIlBAcHo1kz09deXD6ZOX/+PEJCQmQXg4iIiBrh3Llz6NSpk8ljXD6Z8fHxAaA0hq+vr+TSEBERkSXKysoQEhKi+x43xeWTmdpbS76+vkxmiIiInIwlXUTYAZiIiIicGpMZIiIicmpMZoiIiMipuXyfGSJyDTU1NaisrJRdDCKyEa1WCw8PD5uci8kMETm8yspK5Ofno6amRnZRiMiG2rRpg8DAwCbPA8dkhogcmhACBQUF8PDwQEhIiNnJs4jI8QkhcPXqVRQVFQEAgoKCmnQ+JjNE5NCuX7+Oq1evIjg4GC1btpRdHCKykRYtWgAAioqK0LFjxybdcuKfOETk0KqrqwEAXl5ekktCRLZW+wdKVVVVk87DZIaInALXViNyPbb6XDOZISIiIqfGPjNE9eTkAHl5QFgYEB4uuzRERGQOr8wQ/U9JCRAXB0RGAiNHAhERyvMLF2SXjIiITGEyQ/Q/8fFAerr+tvR0YMIEOeUh55aYmAiNRgONRgOtVoubb74Z8+bNw5UrVwAAp0+fhkajQXZ2tt7zhh779+83OH9KSorR42sff//73xEUFISSkhK91/7444/w8vLCf/7zHwDQe42Pjw/69euH9evX645ftGhRg+e/5ZZbVGo9IuswmSGCcmspLQ3438AZnepqZXturpxykY3l5ABbt9rtBxoXF4eCggKcOnUKr776Kj766CPMmzfP5GvS09NRUFCg97jjjjsMjhs3bpzeMVFRUXjiiSf0ts2fPx8hISF4+umnda+rqqpCYmIi4uPjMXr0aN325ORkFBQUICsrC3369MFDDz2Effv26fb36NHDoFzfffedDVqJqOnYZ4YISh8ZU06eZP8Zp1ZSolx6S0u7sS02FkhNBdq2VS2st7c3AgMDAQDx8fHYvXs3Nm7ciGXLlhl9jb+/v+41prRo0UI3TwegDF1v2bKlwWtXr16N22+/HV999RX+/Oc/Y/HixSgpKcF7772nd1ztTKyBgYFYvnw51qxZg02bNiEqKgoA4OnpaVG5iGRgMkMEoFs30/vDwuxTDlKJqXuI27bZrRgtWrRo8nwa1rrllluwZMkSPPnkk/Dx8UFSUhK2bt0KX19fo6/RarXw9PS0e1mJGou3mYigdPaNjQXqT0Dp4aFs51UZJ+Yg9xAPHjyIzz//HMOGDTN53IABA9C6dWu9R3X9sltp5syZ6NmzJ0aOHIknn3wSQ4cONXpsRUUFXn31VZSVlemV9ejRowblevzxx5tULiJb4ZUZov9JTVX+UK97JyI6WtlOTkziPcTNmzejdevWuH79OqqqqjB69Gi8//77Jl/zxRdf4NZbb9Xb1tSVhTUaDZ5//nlkZGTghRdeaPCYCRMmwMPDA9euXYOfnx/efPNNjBgxQrc/MjISmzZt0nuNj49Pk8pFZCtMZoj+p21b5Y5Dbq7y/WbPeWY4t42KJN5DHDJkCJYtWwatVovg4GBotVqzrwkJCUGYCmXy9PTU+7e+d955B9HR0fD19UXHjh0N9nt5ealSLiJbYDJDVE94uP0SCkn9Ut1L7T3E9HT9W00eHsqlNxV/2K1atXKaBCAwMNBpykpUH5MZIokcpF+q63Oie4jFxcUoLCzU29amTRs0b95cUokU169fNyiXRqNBQECApBIR3cBkhkiS2n6p9dXtl8pbTjYi8x6ilaKjow22paamYvz48RJKc8OxY8cQFBSkt83b2xvl5eWSSkR0g0YIIWQXQk1lZWXw8/NDaWmpyaGIRPa2dauybIIxW7YAdfpfuq3y8nLk5+cjNDRU+tUJIrItU59va76/pQ7NbmiK7LqTMgkhsGjRIgQHB6NFixYYPHgwjh07JrHERLbDuW2IiGxD+jwz9afIPnr0qG7f66+/jrfffhsffPABsrKyEBgYiJiYGFy6dEliiYlsg3PbEBHZhvRkpnaK7NpHhw4dAChXZd599108//zzePDBB9GzZ0+sWrUKV69exeeff270fBUVFSgrK9N7EDmq1FSlH2pdDtovlYjIYUlPZnJzcxEcHIzQ0FCMHz8ep06dAgDk5+ejsLAQw4cP1x3r7e2NQYMGYe/evUbPl5SUBD8/P90jJCRE9ToQNVZtv9ScHKWPTE6O8pzDsomILCc1menfvz9Wr16NtLQ0fPrppygsLMSAAQP0hibWH/YXEBBgMDywrgULFqC0tFT3OHfunKp1ILKF8HClsy9vLRERWU/q0Oy6U2X36tULUVFR6NatG1atWoW7774bgDKPQV1CCINtdXl7e8Pb21udAhMREZHDkX6bqa5WrVqhV69eyM3N1Y1qqn8VpqioiJM0ERERkY5DJTMVFRU4fvw4goKCEBoaisDAQOzYsUO3v7KyEpmZmRgwYIDEUhIREZEjkXqbad68ebjvvvvQuXNnFBUV6ZadT0hIgEajwaxZs7BkyRKEh4cjPDwcS5YsQcuWLREfHy+z2ERERORApF6Z+fXXXzFhwgRERkbiwQcfhJeXF/bv348uXboAAJ599lnMmjULTz31FPr164fffvsN27dv57LzROQUEhMTodFoMG3aNIN9Tz31FDQaDRITE3XHjhkzxuC1S5cu1Xvdxo0b9foNZmRkQKPR4OLFiwYxIiMj4eXlhd9++03vWFOPlJQUpKSkoE2bNg3WqU2bNkhJSdE9r/va1q1bo0+fPnr7zcU1NaCDyFJSk5k1a9bg/PnzqKysxG+//YZ169ahe/fuuv0ajQaLFi1CQUEBysvLkZmZiZ49e0osMRGRdUJCQrBmzRpcu3ZNt628vBypqano3Lmzydc2b94cr732Gi5cuGB13O+++w7l5eV46KGHdMnFgAED9CYpffjhhxEXF6e3bdy4cVbHSk5ORkFBAX788UeMGzcOjz32GNIaWHjsxIkTerEKCgrQsWNHq+MR1edQfWaIiNSUk6OsiZWba7+Yt99+Ozp37oz169frtq1fvx4hISG47bbbTL42OjoagYGBSEpKsjruihUrEB8fj0mTJmHlypUQQsDLy0tvktIWLVrA29vbYJu12rRpg8DAQHTr1g0LFy5Eu3btsH37doPjOnbsqBcrMDAQzZrxa4iaju8iInJ5JSVAXBwQGaks7hkRoTxvxAWPRnnssceQnJyse75y5UpMnjzZ7Os8PDywZMkSvP/++/j1118tjnfp0iV8+eWXmDhxImJiYnDlyhVkZGQ0puhWqa6uxtq1a1FSUgKtVqt6PKJaTGaIyOXFxwPp6frb0tOBCRPsE3/SpEn47rvvcPr0aZw5cwbff/89Jk6caNFrH3jgAfTt2xcvvfSSxfHWrFmD8PBw9OjRAx4eHhg/fjxWrFjR2OKbNWHCBLRu3Rre3t4YN24c2rVrh8cff9zguE6dOqF169a6R2RkpGplIvcidTQTEZHacnKABrpvoLpa2Z6bq/7My+3bt8eoUaOwatUqCCEwatQotG/f3uLXv/baaxg6dCjmzp1r0fErVqzQS5YmTpyIgQMH4uLFi0Y79jbFO++8g+joaJw7dw5z5szB7NmzEdbAsu/ffvut3gAOT09+BZFt8J1ERC4tL8/0/pMn7bOMxOTJkzF9+nQAwIcffmjVawcOHIjY2FgsXLhQN/rJmF9++QUHDhxAVlYW5s+fr9teXV2N1NRUPPnkk2bj+fr64vLly6iuroZHnWXdq6urcfnyZfj5+ekdHxgYiLCwMISFheHLL7/Ebbfdhn79+ukN6ACA0NBQVZIpIt5mIovI6DgpI6bMuKSObt1M72/gAoIq4uLiUFlZicrKSsTGxlr9+qVLl+Lrr782udAuoFyVGThwIH788UdkZ2frHs8++6zFt5puueUWVFdX44cfftDbfuTIEVRXV5u8PRQWFoaxY8diwYIFFsUisgUmM2SSjI6Tsjpryu4kSuqIiABiY4E6FxgAKM9jY+23uKeHhweOHz+O48eP613tsFSvXr3wyCOP4P333zd6TFVVFf71r39hwoQJ6Nmzp97j8ccfx+HDh/Hjjz+ajdW9e3eMGDECkydPRnp6OvLz85Geno4pU6ZgxIgRBldc6ps7dy6+/vprHDp0SG97UVERCgsL9R5VVVWWNQCRCUxmyCQZHSdlddaU3UmU1JOaCkRH62+Ljla225Ovry98fX0b/fq///3vEEIY3b9p0yYUFxfjgQceMNgXHh6OXr16WXx1Zs2aNYiOjsaTTz6J7t2748knn8SwYcOQakGj9erVC9HR0XjxxRf1tkdGRiIoKEjvcfjwYYvKQ2SKRpj6ZLiAsrIy+Pn5obS0tEm/RNxRTo5ylcLUflv/VSsjpsy4ZF55eTny8/MRGhqK5s2bN+lcublKH5mwMP48iRyBqc+3Nd/fvDJDRlnScdIVYsqMS/YVHg6MGMFEhsjVMJkho2R0nJTVWdNROokSEZH1mMyQUTI6TsrqrOkonUSJiMh6TGbIJBkdJ2V11nSUTqJERGQdTppHJrVtC2zbZt+OkzJiyowLKB2M8/LYMZWIqDGYzJBFwsPt/yUrI6a945aUKEPC6063HxurXA1q29Y+ZSAicna8zUQkEee2ISJqOiYzRJLULoBYXa2/ve4CiEREZB6TGSJJOLdNw2StjVVaCpw/D5SV2S9mebkSt7zcfjGJXBH7zBBJwrlt9BnrP7Rqlbpxy8uB//s/4Pr1G9s8PYFbbwW8vdWJef06cOqUfuLk6wvcfLMSm4iswyszRJJwbht9xvoPzZ2rbtz6iQygPD9+vOnnTkxMhEajMXjs3n0SZWVAcvIS9O/vgZSUpSgrUxKcWikpKWjTpo3e87rnCAgIwH333Ydjx44BAKZMmYJevXqhsrJSrwxbtmyBVqvFn/70pwbLUvdx+vRpAMDevXvh4eGBuLg4gzqdPn0aGo0G2dnZDda5fjlrH01dioLIFCYzRBJxbhuFqf5D330HqLWwcmmpYSJT6/p129xyiouLQ0FBge6Rn18AX99QAMDXXydj0qRn8fXXKwEo8UzdcvL19UVBQQHOnz+Pb775BleuXMGoUaNQWVmJd999F5cuXcJLL72kO/7ixYv4y1/+gueffx5r1qzRK0enTp3wyiuv6G0LCQkBAKxcuRIzZszAd999h7Nnz1pd59py1n2cOXPG6vMQWYoXNImslJYGHDgAREUBMTFNO5fMuW0cibn+Q8YSDmvlFOcgryQPYe3CEO4fjitXTB9/+bJy+6cpvL29ERgYqHteWgr88Qdw+HAmKiquYdq0V7Bly2ocObIHt98+EBUVgLGLGBqNRneuoKAgzJ49G/fffz9OnDiBXr16ISUlBcOHD8eYMWPQv39/zJo1C0FBQXjhhRfg6emJ1q1b687l4eEBHx8fvbIBwJUrV7B27VpkZWWhsLAQKSkpBqtfm1O3nET2wGSGyEJ5eUD//kBx8Y1t/v5AVhYQGtq0c8uaU8dRmOs/1NR+JCXXShC/Lh5peTc65MR2i8XHMakAjE/oU+e732Zq++Fs2rQCw4dPgKenFsOHT8CmTStw++0DLe6nc/HiRXz++ecAAK1WCwAYPHgwnnrqKSQkJODvf/871q5di8OHD8PTigb84osvEBkZicjISEycOBEzZszA3/72N2g0GqvqSWRPvM1EZKH6iQygPL/zTjnlcSWm+g/dcw/wv+/qRotfF4/0U/odctJPpWPqjglGEyVPz6ZflQGAzZs3o3Xr1rrHpEkPoVmzMuzatQ4jRkwEAIwYMRE7d36FZs3KjF6VAYDS0lK0bt0arVq1Qtu2bbFmzRrcf//9uOWWW3THJCUlQaPRYPz48ViyZAluvfVWq8q7YsUKTJyolCsuLg6XL1/Gzp07rTpHbTnrPoYPH27VOYiswWSGyAJpaYaJTK3iYmDHDvuWxxUZ6z/01ltNO29OcQ7S8tJQLfQ75FSLaqTlpUEbkGuQ0NSOZrKFIUOGIDs7W/d47733sG/f5+jc+WZERPQBAERG9kXnzjfjwIE1Js/l4+OD7OxsHD58GMuXL0e3bt2wfPlyvWNatGiBuXPnomXLlpg5c6ZVZT1x4gQOHjyI8ePHAwA8PT0xbtw4rFy50qrz1Jaz7iM5OdmqcxBZg7eZiCxw4IDp/fv2Nb3/jLsz1n+ovFwZtt1YeSWmO+ScvXwSI/qGo6xM6SPTurVtrsjUatWqFcLqjbNftWolcnOP4e67b/wKrqmpQUrKCjz55F+MnqtZs2a6c91yyy0oLCzEuHHjsGfPHr3jPD094eHhYfWtoRUrVuD69eu46aabdNuEENBqtbhw4QLaWrjGRt1yEtkDr8wQWaB/f9P7o6LsUw53EB4OjBhhuz5E3dqZ7pAT1k750vX1BYKDbZvINOTo0aM4dOgQMjIy9K5c7NmzB1lZWfj5558tPtfs2bPx448/YsOGDU0u1/Xr17F69Wq89dZbeuX68ccf0aVLF/z73/9ucgwitfDKDJEFYmOVzr4N3Wry9+dVGUcW4R+B2G6xSD+VrneryUPjgeiboxHub9+e1ytWrMBdd92FgQMHGuyLiorCihUr8M4771h0Ll9fXzz++ON46aWXMGbMmCZ10t28eTMuXLiAKVOmwM/PT2/fn//8Z6xYsQLTp0/XbTtx4oTBObp37w5AuZpTWFhosL9jx45o1ox/Q5Pt8V1FZKGsLCVxqat2NBM5ttSxqYi+Wb9DTvTN0Ugda98JfSorK/HZZ59h7NixDe4fO3YsPvvsM4OJ70yZOXMmjh8/ji+//LJJZVuxYgWio6MNEpnacmVnZ+PIkSO6bePHj8dtt92m9zh//jwAoKysDEFBQQaPoqKiJpWRyBiNEELILoSaysrK4Ofnh9LSUviqff2Y3MKOHUofGVvMM0PmlZeXIz8/H6GhoU2eRTa3OBcnS07q5pkhIrlMfb6t+f52mCsztcMJZ82apdvW0FTgd999t7xCEkFJYF58kYmMMwr3D8eI8BFun8hwgUtyNQ7RZyYrKwuffPIJevfubbAvLi5Ob0ifl5eXPYtGROQyuMAluSrpV2YuX76MRx55BJ9++mmDw/5qpwKvfbRr105CKYmInF/9RAaAwQKXRM5IejLz9NNPY9SoUYiuP1vW/2RkZKBjx46IiIjAE088YbYDWUVFBcrKyvQeRETurrzc+MKZ5ha4JHJ0Ui8srlmzBkeOHEGWkeEgI0aMwEMPPYQuXbogPz8ff/vb3zB06FAcPnwY3kYWMElKSsLLL7+sZrGJiJxORYX5/U3sX00kjbRk5ty5c5g5cya2b99udITCuHHjdP/v2bMn+vXrhy5duuCbb77Bgw8+2OBrFixYgDlz5uiel5WV6Za1JyJyV+YWsLR0gUsiRyQtmTl8+DCKiopwxx136LZVV1djz549+OCDD1BRUQGPeqvOBQUFoUuXLsjNzTV6Xm9vb6NXbYiI3FXz5kpn34ZuNfn68qoMOTdpycywYcNw9OhRvW2PPfYYbrnlFsyfP98gkQGA4uJinDt3DkFBQfYqJhGRy7j5ZuOjmYicmbRkxsfHBz179tTb1qpVK/j7+6Nnz564fPkyFi1ahLFjxyIoKAinT5/GwoUL0b59ezzwwAOSSk1E5Lw8PYGICKWzb0WFcmuJV2TIFTjszAIeHh44evQoVq9ejYsXLyIoKAhDhgzBF198AR8fH9nFIyJyWs2bM4kh1yJ9aHZdGRkZePfddwEALVq0QFpaGoqKilBZWYkzZ84gJSWFnXmJyCnUzmA+bdo0g31PPfUUNBoNEhMT9Y4fM2YMAKX/4IABAwzWcCotLUVISAheeOEFAMDp06cNZkmvfezfvx8AkJKSotvm4eGBtm3bon///njllVdQWlpqsg4ZGRnQaDS4ePGiwb6uXbvqfl/XPm+oHEuXLtUra3Z2doNl9/Pzw913342vv/5aL05KSgratGljspyrVq3CXXfdhVatWsHHxwcDBw7E5s2bDY6rrq7GO++8g969e6N58+Zo06YNRowYge+//94gpkajwa233mpwjrVr10Kj0aBr164NluW///0vtFotPvvsswb3T506VW+C2GvXrqFt27Zo164drl27ZnB83XZt0aIFbrnlFrzxxhuouxKRubat/74YPHiw0f116zZ48GC9WfmNva7ue3z37t0YMmQI2rVrh5YtWyI8PBwJCQm4fv16g+1hKw6VzBARuZKQkBCsWbNG70uqvLwcqamp6Ny5s9HXeXh4YNWqVdi2bRv+/e9/67bPmDED7dq1w4svvqh3fHp6OgoKCvQedQdX+Pr6oqCgAL/++iv27t2Lv/zlL1i9ejX69u2rWxzSFl555RWDcsyYMcPka2rLfuDAAdx1110YO3Ysfv75Z4tjzps3D1OnTsXDDz+MH3/8EQcPHsS9996L0aNH44MPPtAdJ4TA+PHj8corr+Cvf/0rjh8/jszMTISEhGDw4MHYuHGj3nlbtWqFoqIi7Nu3T2/7ypUrTf7sAgICMGrUKL2Z62tdu3YNa9aswZQpU3Tb1q1bh549e6J79+5Yv359g+esbdfjx49j3rx5WLhwIT755BOzbWPsfbF+/Xrd84MHDxoca2y6FAB44oknDM75+uuvAwCOHTuGESNG4M4778SePXtw9OhRvP/++9BqtaipqTFb3qZw2NtMRES2lwMgD0AYAPXXZ7r99ttx6tQprF+/Ho888ggAYP369QgJCcHNZnrdhoeHIykpCTNmzMCQIUOQlZWFNWvW4ODBgwbLuvj7+yMwMNDouTQajW5/UFAQbr31Vtx3333o0aMHnn32WaNXEazl4+NjshwNqS17YGAgFi9ejPfffx+7d+826FPZkP379+Ott97Ce++9p5c0LV68GOXl5ZgzZw5Gjx6NkJAQrF27Fl999RU2bdqE++67T3fsJ598guLiYjz++OOIiYlBq1atAACenp6Ij4/HypUrERUVBQD49ddfkZGRgdmzZyM11fiK61OmTMHo0aNx+vRpvSs4X331FcrLyzFx4kTdthUrVmDixIkQQmDFihW690ldddv18ccfx7Jly7B9+3ZMnTrVZPsYe1/UnUm//H+zJZp7D9Vq2bKl0eN27NiBoKAgXXIDAN26dUNcXJzZ8zYVr8wQkRsoARAHIBLASAAR/3t+QfXIjz32mN5f6StXrsTkyZMteu2MGTPQp08fPProo/jLX/6CF198EX379rVJuTp27IhHHnkEmzZtQnV1tU3O2RRVVVX49NNPAQBardai16SmpqJ169YNfqnPnTsXVVVVWLduHQDg888/R0REhF4iU/fY4uJi7NixQ2/7lClT8MUXX+Dq1asAlNtPcXFxCAgIMFmukSNHIjAwECkpKXrbV65ciTFjxsDf3x8AkJeXh3379uHhhx/Gww8/jL179+KUibUlhBDIyMjA8ePHLW4jewoMDERBQQH27Nlj99hMZojIDcQDSK+3LR3ABNUjT5o0Cd999x1Onz6NM2fO4Pvvv9f7y9wUjUaDZcuWYefOnQgICMBzzz3X4HEDBgxA69at9R6WJCi33HILLl26hOLiYpPHderUyeD8Z8+eNThu/vz5BsdlZGSYPHdt2Zs3b465c+eia9euePjhh82WHQBycnLQrVu3BhcgDg4Ohp+fH3JycnTHNtQHBoBue+2xtfr27Ytu3brhq6++ghACKSkpFiWiHh4eePTRR5GSkqLr25Kfn4/MzEy9W0wrV67EiBEjdH1m4uLisHLlSoPz1bart7c3hgwZAiEE/vrXv5otR2PfF6Z89NFHBudctWoVAOChhx7ChAkTMGjQIAQFBeGBBx7ABx98YJdlhXibiYhcXA6AtAa2V/9vey7UvOXUvn17jBo1CqtWrYIQAqNGjUL79u0tfv3KlSvRsmVL5Ofn49dff22w4+kXX3xh8EXd0Fxd9dV+0Wo0GpPHffvttwajSAcPHmxw3DPPPKPXqRkAbrrpJpPn/uKLL3DLLbcgJycHs2bNwvLly222oLAQwmzd6mro2MmTJyM5ORmdO3fG5cuXMXLkSL2+OMZMmTIFr732Gnbt2oVhw4Zh5cqV6NSpk24dwurqaqxatQr/+Mc/dK+ZOHEiZs+ejZdfflnv51fbrr///juef/55DB06FAMGDDBbhsa+L0x55JFH8Pzzz+tt69ixo+7cycnJePXVV7Fr1y7s378fixcvxmuvvYaDBw+qOkcckxkicnF5ZvafhNr9ZyZPnozp06cDAD788EOLX7dv3z6888472Lp1K15//XVMmTIF6enpBl+6ISEhCAsLs7pcx48fh6+vr+62hzGhoaEGI4o8PQ2/Ptq3b291OUJCQhAeHo7w8HC0bt0aY8eOxS+//KL7gjQlIiIC3333HSorKw2uzpw/fx5lZWUIDw/XHfvLL780eJ7jx48DgO7Yuh555BE8++yzWLRoER599NEG692Q8PBw3HvvvUhOTsaQIUOwatUqPPbYY2jWTLkhkpaWht9++01v2R5ASXK2b9+OESNG6LbVtmtYWBjWrVuHsLAw3H333UYXaK7V2PeFKX5+fmbPedNNN2HSpEmYNGkSXn31VURERGD58uWqrpvI20xE5OK6mdlv21/2DYmLi0NlZSUqKysRGxtr0WuuXbuGhIQETJ06FdHR0fjnP/+JrKwsfPzxxzYpU1FRET7//HOMGTNG9wUr26BBg9CzZ08sXrzYouPHjx+Py5cvN9gmb775JrRarW54+/jx45Gbm2sw9BsA3nrrLfj7+yMmJsZgX7t27XD//fcjMzPT4r5OtaZMmYL169dj3bp1+PXXX/HYY4/p9q1YsQLjx49Hdna23uORRx7BihUrjJ6zbdu2mDFjBubNm6c3PNtRtW3bFkFBQbhy5YqqcXhlhohcXASAWCh9ZOr2F/AAEA17jGry8PDQ/fVv6WX+5557DjU1NXjttdcAAJ07d8Zbb72FOXPmIC4uTu92U3FxMQoLC/Ve36ZNG90ivkIIFBYWQgiBixcvYt++fViyZAn8/Px088DYwqVLlwzK0bJlS/j6+lp8jrlz5+Khhx7Cs88+q7tFVV1drZtDpZaXlxeioqIwc+ZMPPPMM6isrMSYMWNQVVWFzz77DP/4xz/w7rvv6uYmGz9+PL788kskJCTgjTfewLBhw1BWVoYPP/wQmzZtwpdffqkbyVRfSkoKPvroI7NXsOp76KGH8Ne//hVTp07FsGHDdD+z33//HV9//TU2bdpkMGorISEBo0aNwu+//44OHTo0eN6nn34ar732GtatW4c///nPRuObe180xtWrVw3O6e3tjbZt2+Ljjz9GdnY2HnjgAXTr1g3l5eVYvXo1jh07hvfff7/RMS3hGOk4EZGqUqEkLnVF/2+7ffj6+lr8pZ6ZmYkPP/wQKSkpel+wTzzxBAYMGIApU6bo/VUeHR2NoKAgvUfdeVPKysoQFBSEm266CVFRUfj444+RkJCAH374wab9GF588UWDcjz77LNWneNPf/oTunbtqnd15vLly7jtttv0HiNHjgQAvPvuu/joo4+wZs0a9OrVC3fccQcyMzOxceNGveHaGo0Ga9euxfPPP4933nkHt9xyC+69916cOXMGu3fv1k1Y2JAWLVpYncgASiI3fvx4XLhwQe+qzurVq9GqVSsMGzbM4DVDhgyBj48P/vWvfxk9b4cOHTBp0iQsWrTI5Pwt5t4XjfHpp58anHPCBKUj/V133YXLly9j2rRp6NGjBwYNGoT9+/dj48aNGDRoUJPimqMRznCdqgnKysrg5+eH0tJSq/46ICLHUF5ejvz8fISGhjbpL0pFLpQ+MvaZZ4aITDP1+bbm+5u3mYjIjYSDSQyR6+FtJiIiInJqTGaIiIjIqTGZISIiIqfGZIbIAeTkAFu3Arm5skuib/FiYOhQwIajdxsdU+2xCgUFwIkTyr/2Ul4OlJYq/9qTjLiy6kqOzVafa45mIpKopASIjwfS6sy2HxsLpKYCbdvKK9euXUB0NFD3t4NGA2RkAAMH2jtmFTp0OKlba8fWysqUZLK+yEig3gz+NnP9OnDqlBK7lq8vcPPNgIUTzDpNXFl1JedQXFyMoqIiREREGMzBZM33N5MZIoni4oD0dKDu2m8eHsqX+rZt8srVrJl+UlFLowFMTGuhSsxmzQROnTqLqqoqBAcH23y22p9/Nr6v3nxmNpOfDzQ0IWqrVkBoqDoxZcWVVVdybEIIXL16FUVFRWjTpk2D8x1xaDaRE8jJ0b8iU6u6Wtmemws0sFSM6hYvbjipAJTtS5cCRhZvViVmTY0G69cHYeTIfJw5c8amcUtLgYsXje/PzgZsfTGoqgo4f77hfX/8oezXam0bU1ZcWXUl59GmTRsEBgY2+TxMZogkyTOz/uHJk3KSmZ07Te/fvt32yYy5mF9/7YWZM8NRWVlp07gJCcCBA8b39+8PrFpl05DYsweYNs34/k8+UedWnoy4supKzkGr1TZ5Fe9aTGaIJOlmZv1DGy92a7Fhw4Ddu43vHz5cTsxmzZrZYAZgfb17A2vXGt8/bRpg45Do2hUwdYEpNNT2MWXFlVVXcj/sM0MkEfvMyI0pK66sn7uMuI76HifHZ833N4dmE0mUmqr8Uq8rOlrZLlNGhvJlXlftaCZXiikrrqyfu4y4jvoeJ9fCKzNEDiA3V+kjExYmp5+MMUuXKn1khg+3fT8ZR4opK66sn7uMuI76HifHxaHZdTCZISIicj68zURERERug8kMEREROTUmM0REROTUmMwQERGRU2MyQ0RERE6NyQwRERE5NYdJZpKSkqDRaDBr1izdNiEEFi1ahODgYLRo0QKDBw/GsWPH5BWSCMoCkVu3KvNm2Mvs2UCfPsC8efaLCcipqztZsQKYNAlISXH9uGlpwCuvADt22C8muRHhAA4ePCi6du0qevfuLWbOnKnbvnTpUuHj4yPWrVsnjh49KsaNGyeCgoJEWVmZxecuLS0VAERpaakKJSd3UlwsRGysEMrk98ojNlaIkhL1Ym7YoB+v9rF5s3oxhZBTV3dy6JAQWq1++2q1Qvzwg+vFPXlSCH9//Zj+/kKcOqVeTHIN1nx/S5807/Lly7j99tvx0Ucf4dVXX0Xfvn3x7rvvQgiB4OBgzJo1C/PnzwcAVFRUICAgAK+99hqmTp1q0fk5aR7Ziow1ZupPs1+Xmp9crqejLi8voKrKcLtWC9h4YXDpcdu3B4qLDbf7+wN//KFOTHINTjVp3tNPP41Ro0Yhut7iHfn5+SgsLMTwOkv0ent7Y9CgQdi7d6/R81VUVKCsrEzvQdRUOTnKZfK6X+6A8jwtTZ3bMLNnm96v1i0nGXV1JytWNJxQAMp2tW79yIibltZwIgMo23nLiWxFajKzZs0aHDlyBElJSQb7CgsLAQABAQF62wMCAnT7GpKUlAQ/Pz/dIyQkxLaFJreUl2d6/8mTto+5a5fp/Wp9Ecioqzsxt4Dlzp2uE/fAAdP79+2zfUxyT9KSmXPnzmHmzJn47LPP0Lx5c6PHaepdZxdCGGyra8GCBSgtLdU9zp07Z7Myk/vq1s30/rAw28ccOtT0/pgY28cE5NTVnQwebHr/sGGuE7d/f9P7o6JsH5Pck7Q+Mxs3bsQDDzwADw8P3bbq6mpoNBo0a9YMJ06cQFhYGI4cOYLbbrtNd8zo0aPRpk0brFq1yqI47DNDtsI+M+wzYyvsM8M+M2SeU/SZGTZsGI4ePYrs7Gzdo1+/fnjkkUeQnZ2Nm2++GYGBgdhR51p6ZWUlMjMzMWDAAFnFJjeWmqp8mdcVHa1sV8vmzdZttxUZdXUnBw8qCURdWq2y3dXiZmUpiUtd/v7KdiJb8ZQV2MfHBz179tTb1qpVK/j7++u2z5o1C0uWLEF4eDjCw8OxZMkStGzZEvHx8TKKTG6ubVvlqkRurtJvJCwMCA9XN+aoUcoVmHnzlD4yMTHAm2+qGxOQU1d30revciUkJUXpqzJsGJCY6JpxQ0OVKzA7dih9ZKKi1LtFSu5L+tDsugYPHqwbmg0o/WNefvllfPzxx7hw4QL69++PDz/80CAJMoW3mYiIiJyPNd/fDpXMqIHJDBERkfNxij4zRERERLbAZIaIiIicGpMZIiIicmpMZoiIiMipMZkhIiIip8ZkhoiIiJyatEnzXEFOjrIonztMKCajrmlpykJ19p5kS0ZcWXV1RLI+V+70eSZyNUxmGqGkBIiPV76AasXGKlO9t20rr1xqkFHXvDxlgbq667nUTn8eGqpOTFlxZdXVEcn6XLnT55nIVXHSvEZwp0X4ZNRV1sJ0MuJyEb4bZH2u3OnzTORMOGmeinJylL/g6v7iA5TnaWnKWjauQkZd09Ia/nIHlO111h11+riy6uqIZH2u3OnzTOTKmMxYKS/P9P6TJ+1TDnuQUdcDB0zv37fP9jFlxZVVV0ck63PlTp9nIlfGZMZK3bqZ3h8WZp9y2IOMuvbvb3p/VJTtY8qKK6uujkjW58qdPs9ErozJjJUiIpTOgR4e+ts9PJTtrjQKQkZdY2OV/iIN8fdXb6SPjLiy6uqIZH2u3OnzTOTKmMw0Qmqq0jmwruhoZburkVHXrCzDL/naET5qkhFXVl0dkazPlTt9nolcFUczNUFurnJP3R3mpZBR1x07lH4j9p57RUZcWXV1RLI+V+70eSZyBtZ8fzOZISIiIofDodlERETkNpjMEBERkVPjcgZEZBTXKyIiZ8ArM0RkoKREmeY/MhIYOVIZwhwXB1y4ILtkRESGmMwQkYH4eGW9orrS04EJE+SUh4jIFCYzRKSH6xURkbNhMkMWyckBtm617xfZ4sXA0KHA0qX2i2lJXDXKZe6cs2cDffoA8+bZLqYxstcrWrECmDQJSElRN44jxJX1HpfxeSZSlXBxpaWlAoAoLS2VXRSnVFwsRGysEMCNR2ysECUl6sXcuVMIjUY/pkYjRGamejEtiatGucydc8MG/X21j82bm1xdo06caDhm7SMnR524hw4JodXqx9JqhfjhB3XiyYwr6z0u4/NM1FjWfH9z0jwyKS5O6StR95aDh4cy3fu2berEbNZM+TVbn0YD1NSoE9OSuGqUy9w5NRrjr1XzkysjrpcXUFVluF2rBSor1YkpK66s97iMzzNRY3HSPLIJGX0nFi82/mUphHqX483FHT7c9uUyF/P//T/Tr1frltPixab3q/EzWLGi4YQCULardetHRlxZ73H2hSJXxmSGjJLRd2LnTtP7t2+3fUxL4h48aHp/Y8plLmZ2tun9O3ZYH9MSMn4GGRmm95srkzPFlfUel90XikhNTGbIqG7dTO8PC7N9zGHDTO8fPtz2MS2Je9ddpvc3plzmYvbta3q/WgtSyvgZDB5ser+5MjlTXFnvcRmfZyJ7YZ8ZMol9Ztyzz4yMnwH7zLDPDFFd7DNDNpOaqvyiqys6WtmulowMwy9xjcb8LQG146pRLnPn3Ly54dcZ224rMn4GBw8qCURdWq35W3zOGFfWe1zG55nIHnhlhiySm6vcU7fnGj1Llyr9B4YPB557zj4xLYmrRrnMnXPePKWPTEwM8Oabtolpi3KpISVF6VcybBiQmGifmLLiynqPy/g8E1nLmu9vqQtNLlu2DMuWLcPp06cBAD169MCLL76IESNGAAASExOxatUqvdf0798f+/fvt3dR3V54uP1/6T33nH1/wVsat0MH4KabgMBA+8WsPwLFXmTE/e035VFYaN+4iYn2TZ4Aee9xGZ9nWYuWcrFU9yD1yszXX38NDw8PhP2v59mqVavwxhtv4IcffkCPHj2QmJiI//73v0hOTta9xsvLC+3atbM4Bq/MkK0cPgxERen3r6i9HWGus25jbdwIPPCA4fbNm4FRo9SJCQC7dim3H+r+dqi9DTJwoOvEJPWVlChrfaWl3dgWG6vc2mrb1vXiku1Y8/3tcLeZ2rVrhzfeeANTpkxBYmIiLl68iI0bNzb6fExmyFZkdBR1pw7AsjrFkrpkdTpmZ2fn55QdgKurq7FmzRpcuXIFUVFRuu0ZGRno2LEjIiIi8MQTT6CoqMjkeSoqKlBWVqb3IGoqGZOrzZ5ter+ak+bZe1I3WRPJkbpkTdTHCQLdj/Rk5ujRo2jdujW8vb0xbdo0bNiwAd27dwcAjBgxAv/+97+xa9cuvPXWW8jKysLQoUNRUVFh9HxJSUnw8/PTPUJCQuxVFXJhMiZX27XL9H5XmjRP1kRypC5ZE/VxgkD3Iz2ZiYyMRHZ2Nvbv348nn3wSCQkJ+OWXXwAA48aNw6hRo9CzZ0/cd9992Lp1K3JycvDNN98YPd+CBQtQWlqqe5w7d85eVSEXJmNytaFDTe93pUnzZE0kR+qSNVEfJwh0Pw7XZyY6OhrdunXDxx9/3OD+8PBwPP7445g/f75F52OfGbIV9plhnxmyHvvMUGM5ZZ+ZWkIIo7eRiouLce7cOQQFBdm5VERyJldzp0nzZE0kR+qSNVEfJwh0L1KTmYULF+Lbb7/F6dOncfToUTz//PPIyMjAI488gsuXL2PevHnYt28fTp8+jYyMDNx3331o3749HmhorKqbWLECmDRJvVWEHS2uDDk5wNathp0E+/ZVrsAkJgJduwKTJyvPbTEs21jMUaOUqxWTJ9+IKYS6w7IBZSh0TQ0wdy7Qu7fyb02NukOkZcSsKzERCA0FpkyxTzxZMQGlc3mfPup1Iq+rbVvlSsjixcCQIUBSkvJc7eHRsuKSJEKiyZMniy5duggvLy/RoUMHMWzYMLF9+3YhhBBXr14Vw4cPFx06dBBarVZ07txZJCQkiLNnz1oVo7S0VAAQpaWlalTBbg4dEkKrFUL5KlMeWq0QP/zgmnFlKC4WIjZWv66xsUKUlCj71WgLczHN7VeLjLiy6rp6tX7M2seaNa4VUwghNmxoOO7mzerF5O8uaixrvr8drs+MrblKnxl3WoRPFnP32NVoC3Mx3am/gay6yuiXJKsvlIy4/N1FjeXUfWbIkIw5TmTGlcHcvBRLlti+LczF3L7dfebokDUviLnlC9S4/SMjJiBn3iL+7iJ7YTLjBGTMcSIzrgzm5qVITze9vzFtYS6muSXIXGmODlnzgmRmmt5vbq4fZ4lpyXnVmLeIv7vIXpjMOAEZc5zIjCuDuXkp6o+KqK8xbWEu5t13m97vSnN0yJoXZNAg0/vNzfXjLDEtOa8a8xbxdxfZC/vMOAned1Yf+8xYXi5XiQmwz4zacfm7ixqLfWZckIw5TmTGlcHcvBRqtIW5mO40R4esuq5ZY912Z40JyJm3iL+7yB54ZcbJpKQo93uHDTPfkdAV4sqQm6v00QgLA8LDDffHxABZWUD//krnVFvYvl3pIxMV1fDl/sceU/oBDB2qdG60lwcfBL7/Hrj3XuCrr+wTU432tURYGJCfr8z7Yq+1e3r0UDo/R0QAx47ZJyYA3H478MsvSvzDh+0Tk7+7yFpWfX+rOkjcAbjKPDMknxpzdJibW0XGvCBCyJkH5aOPGo75z3+qF1MIOXWV9XOVFZeoMTjPTB2udmWG5FGjv4G5fiLu1LeCdVU3psy4RI3BPjNENqbGHB3m5lYxd0lcranoZcyDYm6Fkj//2fYxATl1lTHfi8y4RPbAZIbIAmrM0WFubhVz85GoMS+IJXHVmAdl717T+7/91vYxATl1lTHfi8y4RPbAZIbIAmrM0WFubhVz85GoMS+IJXHVmAdlwADT+++91/YxATl1lTHfi8y4RPbQpD4zJ0+eRF5eHgYOHIgWLVpACAGNqZuyErDPDNkK+8yoG5d1VTemzLhEjaF6n5ni4mJER0cjIiICI0eOREFBAQDg8ccfx9y5cxtzSiKHp8YcHebmVpExLwggZx6Uf/7Tuu22IqOusn6usuISqa4xw6UmTZokYmNjxblz50Tr1q1FXl6eEEKItLQ00b1798acUjUcmm0bJ04IsWWLEDk5rh3TEnPnCtG7t/KvreTkmK6rGjEtMXmyEF27Kv/ay9ixQnTsqPxrTzLqKiOmEPLeT0TWUH1odmBgINLS0tCnTx/4+Pjgxx9/xM0334z8/Hz06tULly9ftn3W1Ui8zdQ0JSVAfLz+5GWxscqVg7ZtXScmkT3xPU5knuq3ma5cuYKWLVsabP/jjz/g7e3dmFOSg4qPN1wxOj0dmDDBtWIS2RPf40S21ahkZuDAgVi9erXuuUajQU1NDd544w0MGTLEZoUjuczNg5Kb6xoxieyJ73Ei2/NszIveeOMNDB48GIcOHUJlZSWeffZZHDt2DCUlJfj+++9tXUaSxNw8KCdPNrx2kbPFtLWcHKUextZ2IvfmCu9xIkfTqCsz3bt3x08//YS77roLMTExuHLlCh588EH88MMP6GZu8gxyGuZ+lGFhrhHTVkpKlKHWkZHAyJHK4oFxccCFC7JLRo7Emd/jRI6KazORSebmQXGVmLbgrOUm++N7hcg81TsAJycn48svvzTY/uWXX2LVqlWNOSU5KHPzoLhKzKZiPwiyhjO+x4kcWaP6zCxduhTLly832N6xY0f85S9/QUJCQpMLRo6hbVvlL8XcXOVevj36gciI2VTsB0HWcMb3OJEja1Qyc+bMGYSGhhps79KlC86ePdvkQpHjCQ+3/y9bGTEB8x14c4pzkFeSh7B2YQj3Vw5Qux9EQzHtQUZcd6rrqdIcZF3Mg2dZGMJhv7qmHcrBgdw8REWEIeYOZlHk/BqVzHTs2BE//fQTunbtqrf9xx9/hL+/vy3KRWR35iYyK7lWgvh18UjLu3FAbLdYpI5NRUREW8TGGu8H0dikzFTMti3Um11NRlx3qmve+RL0fz0exW3/FzMH8P9XLLLmpyI0SL26yopLpLZG9ZkZP348/vrXv2L37t2orq5GdXU1du3ahZkzZ2L8+PG2LiORXZibyCx+XTzST+kfkH4qHRPWKQeY6geRkwNs3Xqj70z950bLZCamWmTEdae69n89HsV++jGL/dJx52vq1lVWXCK1NWo0U2VlJSZNmoQvv/wSnp7KxZ2amho8+uijWL58Oby8vGxe0MbiaCayRE6OMqTamLRDOYjdbPyAnOk5ulsTdftB+PsbXu3x9weKi288NzaNfU5xDiI/sCymLcmI6051TTuUg7hvjMfc/qccVW79yIpL1Fiqj2by8vLCF198gf/7v//Dv//9b6xfvx55eXlYuXKlQyUyRJYy14F3f47pA06WnNT9PzwcGDFC+behqz11ExnA+DT2eSWWx7QlGXHdqa4Hck3H3JejTl1lxSWyh0b1makVERGBiIgIW5WFSBpzHXjvjugG5BjfH9bOsIdv7XBtc+oO367bt6ZbO9OFaiimLciI60517R9u+r0UFaFOXWXFJbIHi5OZOXPm4O9//ztatWqFOXPmmDz27bffbnLBiOwpIgImO/AOvyMCscdjkX4qHdXixgEeGg9E3xzd4K0Ic1d76qs/fDvCPwKx3ayLaQsy4rpTXWP7RcD/s1il70qzOm+2Gg/4l0ardqtHVlwie7D4NtMPP/yAqqoqAMCRI0fwww8/NPjIzs5Wq6xEqjI3kVnq2FRE36x/QPTN0Ugd2/BMZ9au7NHQ8G1rY9qKjLjuVNes+anwL9WP6V8ajaz56tZVVlwitUldzmDZsmVYtmwZTp8+DQDo0aMHXnzxRYwYMQIAIITAyy+/jE8++QQXLlxA//798eGHH6JHjx4Wx2AHYOe1YgWQkQEMGwYkJtov7pIlyhWa4cOB554z3L9yYy52/3gSw24LQ+L9hn/NpqUBBw4AUVHAW28ZXu2pz8MDuOkmoE0bICYGePNN/XPExABLluci/YeTGN4vDM89Yb+/oGe/motd2ScRc3sY3lxon7iy6iojbtLHudhxxP513XE4F/tyTnKeGXJoVn1/CytVVVUJDw8PcfToUWtfamDTpk3im2++ESdOnBAnTpwQCxcuFFqtVvz8889CCCGWLl0qfHx8xLp168TRo0fFuHHjRFBQkCgrK7M4RmlpqQAgSktLm1xeso9Dh4TQaoUAbjy0WiF++EFu3JMnhfD319/v7y/EqVPG9/v6CuHpqb+tqQ97tMXOnUJoNPpxNRohMjPVi+moP3c1mHsvEZF139+NujLTrVs3rF+/Hn369GlcumVCu3bt8MYbb2Dy5MkIDg7GrFmzMH/+fABARUUFAgIC8Nprr2Hq1KkWnY9XZpyPlxfwvzuaerRaoLJSXtz27Q1HIgHKUOs//jC+Xw1qt0WzZspXbH0aDVBTo05MR/25q8Hce4mI7DA0+4UXXsCCBQtQUlLSqAI2pLq6GmvWrMGVK1cQFRWF/Px8FBYWYvjw4bpjvL29MWjQIOzdu9foeSoqKlBWVqb3IOexYkXDXyyAsj0lRU7cuXONJyrFxUBSkv0SmdoyqdUWixc3nMgAyvalS20f01F/7mrETUsz/V7ascP2MYlcXaOSmffeew/ffvstgoODERkZidtvv13vYY2jR4+idevW8Pb2xrRp07BhwwZ0794dhYWFAICAgAC94wMCAnT7GpKUlAQ/Pz/dIyQkxPoKkjQZGab379wpJ279uWLqk/EFpFZbmDvv9u22j+moP3c14h44YHr/vn22j0nk6ho1z8yYMWOg0WjQiDtUBiIjI5GdnY2LFy9i3bp1SEhIQGZmpm6/RqPRO14IYbCtrgULFugNHS8rK2NC40QGDwY++8z4/mHD5MSNjgZ++sn4/pgYYPdumxfLJLXaYtgw03Wpc7HUZhz1565G3P79Te+PirJ9TCJXZ1WfmatXr+KZZ57Bxo0bUVVVhWHDhuH9999H+/btbVag6OhodOvWDfPnz0e3bt1w5MgR3Hbbbbr9o0ePRps2bbBq1SqLzsc+M87HUftOsM8M+8zYCvvMEJmnWp+Zl156CSkpKRg1ahQmTJiA9PR0PPnkk00qbH1CCFRUVCA0NBSBgYHYUef6fWVlJTIzMzFgwACbxiTHcvCg8kVSl1arbJcZNytL+bKpy99f2W5sv5+f4ZpLJi4sWsQebZGRYVhOjcb8bZmmcNSfuxrMvZeIyErWDJO6+eabRWpqqu75gQMHhKenp7h+/bo1p9FZsGCB2LNnj8jPzxc//fSTWLhwoWjWrJnYvn27EEIZmu3n5yfWr18vjh49KiZMmMCh2W4kOVmIiROVfx0p7vbtQrz8svKvpfvrb6v/fO5cIXr3Vv5taL+stkhKEmLIEOVfe3HUn7sazL2XiNyZakOzvby8kJ+fj5tuukm3rUWLFsjJyWlUv5QpU6Zg586dKCgogJ+fH3r37o358+cjJiamNtHCyy+/jI8//lhv0ryePXtaHIO3mYiIiJyPNd/fViUzHh4eKCwsRIcOHXTbfHx88NNPPyE0NLTxJVYRkxlyBjk5ylpOYWHK+kz1nxMRuRtrvr+tGs0khEBiYiK8vb1128rLyzFt2jS0atVKt239+vVWFpnIPZWUAPHx+qtr+/vrdw6NjVXWh6rf94aIiBRWJTMJCQkG2yZOnGizwhC5m/h4wzls6o9ySU8HJkwAtm2zX7mIiJyJ1IUm7YG3mchR5eQAkZHWHc9bTkTkLlRfzoCImi4vz7rjT55UpxxERM6OyQyRJN26WXd8WJg65SAicnZMZoislJYGvPJK09djiohQOvd6eJg+zsMDuOMOpROwvdeAslVdHT2mrLiy6krkahq1NhORO8rLU9bVqdtBt3bW1sbOTJCaqnTuNTWaqVkz4PBh5WGLmJZQo66OGFNWXFl1JXJV7ABMZCE119PJzVX6xNTOK1P7PD4euHhRnZimyFg7SNZ6Re5UVyJnwg7ARDaWlmZ8Ecni4qbfJggPB0aMuDFaKTxcuSLTUCJjq5jGqF1XR4kpK66suhK5MiYzRBY4cMD0/n37XCOmrLisq7oxiVwd+8wQWaB/f9P7o6JM72/McgVNjdlYMuKyrurGJHJ17DNDZKHG9HNo6nIF7EfCuhK5K/aZIVJBVpbyZVNX7QgUY6xZrsBWMW1BRlzWVf26ErkqXpkhstKOHUq/hqgoICbG+HG2XK7A0pi2JiMu60pEgHXf30xmiFSydSswcqTlx2/ZooxoIiIi3mYicghcroCIyD6YzBCpxJrlCmJjuSI2EVFjMZkhUlFqKhAdrb+tfsfP6GjlOCIiahzOM0OkorZtgW3bjC9XYGqeGSIisgyTGSI7CA/XT1rqPyciosbjbSYiIiJyakxmiIiIyKkxmSEiIiKnxmSGiIiInBo7ABM5gMasqk1ERAomM0QSNXVVbSIi4m0mIqmauqo2ERExmSGSJidHuSJTXW36uOpq5bjcXPuUi4jI2TCZIZIkL8+640+eVKccRETOjskMkSRcVZuIyDaYzBDZWE4OsHWr+dtCXFWbiMg2pCYzSUlJuPPOO+Hj44OOHTtizJgxOHHihN4xiYmJ0Gg0eo+7775bUomJjCspAeLigMhIYORIJVmJiwMuXDD+Gq6qTUTUdFKTmczMTDz99NPYv38/duzYgevXr2P48OG4cuWK3nFxcXEoKCjQPbZs2SKpxETGNTQyydxIpNpVtXNygC1blH//+EP/+bZtHJZNRGSKRgghZBei1u+//46OHTsiMzMTAwcOBKBcmbl48SI2btxo0TkqKipQUVGhe15WVoaQkBCUlpbC19dXjWITISdHuSJjaj9vExERWa6srAx+fn4WfX87VJ+Z0tJSAEC7du30tmdkZKBjx46IiIjAE088gaKiIqPnSEpKgp+fn+4REhKiapmJAPMjkzgSiYhIPQ5zZUYIgdGjR+PChQv49ttvddu/+OILtG7dGl26dEF+fj7+9re/4fr16zh8+DC8vb0NzsMrMyQDr8wQEdmWNVdmHGY5g+nTp+Onn37Cd999p7d93Lhxuv/37NkT/fr1Q5cuXfDNN9/gwQcfNDiPt7d3g0kOkZpqRyalp+tPgufhoXTgZSJDRKQeh7jNNGPGDGzatAm7d+9Gp06dTB4bFBSELl26IJfToZKDaWhkEkciERGpT+qVGSEEZsyYgQ0bNiAjIwOhoaFmX1NcXIxz584hKCjIDiUkslztyKTcXKWPDFe8JiKyD6nJzNNPP43PP/8c//nPf+Dj44PCwkIAgJ+fH1q0aIHLly9j0aJFGDt2LIKCgnD69GksXLgQ7du3xwMPPCCz6ERGhYdbn8SkpQEHDgBRUUBMjOFzIiIyTmoHYI1G0+D25ORkJCYm4tq1axgzZgx++OEHXLx4EUFBQRgyZAj+/ve/WzxKyZoORET2lpcH9O+vv1K2RgPU/VT6+wNZWYAFFy6JiFyGNd/fDjOaSS1MZsiRtW+vn8gY4++vTKZHROQunHaeGSJ3kpZmWSIDKMft2KFueYiInBWTGSJJDhyw7vh9+9QpBxGRs2MyQyRJ//7WHR8VpU45iIicHZMZIkliYw1XyDbG35+jmoiIjGEyQyRRVpZhQlN/kF/taCYiImqYwyxnQOSOQkOVUUo7dih9Ymrnlan/nIiIjOPQbCIiInI4HJpNREREboPJDBERETk1JjNERETk1JjMEBERkVNjMkNEREROjckMEREROTXOM0NUX04OkJcHhIUB4eGuG5OIyEXwygxRrZISIC4OiIwERo4EIiKU5xcuuFZMIiIXw2SGqFZ8PJCerr8tPR2YMMG1YhIRuRgmM0SAcpsnLQ2ortbfXl2tbM/NdY2YREQuiMkMEaD0VzHl5EnXiElE5IKYzBABQLdupveHhblGTCIiF8RkhghQOt7GxgIeHvrbPTyU7WqMMJIRk4jIBTGZIaqVmgpER+tvi45WtrtSTCIiF6MRQgjZhVCTNUuIEwFQOt6ePGnfOV9kxCQicmDWfH9z0jyi+sLD7Z9QyIhJROQieJuJiIiInBqTGSIiInJqTGaIiIjIqbHPDBEZ506LbrKurhdTVlx3qqujEC6utLRUABClpaWyi0LkPIqLhYiNFQK48YiNFaKkxLViyorLurKurhBXZdZ8fzOZISJDsbFCeHjo/3L08FC2u1JMWXFZV9bVFeKqzJrvb84zQ0T6cnKAyEjT+219CVtGTFlxWVd1Y8qK6051tRNrvr+ldgBOSkrCnXfeCR8fH3Ts2BFjxozBiRMn9I4RQmDRokUIDg5GixYtMHjwYBw7dkxSid1YTg6wdat9V3KWEVNmXEfhTotusq7qxmVd1Y0pM66DkZrMZGZm4umnn8b+/fuxY8cOXL9+HcOHD8eVK1d0x7z++ut4++238cEHHyArKwuBgYGIiYnBpUuXJJbcjZSUAHFxSuY/cqSynlBcHHDhgmvFlBnX0bjTopusq7pxWVd1Y8qM62jUvudljaKiIgFAZGZmCiGEqKmpEYGBgWLp0qW6Y8rLy4Wfn59Yvny5Redkn5km4n1np7/v3Cj8ubOuzhxTVlx3qqsdOG0H4NzcXAFAHD16VAghRF5engAgjhw5onfc/fffLx599NEGz1FeXi5KS0t1j3PnzjGZaawTJ/Q/HPUfOTmuEVNmXEdVUmL/0REyYsqKy7qyrq4QV2XWJDMOM8+MEAJz5szBPffcg549ewIACgsLAQABAQF6xwYEBODMmTMNnicpKQkvv/yyuoV1F5bci7V1xzIZMWXGdVRt2wLbttl3AUwZMWXFZV3VjSkrrjvV1cE4TDIzffp0/PTTT/juu+8M9mk0Gr3nQgiDbbUWLFiAOXPm6J6XlZUhJCTEtoV1F7zvrH5cR+dOi26yrq4XU1Zcd6qrg3CI5QxmzJiBTZs2Yffu3ejUqZNue2BgIIAbV2hqFRUVGVytqeXt7Q1fX1+9BzVSRAQQGwt4eOhv9/BQtqvxoZERU2ZcIiJqMqnJjBAC06dPx/r167Fr1y6Ehobq7Q8NDUVgYCB27Nih21ZZWYnMzEwMGDDA3sV1T6mpQHS0/rboaGW7K8WUGZeIiJpE6m2mp59+Gp9//jn+85//wMfHR3cFxs/PDy1atIBGo8GsWbOwZMkShIeHIzw8HEuWLEHLli0RHx8vs+juo/Ze7PbtwP79QFQUEBPjejGtiavG+ifmzulOa72wrq4bl0gtavdGNgVAg4/k5GTdMTU1NeKll14SgYGBwtvbWwwcOFA32skSHJrdRFzf5EZcNcolI6YtyuUqMWXFdae6EjWS0w7NVgOTmSbiXA034qpRLhkxbVEuV4kpK6471ZWokbg2Ux1cm6kJuL7JDWlpSkdgW5ZLRkxblMudfu6sK5E0TrM2Ezk4rm9yw/79pvc3plwyYlqCP3d147pTXYnshMkMGcd5Zm64+27T+xtTLhkxLcGfu7px3amuRHbCZIaM4zwzN+IOH277csmIaYtyudPPnXUlcg6q9+CRTNUOwCdOCLFli33X7bF3TJnrm4RDiDgIEWbn9U2MxVWjLcyd053WemFdXbOuRI3klGszOZWSEiA+XumgWSs2VplcrW1b14kJKL/u7K4EePUQ0K/OpkOHAM1FACrW1VxcNdY/Mde+7rTWC+uqbkyZcYlUxtFMjREXB6SnA9XVN7Z5eCizxW7bZpsYjhBTVtxD7YG+xfpTOl4HkO0P9PtDnZiy4sr6uRIROThrvr+ZzFiLQznVjZufBoTGGd9/ejvQVYXZgGXE5VBZIiKjODRbTRzKqW7cPw6Y3v/7PtvHlBWXQ2WJiGyCyYy1OJRT3bjt+5ve3yHK9jFlxeVQWSIim2AyYy0O5VQ3bmgscMhf6atS13Uo29W4xSQrLofKEhHZBJOZxkhNVTpo1hUdrWx3pZiy4oZlKZ1u68r2V7arSUZcWT9XIiIXwqHZjSGjz7Q7DeWs8QNe6AfkpQFhAE4C6NYPSG3jenFdu/89EZFdcDRTY3A4rbrcaRg630tERA3i0Ow6XGJotjtxp2HofC81LCdHGell7wndZMRlXYmM4tBsNXE4rbrcaRg630v6SkqUK1WRkcDIkUoH6bg44MIF14vLuqpfV3IrTGasxeG06nKnYeh8L+mLj1duudWVng5MmOB6cVlX9etKboXJjLU4nFZd7jQMne+lG3JylHXH6vYdApTnaWlKB3RXicu6ql9XcjtMZhqDw2nV5U7D0PleUvD2oroxZcXlrVSyEw7NbgyZK8+6Q4fB2vZdskS5HD18OPDcc+rGrBt35Upg925g2DAgMdHwOFu2RW3M7duB/fuBqCggpoEJ+iT93BfPLsLOXc0wPEbguTc7qBeozi23NMTgAPojCvsQg53KRt5edM64vJVK9iJcXGlpqQAgSktLZRelaYqLhYiNFUKZmUR5xMYKUVLienEPHRJCq9WPqdUK8cMP6sUUQoiTJ4Xw99eP6+8vxKlTyn412sLcOSX93HduKBUaVOuF1aBaZG6+qFrMk/ckCH8U6Tc/isSpex9VLaYQQmlPDw/9NvbwULa7UkxZcWXVlZyeNd/fTGachTv98qufyNRNaNRUP5Gpm9AIoU5bmDunpJ+7ksjU1GuKGqFBtWox/ds1HNO/nXoxhRBKYmjvhFFGTFlxZdWVnB6TmTpcIpk5caLhL9naR06O68T95z9Nx0xOtn1MIYTYts103JUrbd8W5to3LU3Kz/3VWf81GTZpbpHNY5pr/u3bbR7yBllXPYVQfoZbtqj3GXakuLLqSk7Lmu9vdgB2Bu7UYTAjw/T+nTttHxMADhwwvX/XLtP7G9MW5tp3/37bx7TAzl2mfy1s36GxeUxzzb9vn81D3iBz6HB4ODBihP1HrsmIK6uu5BaYzDgDd+owOHiw6f3Dhtk+JgD07296/9Chpvc3pi3qtG8OwrEVcchFnfPcfbftY1pg2G3FJvcP71di85j9bzpncn9U519tHhMAhw4TuQgmM87AneZemTIF0Gob3qfVNjy6yBZiY03HfewxwN+/4f3+/o1ri4gIlAwZizhsQyRyMBJbEYFcxGEbLgwdq4zi0hi5CqLRqPZzf37cKWhQDUDU2yOgQTWe+7OZK0qNEBv8M/zxe4Mx/fE7YgKO2jwmAA4dJnIRTGachTvNvXLwoGFiodUq29WSkwNUVTW8r6pKGTpdbOSKRXFxo/+Cj9ekIl2j377pmmhMQCqwYoXSg6MhQgApKY2KadaRI8jAQGhQo7dZgxpkYCDw44+qxMzCnfDHH3qb/fEHsnCnOjEB+UOHZ88G+vQB5s1TN059ixcrVxuXLrVfzBUrgEmT1HvfOlLctDTglVeAHTvsF1NWXFl1rc8OfXikcokOwHW5U4fB5GQhJk5Ur9NvXVu2mO6B+vLLpvdv2WJ1SLP9q++fa/qAiRNVaAihnPd/MZLwjBiCHSIJz6gbt07M7RgmXsYLYjuGqV9XIYTo27fh9u3bV72YGzY0HHPzZvViCiHEzp1CaDT6MTUaITIz1Yspa6oFGXHNTe/gSnHtEJOjmepwuWSG1CFhZJG5/GnLTDNDfNRK8mSMKJM1ik0I03FdKaYQholM3YRGLbKmWpAR19z0Dq4U1w4xOZqJyFrm+gedM91BFd9/b3VIs3c4erc0fcBNN1kd0yKdOtk/royYgHKbxxQ1bv/IiAkot5aEaHifEOrcclqxwvTtW7Vu/ciIm5Zm+la0WrdhZMSVVVcTpCYze/bswX333Yfg4GBoNBps3LhRb39iYiI0Go3e425zIzyIGstU/yAVhoyb7V/9627TJ1BrvLKMcdKyxmabG3Kvxi9lGTEB8+/R7dttH1PWVAsy4sp6D7vT59UEqcnMlStX0KdPH3zwwQdGj4mLi0NBQYHusWXLFjuWkNxK7TpJOTnAli3Kv9u2KdtVGjJusn+1ueHiUVGNimmWjLiy6mpuyH1D62Q5Y0zA/Ht0+HDbx5Q11YKMuPy8qhvTDI0Qxq472pdGo8GGDRswZswY3bbExERcvHjR4IqNNcrKyuDn54fS0lL4+vo2vaDkvry8Gr50rdUClZVNOrXRNUvbt2/4cq6/P/DHH4bbbUVGXFl1NTb8HTB+W8YZYwJAs2YNn1+jAWpqDLfbgoqfG4eLy8+rTWNa8/3t8H1mMjIy0LFjR0REROCJJ55AUVGRyeMrKipQVlam9yCyCRWHjBudHDUry3B+G39/ZbuaZMSVVdfNm63b3lQ5Oab3qzlRX0aGYSKl0Zi/LdMUMqZakBWXn1f162qEQ1+Z+eKLL9C6dWt06dIF+fn5+Nvf/obr16/j8OHD8Pb2bvA8ixYtwssvv2ywnVdmyGJpaco94aiohi/5p6Qo99yHDbPdJH7mYu7YodyHNrZfLTLiyqrrvHlK7JgY4M031YuzdSswcqTx/Vu2KJmtmpYuVfrIDB8OPPecurFqqfG5cdS4/LzahDVXZhw6mamvoKAAXbp0wZo1a/Dggw82eExFRQUqKip0z8vKyhASEsJkhszLy1PuBde9dFr7l0ZoqOvEJLlycoDISNP7uX4RkWvdZqorKCgIXbp0Qa6Jy7De3t7w9fXVexBZpH5SASjP77zTtWKSXLKWJyFyYU6VzBQXF+PcuXMICgqSXRRyNZyrgexJ1vIkRC7KU2bwy5cv42Sdhdzy8/ORnZ2Ndu3aoV27dli0aBHGjh2LoKAgnD59GgsXLkT79u3xwAMPSCw1uSRL5k2ovR9srn+LI8e0loy47lDX2mkAVq4Edu+2fz8S/lxdL6asuDk5yu1yg6GYdmazeYcbYffu3QKAwSMhIUFcvXpVDB8+XHTo0EFotVrRuXNnkZCQIM6ePWtVDC5nQBbZZmbpgO3bbb8WiYyYlnLRtV4cJi7ryro6e9ziYiFiY/VjxsYKUVJisxDWfH87TAdgtXCeGbKYuXkT1JhXQUZMW5TLVWLKisu6sq7OHjcuDkhPB6qrb2zz8FBul27bZpMQLtsBmEhVpuZNUKt/i4yY5rhT/yHWVd24rKu6MWXFzclR4tZNZADleVqaunMlGcFkhqhWaKjyV8z27cDLLyv//vGHsl2ttUj8/IB+/fS39esHtGnDtV7UjikrLuuqbkxZcd2prnl5pvfX6QtrL1I7ABM5pJgYw85zaq1FEh+vXKqtKz0dmDDB/OrKXOvFOeOyrurGlBXXnerarZvp/WFhto9pBq/MEFnC3AR2Xbtaf05zl2pvvtnwFlQtf3/1RiuoUVdHjCkrrqy6xsYC7do1vK9dO3XeTzJiyoors672/j3hgHMlMZkhsoQal1UtOaeM9U9kXEKWddnaneoKAD16WLfdWWPKiiurrjJ+TzjYXEm8zURkCTUuq1pyztp+PPZcc0XGJWRZl63dqa45OcC33za879tvlU6btv6LWkZMWXFl1RWQ83uidq6k3FwlAZc8zwyvzBBZQo3LqtacMyYGePFF+0yEJeMSsqzL1u5UV3e6CuVOda3Lnr8naoWHKwujSl6Gg8kMkaXUuKzqYJdqdWSUS1ZbuEtd3ekqlDvVlQA40KrZauGkeWRzalxWdZBLtQZklEtWW7hDXe0w0ZlDxJQVV1ZdXZQ1399MZoiI3MWFC8qw/7S0G9tiY5UrQm3buk5MWXFl1dVFMZmpg8kMEVE97nAVSmZcR73S6mSYzNTBZIaIiMj5cG0mIiIichtMZoiIiMipMZkhIiIip8YZgInIseTkKBOQsfMkEVmIV2aIyDGUlCjzdERGAiNHKjPlxsUpw12JiExgMkNEjiE+XplwrK70dGXeDiIiE5jMEJF8OTnKRGN1Z04FlOdpacq8HURERjCZISL5HGGRPiJyWkxmiEg+LtJHRE3AZIaI5IuIUNaw8fDQ3+7hoWznqCYiMoHJDBE5htRUZXXhuqKjle1ERCZwnhkicgxt2wLbtslbpC8tDThwAIiKAmJiXDcmwLl8yOUwmSEixxIebt8v2Lw8oH9/oLj4xjZ/fyArCwgNdZ2YgDKXT3y8kkTVio1Vrn61bateXCKV8TYTEbm3+kkFoDy/807XiglwLh9yWUxmiMh9paUZJhW1iouBHTtcIybAuXzIpTGZISL3deCA6f379rlGTIBz+ZBLYzJDRO6rf3/T+6OiXCMmwLl8yKUxmSEi9xUbq3S8bYi/vzojjGTEBDiXD7k0JjNE5N6ysgyTi9qRRa4UE+BcPuSypCYze/bswX333Yfg4GBoNBps3LhRb78QAosWLUJwcDBatGiBwYMH49ixY3IKS0SuKTQU+OMPYPt24OWXlX//+EPdIdIyYgI35vLJyQG2bFH+3baNw7LJ6UmdZ+bKlSvo06cPHnvsMYwdO9Zg/+uvv463334bKSkpiIiIwKuvvoqYmBicOHECPj4+EkpMRC4rJsa+E9fJignYfy4fIpVJTWZGjBiBESNGNLhPCIF3330Xzz//PB588EEAwKpVqxAQEIDPP/8cU6dObfB1FRUVqKio0D0vKyuzfcGJiIjIYThsn5n8/HwUFhZi+PDhum3e3t4YNGgQ9u7da/R1SUlJ8PPz0z1CQkLsUVwiIiKSxGGTmcLCQgBAQECA3vaAgADdvoYsWLAApaWluse5c+dULScRERHJ5fBrM2k0Gr3nQgiDbXV5e3vD29tb7WIRERGRg3DYKzOBgYEAYHAVpqioyOBqDREREbkvh01mQkNDERgYiB111imprKxEZmYmBgwYILFkRERE5Eik3ma6fPkyTtZZDyQ/Px/Z2dlo164dOnfujFmzZmHJkiUIDw9HeHg4lixZgpYtWyI+Pl5iqYmIiMiRSE1mDh06hCFDhuiez5kzBwCQkJCAlJQUPPvss7h27RqeeuopXLhwAf3798f27ds5xwwRERHpaIQQQnYh1FRWVgY/Pz+UlpbC19dXdnGIiIjIAtZ8fztsnxkiIiIiSzj80Oymqr3wxJmAiYiInEft97YlN5BcPpm5dOkSAHAmYCIiIid06dIl+Pn5mTzG5fvM1NTU4Pz58/Dx8TE52Z49lJWVISQkBOfOnWP/HRtge9oe29S22J62xfa0PUduUyEELl26hODgYDRrZrpXjMtfmWnWrBk6deokuxh6fH19He5N48zYnrbHNrUttqdtsT1tz1Hb1NwVmVrsAExEREROjckMEREROTUmM3bk7e2Nl156iQth2gjb0/bYprbF9rQttqftuUqbunwHYCIiInJtvDJDRERETo3JDBERETk1JjNERETk1JjMEBERkVNjMtNESUlJuPPOO+Hj44OOHTtizJgxOHHihNHjp06dCo1Gg3fffVdve0VFBWbMmIH27dujVatWuP/++/Hrr7+qXHrHY2l7Hj9+HPfffz/8/Pzg4+ODu+++G2fPntXtZ3sqLGnPy5cvY/r06ejUqRNatGiBW2+9FcuWLdM7hu15w7Jly9C7d2/dJGNRUVHYunWrbr8QAosWLUJwcDBatGiBwYMH49ixY3rnYHveYKo9q6qqMH/+fPTq1QutWrVCcHAwHn30UZw/f17vHGzPG8y9P+tyqe8jQU0SGxsrkpOTxc8//yyys7PFqFGjROfOncXly5cNjt2wYYPo06ePCA4OFu+8847evmnTpombbrpJ7NixQxw5ckQMGTJE9OnTR1y/ft1ONXEMlrTnyZMnRbt27cQzzzwjjhw5IvLy8sTmzZvFf//7X90xbE+FJe35+OOPi27duondu3eL/Px88fHHHwsPDw+xceNG3TFszxs2bdokvvnmG3HixAlx4sQJsXDhQqHVasXPP/8shBBi6dKlwsfHR6xbt04cPXpUjBs3TgQFBYmysjLdOdieN5hqz4sXL4ro6GjxxRdfiP/7v/8T+/btE/379xd33HGH3jnYnjeYe3/WcrXvIyYzNlZUVCQAiMzMTL3tv/76q7jpppvEzz//LLp06aL35rl48aLQarVizZo1um2//fabaNasmdi2bZu9iu6QGmrPcePGiYkTJxp9DdvTuIbas0ePHuKVV17RO+72228XL7zwghCC7WmJtm3bin/+85+ipqZGBAYGiqVLl+r2lZeXCz8/P7F8+XIhBNvTErXt2ZCDBw8KAOLMmTNCCLanJeq3pyt+H/E2k42VlpYCANq1a6fbVlNTg0mTJuGZZ55Bjx49DF5z+PBhVFVVYfjw4bptwcHB6NmzJ/bu3at+oR1Y/fasqanBN998g4iICMTGxqJjx47o378/Nm7cqHsN29O4ht6f99xzDzZt2oTffvsNQgjs3r0bOTk5iI2NBcD2NKW6uhpr1qzBlStXEBUVhfz8fBQWFuq1lbe3NwYNGqRrK7ancfXbsyGlpaXQaDRo06YNALanKQ21p6t+HzGZsSEhBObMmYN77rkHPXv21G1/7bXX4Onpib/+9a8Nvq6wsBBeXl5o27at3vaAgAAUFhaqWmZH1lB7FhUV4fLly1i6dCni4uKwfft2PPDAA3jwwQeRmZkJgO1pjLH353vvvYfu3bujU6dO8PLyQlxcHD766CPcc889ANieDTl69Chat24Nb29vTJs2DRs2bED37t117REQEKB3fN22YnsaMtae9ZWXl+O5555DfHy8blFEtqchU+3pqt9HLr9qtj1Nnz4dP/30E7777jvdtsOHD+Mf//gHjhw5Ao1GY9X5hBBWv8aVNNSeNTU1AIDRo0dj9uzZAIC+ffti7969WL58OQYNGmT0fGxPw/YElGRm//792LRpE7p06YI9e/bgqaeeQlBQEKKjo42ez53bMzIyEtnZ2bh48SLWrVuHhIQEXTINwKBdLGkrtqdhe9ZNaKqqqjB+/HjU1NTgo48+MntOtqdhe167ds1lv494ZcZGZsyYgU2bNmH37t3o1KmTbvu3336LoqIidO7cGZ6envD09MSZM2cwd+5cdO3aFQAQGBiIyspKXLhwQe+cRUVFBn/huQtj7dm+fXt4enoa/NV266236kYzsT0NGWvPa9euYeHChXj77bdx3333oXfv3pg+fTrGjRuHN998EwDbsyFeXl4ICwtDv379kJSUhD59+uAf//gHAgMDAcDgL9i6bcX2NGSsPWtVVVXh4YcfRn5+Pnbs2KG7KgOwPRtirD1d+fuIyUwTCSEwffp0rF+/Hrt27UJoaKje/kmTJuGnn35Cdna27hEcHIxnnnkGaWlpAIA77rgDWq0WO3bs0L2uoKAAP//8MwYMGGDX+shmrj29vLxw5513GgwvzsnJQZcuXQCwPesy155VVVWoqqpCs2b6vwo8PDx0V8HYnuYJIVBRUYHQ0FAEBgbqtVVlZSUyMzN1bcX2NK+2PYEbiUxubi7S09Ph7++vdyzb07za9nTp7yO7dzl2MU8++aTw8/MTGRkZoqCgQPe4evWq0dfU7z0uhDIUrlOnTiI9PV0cOXJEDB061OGHwqnBkvZcv3690Gq14pNPPhG5ubni/fffFx4eHuLbb7/VHcP2VFjSnoMGDRI9evQQu3fvFqdOnRLJycmiefPm4qOPPtIdw/a8YcGCBWLPnj0iPz9f/PTTT2LhwoWiWbNmYvv27UIIZWi2n5+fWL9+vTh69KiYMGFCg0Oz2Z4KU+1ZVVUl7r//ftGpUyeRnZ2t9x6uqKjQnYPteYO592d9rvJ9xGSmiQA0+EhOTjb6mobePNeuXRPTp08X7dq1Ey1atBB/+tOfxNmzZ9UtvAOytD1XrFghwsLCRPPmzUWfPn305kQRgu1Zy5L2LCgoEImJiSI4OFg0b95cREZGirfeekvU1NTojmF73jB58mTRpUsX4eXlJTp06CCGDRum90VRU1MjXnrpJREYGCi8vb3FwIEDxdGjR/XOwfa8wVR75ufnG30P7969W3cOtucN5t6f9bnK95FGCCHsdx2IiIiIyLbYZ4aIiIicGpMZIiIicmpMZoiIiMipMZkhIiIip8ZkhoiIiJwakxkiIiJyakxmiIiIyKkxmSEiIiKnxmSGiFzG999/j169ekGr1WLMmDGyi0NEdsJkhoislpiYCI1GA41GA61Wi4CAAMTExGDlypW6BSptafPmzRg8eDB8fHzQsmVL3HnnnUhJSTE4bs6cOejbty/y8/Mxe/ZsaLVafPfdd3rHXLlyBTfffDNmz55t83ISkRxMZoioUeLi4lBQUIDTp09j69atGDJkCGbOnIk//elPuH79us3ivP/++xg9ejQGDBiAAwcO4KeffsL48eMxbdo0zJs3T+/YvLw8DB06FJ06dcKgQYMwY8YMJCYm4sqVK7pjnn32WXh7eyMpKclmZaxVWVlp83MSkQVkLw5FRM4nISFBjB492mD7zp07BQDx6aefCiGEeOutt0TPnj1Fy5YtRadOncSTTz4pLl26JIQQ4vLly8LHx0d8+eWXeufYtGmTaNmypSgrKxNnz54VWq1WzJkzxyDWe++9JwCI/fv3N7ggYXJysrh27Zq49dZbxdNPPy2EEGLXrl1Cq9WKrKwsUVNTI1577TURGhoqmjdvLnr37q1XluvXr4vJkyeLrl27iubNm4uIiAjx7rvvNtgOS5YsEUFBQaJLly5NaVYiaiRemSEimxk6dCj69OmD9evXAwCaNWuG9957Dz///DNWrVqFXbt24dlnnwUAtGrVCuPHj0dycrLeOZKTk/HnP/8ZPj4++Oqrr1BVVWVwBQYApk6ditatWyM1NRUhISEoKCiAr68v3n33XRQUFGDcuHFo3rw5Vq9ejU8++QQbN27E5MmTsXDhQvTr1w8vvPACkpOTsWzZMhw7dgyzZ8/GxIkTkZmZCQCoqalBp06dsHbtWvzyyy948cUXsXDhQqxdu1avHDt37sTx48exY8cObN68WY1mJSJzZGdTROR8jF2ZEUKIcePGiVtvvbXBfWvXrhX+/v665wcOHBAeHh7it99+E0II8fvvvwutVisyMjKEEEJMmzZN+Pn5GS1H7969xYgRI3TP/fz8RHJyssFxL774omjWrJm44447RFVVlbh8+bJo3ry52Lt3r95xU6ZMERMmTDAa76mnnhJjx47VPU9ISBABAQGioqLC6GuISH28MkNENiWEgEajAQDs3r0bMTExuOmmm+Dj44NHH30UxcXFuj4sd911F3r06IHVq1cDAP71r3+hc+fOGDhwoNWxTHnhhRdQU1OD5557Dp6envjll19QXl6OmJgYtG7dWvdYvXo18vLydK9bvnw5+vXrhw4dOqB169b49NNPcfbsWb1z9+rVC15eXhaVl4jUwWSGiGzq+PHjCA0NxZkzZzBy5Ej07NkT69atw+HDh/Hhhx8CAKqqqnTHP/7447pbTcnJyXjsscd0CUpERARKS0tx/vx5gziVlZU4deoUwsPDzZZJq9UCADw9PQFAN+Lqm2++QXZ2tu7xyy+/4KuvvgIArF27FrNnz8bkyZOxfft2ZGdn47HHHjPo5NuqVSur2oeIbI/JDBHZzK5du3D06FGMHTsWhw4dwvXr1/HWW2/h7rvvRkRERINJycSJE3H27Fm89957OHbsGBISEnT7xo4dC09PT7z11lsGr1u+fDmuXLmCCRMmWF3O7t27w9vbG2fPnkVYWJjeIyQkBADw7bffYsCAAXjqqadw2223ISwsTO+qDRE5Dk/ZBSAi51RRUYHCwkJUV1fjv//9L7Zt24akpCT86U9/wqOPPoqjR4/i+vXreP/993Hffffh+++/x/Llyw3O07ZtWzz44IN45plnMHz4cHTq1Em3r3Pnznj99dcxb948NG/eHJMmTYJWq8V//vMfLFy4EHPnzkX//v2tLruPjw/mzZuH2bNno6amBvfccw/Kysqwd+9etG7dGgkJCQgLC8Pq1auRlpaG0NBQ/Otf/0JWVhZCQ0Ob1G5EpALZnXaIyPkkJCTohkB7enqKDh06iOjoaLFy5UpRXV2tO+7tt98WQUFBokWLFiI2NlasXr1aABAXLlzQO1/tkO61a9c2GO8///mPuPfee0WrVq1E8+bNxR133CFWrlxpcJyxDsBCCAFAbNiwQfe8pqZG/OMf/xCRkZFCq9WKDh06iNjYWJGZmSmEEKK8vFwkJiYKPz8/0aZNG/Hkk0+K5557TvTp00evHYx1hCYi+9EIIYTUbIqI3N6///1vzJw5E+fPn2dnWiKyGm8zEZE0V69eRX5+PpKSkjB16lQmMkTUKOwATETSvP766+jbty8CAgKwYMEC2cUhIifF20xERETk1HhlhoiIiJwakxkiIiJyakxmiIiIyKkxmSEiIiKnxmSGiIiInBqTGSIiInJqTGaIiIjIqTGZISIiIqf2/wEUZrSdzlqZdgAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ax=None\n", "colors = ['red','blue','green','yellow']\n", "for i,var in enumerate(new_pumpkins['Variety'].unique()):\n", " df = new_pumpkins[new_pumpkins['Variety']==var]\n", " ax = df.plot.scatter('DayOfYear','Price',ax=ax,c=colors[i],label=var)" ] }, { "cell_type": "markdown", "id": "f6bb5fa9", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Our investigation suggests that variety has more effect on the overall price than the actual selling date. So let us focus for the moment only on one pumpkin variety, and see what effect the date has on the price:" ] }, { "cell_type": "code", "execution_count": 9, "id": "115f94ab", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjMAAAGwCAYAAABcnuQpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA3aklEQVR4nO3dfVyUdb7/8fcIiIAwyo0agYiBPioLrbQ2DbU7cyvz1J6UTt5Ue7YytXSz0trd2tOJardOSXe7bdG6Z2OrTY1Tv25XAsPtxpCj3ayBIlrqCt5AQCLC9/eHh9ERGG6c4Zpr5vV8POYPruua+Xy/Xwbn7TXX9/o6jDFGAAAANtXH6gYAAACcCMIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwtVCrG+BrLS0t2rlzp6Kjo+VwOKxuDgAA6AJjjL7//nslJiaqTx/P514CPszs3LlTycnJVjcDAAD0wI4dO5SUlOTxmIAPM9HR0ZKODEZMTIzFrQEAAF1RW1ur5ORk1+e4JwEfZlq/WoqJiSHMAABgM125RIQLgAEAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK1ZGmays7M1duxYRUdHa9CgQZo+fbo2b97s2t/U1KS7775bZ5xxhqKiopSYmKjZs2dr586dFrYagW5rVZ0KNu9RRXV9QNcEgEBh6dpMhYWFuu222zR27FgdPnxY9957ry699FJ99dVXioqKUkNDg0pKSvSLX/xCGRkZ2r9/v+644w5NmzZN69evt7LpCEAHGg5pYV6pisqqXNsy0xOUkzVGzsiwgKkJAIHGYYwxVjeiVVVVlQYNGqTCwkJlZma2e8xnn32mcePGqbKyUkOHDu30NWtra+V0OlVTU8NCk/Bo9gufqri8Ws3H/EmEOBwanxavFTeNC5iaAGAH3fn89qtrZmpqaiRJsbGxHo9xOBwaMGBAu/sbGxtVW1vr9gA6s7WqTkVlVW6hQpKajVFRWZVPvv6xoiYABCK/CTPGGC1evFgTJkzQqFGj2j3m4MGDuueee3Tdddd1mNKys7PldDpdj+TkZF82GwGicl+Dx/3b9no/WFhREwACkd+Emfnz52vjxo3Ky8trd39TU5NmzpyplpYWPfPMMx2+ztKlS1VTU+N67Nixw1dNRgBJiY30uH9YXFRA1ASAQOQXYWbBggXKz89XQUGBkpKS2uxvamrStddeq4qKCr3//vsevzsLDw9XTEyM2wPozPCE/spMT1CIw+G2PcThUGZ6glLjvR8srKgJAIHI0jBjjNH8+fO1cuVKrVmzRqmpqW2OaQ0yZWVl+uCDDxQXF2dBSxEMcrLGaHxavNu28WnxyskaE1A1ASDQWDqbad68eXr55Zf1xhtvaOTIka7tTqdTEREROnz4sK655hqVlJTozTff1ODBg13HxMbGqm/fvp3WYDYTuquiul7b9tZrWFxUr50dsaImAPiz7nx+WxpmHMedXm+Vm5uruXPnatu2be2erZGkgoICTZo0qdMahBkAAOynO5/flt40r7McNWzYsE6PAQAAwc0vLgAGAADoKcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwNcIMAACwtVCrGwDAf22tqlPlvgYNi4tSanxUwNa0qi59Dcy6wdRXf0GYAdDGgYZDWphXqqKyKte2zPQE5WSNkTMyLGBqWlWXvtLXQKjrT/iaCUAbC/NKVVxe7batuLxaC/I2BFRNq+rSV/oaCHX9CWEGgJutVXUqKqtSszFu25uNUVFZlSqq6wOiplV16St9DYS6/oYwA8BN5b4Gj/u37fX+P45W1LSqLn31bU2r6gZTX/0RYQaAm5TYSI/7h8V5/8JCK2paVZe++ramVXWDqa/+iDADwM3whP7KTE9QiMPhtj3E4VBmeoJPZklYUdOquvSVvgZCXX9DmAHQRk7WGI1Pi3fbNj4tXjlZYwKqplV16St9DYS6/sRhzHFXDQWY2tpaOZ1O1dTUKCYmxurmALZSUV2vbXvre/W+FVbUtKoufQ3MusHUV1/qzuc3YQYAAPid7nx+8zUTAACwNcIMAACwNcIMAACwNcIMAACwNUvDTHZ2tsaOHavo6GgNGjRI06dP1+bNm92OMcbo/vvvV2JioiIiIjRp0iR9+eWXFrU4eG2tqlPB5j29emtsK2paWRcA0DOWrppdWFio2267TWPHjtXhw4d177336tJLL9VXX32lqKgj08oeffRRPf7443rppZc0YsQIPfjgg7rkkku0efNmRUdHW9n8oMDKs8G18iwA2JFfTc2uqqrSoEGDVFhYqMzMTBljlJiYqDvuuEN33323JKmxsVGDBw/WI488optvvrnT12Rq9omZ/cKnKi6vdlvELMTh0Pi0eK24aVzA1LSyLgCgLdtOza6pqZEkxcbGSpIqKiq0e/duXXrppa5jwsPDNXHiRK1bt67d12hsbFRtba3bAz3DyrPBt/IsANiR34QZY4wWL16sCRMmaNSoUZKk3bt3S5IGDx7sduzgwYNd+46XnZ0tp9PpeiQnJ/u24QGMlWd9XxcAcOL8JszMnz9fGzduVF5eXpt9juMW0DLGtNnWaunSpaqpqXE9duzY4ZP2BgNWnvV9XQDAifOLMLNgwQLl5+eroKBASUlJru1DhgyRpDZnYfbs2dPmbE2r8PBwxcTEuD3QM6w8G3wrzwKAHVkaZowxmj9/vlauXKk1a9YoNTXVbX9qaqqGDBmi999/37Xt0KFDKiws1Pnnn9/bzQ1KrDwbXCvPAoAdWTqbad68eXr55Zf1xhtvaOTIka7tTqdTERERkqRHHnlE2dnZys3NVXp6uh566CF9+OGHXZ6azWwm72DlWQBAb7LNqtkdXfeSm5uruXPnSjpy9uaBBx7Q7373O+3fv1/nnnuunn76addFwp0hzAAAYD+2CTO9gTADAID92PY+MwAAAN1FmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALZGmAEAALYWanUDYA9bq+pUua9Bw+KilBofFbA1u1LXF+2yoqY32hUoNa2qG0x9BXyJMAOPDjQc0sK8UhWVVbm2ZaYnKCdrjJyRYQFTsyt1fdEuK2p6o12BUtOqusHUV6A38DUTPFqYV6ri8mq3bcXl1VqQtyGganalri/aZUVNb7QrUGpaVTeY+gr0BsIMOrS1qk5FZVVqNsZte7MxKiqrUkV1fUDU7Erdom/2eL1dVtT0RruC6fdOXwF7IMygQ5X7Gjzu37bX+//4WVGzK3U37DjgcX9P2mVFza7g9+7busHUV6C3EGbQoZTYSI/7h8V5/8JBK2p2pe6Y5AEe9/ekXVbU7Ap+776tG0x9BXoLYQYdGp7QX5npCQpxONy2hzgcykxP8MksCCtqdqVu5ohBXm+XFTW90a5g+r3TV8AeCDPwKCdrjManxbttG58Wr5ysMQFVsyt1fdEuK2p6o12BUtOqusHUV6A3OIw57mqwAFNbWyun06mamhrFxMRY3Rzbqqiu17a99b16Xworanalri/aZUVNb7QrUGpaVTeY+gp0V3c+vwkzAADA73Tn85uvmQAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK0RZgAAgK2FWt0AO9taVafKfQ0aFhel1PiogK1pVV36CgDoCsJMDxxoOKSFeaUqKqtybctMT1BO1hg5I8MCpqZVdemr7/sKAIGEr5l6YGFeqYrLq922FZdXa0HehoCqaVVd+ur7vgJAICHMdNPWqjoVlVWp2Ri37c3GqKisShXV9QFR06q69NX3fQWAQEOY6abKfQ0e92/b6/0PICtqWlWXvvq2JgAEIkvDTFFRka688kolJibK4XBo9erVbvvr6uo0f/58JSUlKSIiQqeeeqqeffZZaxr7f1JiIz3uHxbn/Ys3rahpVV366tuaABCILA0z9fX1ysjI0FNPPdXu/kWLFumdd97Rf//3f+vrr7/WokWLtGDBAr3xxhu93NKjhif0V2Z6gkIcDrftIQ6HMtMTfDITxYqaVtWlr77vKwAEGkvDzNSpU/Xggw/q6quvbnf/3//+d82ZM0eTJk3SsGHD9LOf/UwZGRlav359L7fUXU7WGI1Pi3fbNj4tXjlZYwKqplV16avv+woAgcRhzHFXH1rE4XBo1apVmj59umvbLbfcos8//1yrV69WYmKiPvzwQ02bNk1vv/22JkyY0O7rNDY2qrGx0fVzbW2tkpOTVVNTo5iYGK+2uaK6Xtv21vfqvUGsqGlVXfoKAMGrtrZWTqezS5/ffh1mDh06pH//93/XihUrFBoaqj59+ugPf/iDZs2a1eHr3H///XrggQfabPdFmAEAAL7RnTDj17OZli9fro8//lj5+fn6/PPP9dhjj2nevHn64IMPOnzO0qVLVVNT43rs2LGjF1sMAAB6m9/eAfiHH37QsmXLtGrVKl1++eWSpDPPPFOlpaX67W9/q4svvrjd54WHhys8PLw3mwoAACzkt2dmmpqa1NTUpD593JsYEhKilpYWi1oFAAD8jaVnZurq6lReXu76uaKiQqWlpYqNjdXQoUM1ceJELVmyRBEREUpJSVFhYaFWrFihxx9/3MJWAwAAf2LpBcAffvihJk+e3Gb7nDlz9NJLL2n37t1aunSp3nvvPe3bt08pKSn62c9+pkWLFslx3L05OtKdC4gAAIB/sOVsJl8hzAAAYD8BM5sJAACgM4QZAABga4QZAABga4QZAABga4QZAABga4QZAABga4QZAABga4QZAABga4QZAABga4QZAABga4QZAABga5aumm13W6vqVLmvQcPiopQaH2V1cwKOVeNrRV3eSwDQc4SZHjjQcEgL80pVVFbl2paZnqCcrDFyRoZZ2LLAYNX4WlGX9xIAnDi+ZuqBhXmlKi6vdttWXF6tBXkbLGpRYLFqfK2oy3sJAE4cYaabtlbVqaisSs3GuG1vNkZFZVWqqK63qGWBwarxtaIu7yUA8A7CTDdV7mvwuH/bXj6AToRV42tFXd5LAOAdhJluSomN9Lh/WBwXb54Iq8bXirq8lwDAOwgz3TQ8ob8y0xMU4nC4bQ9xOJSZnsBMlBNk1fhaUZf3EgB4B2GmB3Kyxmh8WrzbtvFp8crJGmNRiwKLVeNrRV3eSwBw4hzGHHf1YYCpra2V0+lUTU2NYmJivPraFdX12ra3nnuD+IhV42tFXd5LR3F/ocCsy72U0F3d+fwmzADwC9xfiL4Cx+rO5zdfMwHwC9xfiL4CPUWYAWA57i9EX4ETQZgBYDnuL+TbmlbV5V5K6C2EGQCW4/5Cvq1pVV3upYTeQpgBYDnuL0RfgRNBmAHgF7i/EH0FeuqEpmaXl5dry5YtyszMVEREhIwxchyXwK3G1GzAXri/UGDW5V5K6C6f32dm7969mjFjhtasWSOHw6GysjINHz5cN910kwYMGKDHHnusx433NsIMAAD24/P7zCxatEihoaHavn27IiOPXuA1Y8YMvfPOOz15SQAAgB4J7cmT3nvvPb377rtKSkpy256enq7KykqvNAwAAKArenRmpr6+3u2MTKvq6mqFh4efcKMAAAC6qkdhJjMzUytWrHD97HA41NLSot/85jeaPHmy1xoHAADQmR59zfSb3/xGkyZN0vr163Xo0CHddddd+vLLL7Vv3z4VFxd7u40AAAAd6tGZmdNOO00bN27UuHHjdMkll6i+vl5XX321NmzYoFNOOcXbbQQAAOjQCd1nxg6Ymg0AgP34fGp2bm6uXnvttTbbX3vtNf3xj3/syUsCAAD0SI/CzMMPP6z4+Pg22wcNGqSHHnrohBsFAADQVT0KM5WVlUpNTW2zPSUlRdu3bz/hRgEAAHRVj8LMoEGDtHHjxjbb//d//1dxcXEn3Ch0bGtVnQo271FFdX3A133l0+2645UNem39jl6rKUmFm/foyb99o7VlVe3u98VYdPaawfR7B4Du6tHU7JkzZ2rhwoWKjo5WZmamJKmwsFC33367Zs6c6dUG4ogDDYe0MK9URcd8wGamJygna4yckWEBVXfTtwf0L8+s0+GWI9emr96wU0tXblL+beN12slOn9SUpMq99Zr+dLH2NzS5tg2MDFP+bROUHBfpk7Ho7DWD6fcOAD3VozMzDz74oM4991xddNFFioiIUEREhC699FJdeOGFXDPjIwvzSlVcXu22rbi8WgvyNgRc3WODTKvDLUbTnvbtPYyODzKStL+hSdOe/kiSb8ais9cMpt87APRUj8JM37599corr+gf//iH/vznP2vlypXasmWLXnzxRfXt29fbbQx6W6vqVFRWpebjZtE3G6OisiqffQVgRd1XPt3eJsi0OtxifPaVU+HmPW2CTKv9DU169bPtXh+Lzsa36Js9QfN7B4AT0aMw02rEiBH613/9V11xxRVKSUnxVptwnMp9DR73b9vrmw8XK+r+vWKvx/3FW6o97u+p0m8PeNy/bqvndvVkLDob3w07Dni9ZldY9X4DgJ7q8jUzixcv1n/8x38oKipKixcv9njs448/fsINw1EpsW0X9TzWsLiogKn7o9Q4rd6ws8P9409pe0sAbxidNMDj/vOHe25XT8ais/Edk+y5TYH0eweAE9HlMLNhwwY1NR05DV9SUiKHw9HucR1tR88NT+ivzPQEFZdXu536D3E4ND4tXqnxvvlwsaLujHFDde/qL9r9qim0j0P/ek6y12tK0sSRgzQwMqzdr5oGRobp2rFD9ebG3V4di87GN3PEoKD5vQPAiWA5A5uoaWjSgrwNvT67xIq6X31Xo2lPF7sFmtA+Dp/PZtqxt0HTnv6ow9lMvhiLzl4zmH7vrbZW1alyX4OGxUX1anCyqi6A9nXn87vbYebw4cPq16+fSktLNWrUqBNqaG8IlDDTqqK6Xtv21vf6P7hW1H1t/Q4Vb6nW+FPifXZGpj1ry6pUsn2/zho6UBekJ7TZ74ux6Ow1e3v8rZiazTR0AMfyaZiRpFNOOUUrV65URkZGjxvZWwItzAC9YfYLn3b4NdOKm8YFTE0r6wLwzOcLTd53331aunSp9u3b16MGAvBfVkzNDqbbDwDwvh6FmeXLl2vt2rVKTEzUyJEjddZZZ7k9uqqoqEhXXnmlEhMT5XA4tHr16jbHfP3115o2bZqcTqeio6N13nnnsf4T4ENWTM0OptsPAPC+Hi1nMH36dDkcDp3otcP19fXKyMjQDTfcoGuuuabN/i1btmjChAm66aab9MADD8jpdOrrr79Wv379TqgugI5ZMTU7mG4/AMD7uhVmGhoatGTJEq1evVpNTU266KKLlJOTo/j4nt37Y+rUqZo6dWqH+++99179+Mc/1qOPPuraNnz48B7VArqqs1ktVsx6eeXT7fp7xd5euRjaiqnZwXT7AQDe160w86tf/UovvfSS/u3f/k0RERF6+eWXdeutt+q1117zesNaWlr01ltv6a677tKUKVO0YcMGpaamaunSpZo+fXqHz2tsbFRjY6Pr59raWq+3DYHJHxd9tGrRzVsnDnfrp3TkOpLbJp0SUDUlKSdrTJtp6OPT4pWTNcandQF4T7dmM51yyin6z//8T9fK2J9++qnGjx+vgwcPKiQk5MQa4nBo1apVrqCye/dunXTSSYqMjNSDDz6oyZMn65133tGyZctUUFCgiRMntvs6999/vx544IE225nNhM50NqvFilkvacv+X4c3ECx/6Mc+qWlVXav62sqq2x4AaJ/PZjPt2LFDF1xwgevncePGKTQ0VDt3dnyb955qaWmRJF111VVatGiRRo8erXvuuUdXXHGFnnvuuQ6ft3TpUtXU1LgeO3b4ZmFCBBZ/XPTRqkU3rahrVV+PlRofpckjBxFkABvqVphpbm5usyp2aGioDh8+7NVGSVJ8fLxCQ0N12mmnuW0/9dRTPc5mCg8PV0xMjNsD6Iw/Lvpo1aKbVtS1qq8AAkO3rpkxxmju3LkKDw93bTt48KBuueUWRUUd/d/MypUrT7hhffv21dixY7V582a37d988w0rdMPr/HHRR6sW3bSirlV9BRAYunVmZs6cORo0aJCcTqfrcf311ysxMdFtW1fV1dWptLRUpaWlkqSKigqVlpa6zrwsWbJEr7zyip5//nmVl5frqaee0v/8z/9o3rx53Wk20KnWWS0hxy2UGuJwKDM9wbXoY0f7fbXoZmif9hdu9eWim1bUtaqvAAKDpQtNfvjhh5o8eXKb7XPmzNFLL70kSXrxxReVnZ2tb7/9ViNHjtQDDzygq666qss1WM4AXeWPiz5ateimFXWt6isA/+TztZnshDCD7vK3RR8l6xbdtKKuVX0F4F8IM8cgzAAAYD8+X2gSAADAXxBmAACArRFmAACArRFmAACArXXrpnkAencF61aFm/eo9NsDOmvoQF2QntArNa2qa1Vff53/pdZtrdaEtATdd8VpnT/BpjUl6am/lal4S7UuSE/QvMlpvVLTir8bq+ry99r7mM0EdNHxK1hLvr8PSuXeek1/ulj7G5pc2wZGhin/tglKjvN812K71bWqr+99sUs/+++SNttfmH2OLjptcMDUlKR15VW67g+fttn+yr+fp3NPifNJTSv+bqyqy9+rd2symwnwgeP/YZSOLII47elin9U8/h8LSdrf0KRpT3/ks5pW1bWqr+2FCkm6acX6gKopqd0gI0kznv/YZzWt+Luxqi5/r77va0cIM0AXWLGqc+HmPW3+sWi1v6FJa4+5E7Hd61rV11/nf+lx/4NvfhUQNaUjXy158kxBuddrBtPK7/y9+rZmZwgzQBdYsapz6bcHPO4v2b7f6zWtqmtVX9dt9fx7+6jc+/8oW1FT6vw96osPoGBa+Z2/V9/W7AxhBuiCH6V6vp7AF6s6j04a4HH/WUMHer2mVXWt6uv5wz3/3iakef+CRitqSp2/R31x8aYVfzdW1eXv1bc1O0OYAbrAilWdJ44cpIEdLGA5MDLMZzMHrKhrVV9/Oe10j/t9McPIipqSNP+idI/7fTGrKZhWfufv1bc1O0OYAboo/7bxbf6BbJ0d4buaE9r8o9E6Y8CXrKhrVV9fmH1Ot7bbtaZ0ZNZSd7Z7gxV/N1bV5e/V933tCFOzgW6yYlXntWVVKtm+v9fv5WBFXav6+uCbX+mj8qpeveeLFTWlIxf7ri2r6tX7zATTyu/8vXoHq2YfgzADAID9cJ8ZAAAQNAgzAADA1ggzAADA1ggzAADA1lg1GwAssrWqTpX7GjQsLkqp8VFWNwewLcIMAPSyAw2HtDCvVEXHLCGQmZ6gnKwxcnZwMzIAHeNrJgDoZQvzSlVc7r4+UHF5tRbkbbCoRYC9EWYAoBdtrapTUVmVmo+7xVezMSoqq1JFdb1FLQPsizADAL2ocl+Dx/3b9hJmgO4izABAL0qJjfS4f1gcFwID3UWYAYBeNDyhvzLTExTicF8EMcThUGZ6ArOagB4gzABAL8vJGqPxafFu28anxSsna4xFLQLsjanZANDLnJFhWnHTOFVU12vb3nruMwOcIMIMAFgkNZ4QA3gDXzMBAABbI8wAAABbI8wAAABbI8wAAABb4wJg4DiFm/eo9NsDOmvoQF2QnhCwNQEgUBBmgP9Tubde058u1v6GJte2gZFhyr9tgpLjPN+11U41ASDQ8DUT8H+ODxWStL+hSdOe/iigagJAoCHMADryNc/xoaLV/oYmrS2rCoiaABCICDOApNJvD3jcX7J9f0DUBIBARJgBJI1OGuBx/1lDBwZETQAIRIQZQNLEkYM0MDKs3X0DI8N8MsPIipoAEIgIM8D/yb9tQptw0TqzKJBqAkCgcRhjjNWN8KXa2lo5nU7V1NQoJibG6ubABtaWValk+/5eveeLFTUBwJ915/ObMAMAAPxOdz6/+ZoJAADYGmEGAADYGmEGAADYGmEGAADYGgtNAt0UTKtq09fAq2lVXfoamHW3VtWpcl+DhsVFKTU+qldqtocwA3RRMK2qTV/pq51rWlU3mPp6oOGQFuaVquiYNeQy0xOUkzVGzg5uBupLfM0EdFEwrapNX+mrnWtaVTeY+rowr1TF5dVu24rLq7Ugb4PPanpiaZgpKirSlVdeqcTERDkcDq1evbrDY2+++WY5HA498cQTvdY+oFUwrapNX31bl776tqZVdYOpr1ur6lRUVqXm425T12yMisqqVFFd7/WanbE0zNTX1ysjI0NPPfWUx+NWr16tTz75RImJib3UMsBdMK2qTV99W5e++ramVXWDqa+V+xo87t+2t/fDjKXXzEydOlVTp071eMx3332n+fPn691339Xll1/e6Ws2NjaqsbHR9XNtbe0JtxMIplW16atv69JX39a0qm4w9TUl1vN1OMPiev9CYL++ZqalpUWzZs3SkiVLdPrpp3fpOdnZ2XI6na5HcnKyj1uJYBBMq2rTV9/Wpa++rWlV3WDq6/CE/spMT1CIw+G2PcThUGZ6giWzmvw6zDzyyCMKDQ3VwoULu/ycpUuXqqamxvXYsWOHD1uIYBJMq2rTV/pq55pW1Q2mvuZkjdH4tHi3bePT4pWTNcZnNT3xm4UmHQ6HVq1apenTp0uSPv/8c11++eUqKSlxXSszbNgw3XHHHbrjjju6/LosNAlvC6ZVtelr4NW0qi59Dcy6FdX12ra33if3mbHlqtnHh5knnnhCixcvVp8+R08eNTc3q0+fPkpOTta2bdu69LqEGQAA7Kc7n99+e9O8WbNm6eKLL3bbNmXKFM2aNUs33HCDRa0CAAD+xtIwU1dXp/LyctfPFRUVKi0tVWxsrIYOHaq4uDi348PCwjRkyBCNHDmyt5sKAAD8lKVhZv369Zo8ebLr58WLF0uS5syZo5deesmiVgEAADuxNMxMmjRJ3blkp6vXyQC+5IuF1fxlsbbjWdEuq8YimPoKBBq/vWYG8De+WFjN3xZrs7JdVo1FMPUVCFR+fZ8ZwJ/4YmE1f1usrZUV7bJqLIKpr0CgIswAXeCLhdX8cbE2q9pl1VgEU1+BQEaYAbrAFwur+eNibZI17bJqLIKpr0AgI8wAXeCLhdX8cbE2yZp2WTUWwdRXIJARZoAu8MXCav64WJtV7bJqLIKpr0AgI8wAXeSLhdX8bbG2Vla0y6qxCKa+ttpaVaeCzXt69focK2paVdeqvgYzv1mbyVdYmwne5ouF1Xy5WNuJsKJdVo1FMPQ1mKahB1NfA5UtF5r0FcIMABwx+4VPVVxe7TaTKsTh0Pi0eK24aVzA1LSqrlV9DVTd+fzmayYACALBNA09mPqKIwgzABAEgmkaejD1FUcQZgAgCATTNPRg6iuOIMwAQBAIpmnowdRXHEGYAYAgEUzT0IOpr2A2EwAEnWCYhm5lXX+91YLdMDX7GIQZAADsh6nZAAAgaBBmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArRFmAACArVkaZoqKinTllVcqMTFRDodDq1evdu1ramrS3XffrTPOOENRUVFKTEzU7NmztXPnTusaDAAA/I6lYaa+vl4ZGRl66qmn2uxraGhQSUmJfvGLX6ikpEQrV67UN998o2nTplnQUgAA4K8cxhhjdSMkyeFwaNWqVZo+fXqHx3z22WcaN26cKisrNXTo0C69bm1trZxOp2pqahQTE+Ol1gIAAF/qzud3aC+1yStqamrkcDg0YMCADo9pbGxUY2Oj6+fa2tpeaBkAALCKbS4APnjwoO655x5dd911HhNadna2nE6n65GcnNyLrQQAAL3NFmGmqalJM2fOVEtLi5555hmPxy5dulQ1NTWux44dO3qplQC8YWtVnQo271FFdb3VTQFgE37/NVNTU5OuvfZaVVRUaM2aNZ1+bxYeHq7w8PBeah0AbznQcEgL80pVVFbl2paZnqCcrDFyRoZZ2DIA/s6vz8y0BpmysjJ98MEHiouLs7pJAHxkYV6pisur3bYVl1drQd4Gi1oEwC4sPTNTV1en8vJy188VFRUqLS1VbGysEhMT9ZOf/EQlJSV688031dzcrN27d0uSYmNj1bdvX6uaDcDLtlbVuZ2RadVsjIrKqlRRXa/U+CgLWgbADiwNM+vXr9fkyZNdPy9evFiSNGfOHN1///3Kz8+XJI0ePdrteQUFBZo0aVJvNROAj1Xua/C4f9tewgyAjlkaZiZNmiRPt7nxk1vgAPCxlNhIj/uHxRFkAHTMr6+ZARAchif0V2Z6gkIcDrftIQ6HMtMTOCsDwCPCDAC/kJM1RuPT4t22jU+LV07WGItaBMAu/H5qNoDg4IwM04qbxqmiul7b9tZrWFwUZ2QAdAlhBoBfSY0nxADoHr5mAgAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAtkaYAQAAthZqdQMAwB8Ubt6j0m8P6KyhA3VBekLA1pSkrVV1qtzXoGFxUUqNj+q1uoCvEGYABLXKvfWa/nSx9jc0ubYNjAxT/m0TlBwXGTA1JelAwyEtzCtVUVmVa1tmeoJyssbIGRnms7qAr/E1E4CgdnyokKT9DU2a9vRHAVVTkhbmlaq4vNptW3F5tRbkbfBpXcDXCDMAglbh5j1tQkWr/Q1NWnvMGQw715SOfLVUVFalZmPctjcbo6KyKlVU1/ukLtAbCDMAglbptwc87i/Zvj8gakpS5b4Gj/u37SXMwL4IMwCC1uikAR73nzV0YEDUlKSUWM/X4gyL40Jg2BdhBkDQmjhykAZ2cOHrwMgwn8wwsqKmJA1P6K/M9ASFOBxu20McDmWmJzCrCbZGmAEQ1PJvm9AmXLTOLAqkmpKUkzVG49Pi3baNT4tXTtYYn9YFfM1hzHFXgwWY2tpaOZ1O1dTUKCYmxurmAPBTa8uqVLJ9f6/e88WKmpJUUV2vbXvruc8M/Fp3Pr8JMwAAwO905/Obr5kAAICtEWYAAICtEWYAAICtEWYAAICtEWYAAICtEWYAAICtEWYAAICtEWYAAICtEWYAAICtEWYAAICthVrdAF9rXa2htrbW4pYAAICuav3c7sqqSwEfZr7//ntJUnJyssUtAQAA3fX999/L6XR6PCbgF5psaWnRzp07FR0dLYfDYWlbamtrlZycrB07drDopRcwnt7HmHoX4+ldjKf3+fOYGmP0/fffKzExUX36eL4qJuDPzPTp00dJSUlWN8NNTEyM371p7Izx9D7G1LsYT+9iPL3PX8e0szMyrbgAGAAA2BphBgAA2BphpheFh4frV7/6lcLDw61uSkBgPL2PMfUuxtO7GE/vC5QxDfgLgAEAQGDjzAwAALA1wgwAALA1wgwAALA1wgwAALA1wswJys7O1tixYxUdHa1BgwZp+vTp2rx5c4fH33zzzXI4HHriiSfctjc2NmrBggWKj49XVFSUpk2bpm+//dbHrfc/XR3Pr7/+WtOmTZPT6VR0dLTOO+88bd++3bWf8TyiK+NZV1en+fPnKykpSRERETr11FP17LPPuh3DeB717LPP6swzz3TdZOxHP/qR3n77bdd+Y4zuv/9+JSYmKiIiQpMmTdKXX37p9hqM51GexrOpqUl33323zjjjDEVFRSkxMVGzZ8/Wzp073V6D8Tyqs/fnsQLq88jghEyZMsXk5uaaL774wpSWlprLL7/cDB061NTV1bU5dtWqVSYjI8MkJiaa//qv/3Lbd8stt5iTTz7ZvP/++6akpMRMnjzZZGRkmMOHD/dST/xDV8azvLzcxMbGmiVLlpiSkhKzZcsW8+abb5p//vOfrmMYzyO6Mp4//elPzSmnnGIKCgpMRUWF+d3vfmdCQkLM6tWrXccwnkfl5+ebt956y2zevNls3rzZLFu2zISFhZkvvvjCGGPMww8/bKKjo83rr79uNm3aZGbMmGFOOukkU1tb63oNxvMoT+N54MABc/HFF5tXXnnF/OMf/zB///vfzbnnnmvOPvtst9dgPI/q7P3ZKtA+jwgzXrZnzx4jyRQWFrpt//bbb83JJ59svvjiC5OSkuL25jlw4IAJCwszf/nLX1zbvvvuO9OnTx/zzjvv9FbT/VJ74zljxgxz/fXXd/gcxrNj7Y3n6aefbn7961+7HXfWWWeZ++67zxjDeHbFwIEDzR/+8AfT0tJihgwZYh5++GHXvoMHDxqn02mee+45Ywzj2RWt49meTz/91EgylZWVxhjGsyuOH89A/DziayYvq6mpkSTFxsa6trW0tGjWrFlasmSJTj/99DbP+fzzz9XU1KRLL73UtS0xMVGjRo3SunXrfN9oP3b8eLa0tOitt97SiBEjNGXKFA0aNEjnnnuuVq9e7XoO49mx9t6fEyZMUH5+vr777jsZY1RQUKBvvvlGU6ZMkcR4etLc3Ky//OUvqq+v149+9CNVVFRo9+7dbmMVHh6uiRMnusaK8ezY8ePZnpqaGjkcDg0YMEAS4+lJe+MZqJ9HhBkvMsZo8eLFmjBhgkaNGuXa/sgjjyg0NFQLFy5s93m7d+9W3759NXDgQLftgwcP1u7du33aZn/W3nju2bNHdXV1evjhh3XZZZfpvffe07/8y7/o6quvVmFhoSTGsyMdvT+XL1+u0047TUlJSerbt68uu+wyPfPMM5owYYIkxrM9mzZtUv/+/RUeHq5bbrlFq1at0mmnneYaj8GDB7sdf+xYMZ5tdTSexzt48KDuueceXXfdda5FERnPtjyNZ6B+HgX8qtm9af78+dq4caM++ugj17bPP/9cTz75pEpKSuRwOLr1esaYbj8nkLQ3ni0tLZKkq666SosWLZIkjR49WuvWrdNzzz2niRMndvh6jGfb8ZSOhJmPP/5Y+fn5SklJUVFRkebNm6eTTjpJF198cYevF8zjOXLkSJWWlurAgQN6/fXXNWfOHFeYltRmXLoyVoxn2/E8NtA0NTVp5syZamlp0TPPPNPpazKebcfzhx9+CNjPI87MeMmCBQuUn5+vgoICJSUlubavXbtWe/bs0dChQxUaGqrQ0FBVVlbq5z//uYYNGyZJGjJkiA4dOqT9+/e7veaePXva/A8vWHQ0nvHx8QoNDW3zv7ZTTz3VNZuJ8Wyro/H84YcftGzZMj3++OO68sordeaZZ2r+/PmaMWOGfvvb30piPNvTt29fpaWl6ZxzzlF2drYyMjL05JNPasiQIZLU5n+wx44V49lWR+PZqqmpSddee60qKir0/vvvu87KSIxnezoaz0D+PCLMnCBjjObPn6+VK1dqzZo1Sk1Ndds/a9Ysbdy4UaWlpa5HYmKilixZonfffVeSdPbZZyssLEzvv/++63m7du3SF198ofPPP79X+2O1zsazb9++Gjt2bJvpxd98841SUlIkMZ7H6mw8m5qa1NTUpD593P8pCAkJcZ0FYzw7Z4xRY2OjUlNTNWTIELexOnTokAoLC11jxXh2rnU8paNBpqysTB988IHi4uLcjmU8O9c6ngH9edTrlxwHmFtvvdU4nU7z4Ycfml27drkeDQ0NHT7n+KvHjTkyFS4pKcl88MEHpqSkxFx44YV+PxXOF7oynitXrjRhYWHm97//vSkrKzM5OTkmJCTErF271nUM43lEV8Zz4sSJ5vTTTzcFBQVm69atJjc31/Tr188888wzrmMYz6OWLl1qioqKTEVFhdm4caNZtmyZ6dOnj3nvvfeMMUemZjudTrNy5UqzadMmk5WV1e7UbMbzCE/j2dTUZKZNm2aSkpJMaWmp23u4sbHR9RqM51GdvT+PFyifR4SZEySp3Udubm6Hz2nvzfPDDz+Y+fPnm9jYWBMREWGuuOIKs337dt823g91dTxfeOEFk5aWZvr162cyMjLc7oliDOPZqivjuWvXLjN37lyTmJho+vXrZ0aOHGkee+wx09LS4jqG8TzqxhtvNCkpKaZv374mISHBXHTRRW4fFC0tLeZXv/qVGTJkiAkPDzeZmZlm06ZNbq/BeB7laTwrKio6fA8XFBS4XoPxPKqz9+fxAuXzyGGMMb13HggAAMC7uGYGAADYGmEGAADYGmEGAADYGmEGAADYGmEGAADYGmEGAADYGmEGAADYGmEGAADYGmEGQMAoLi7WGWecobCwME2fPt3q5gDoJYQZAN02d+5cORwOORwOhYWFafDgwbrkkkv04osvuhao9KY333xTkyZNUnR0tCIjIzV27Fi99NJLbY5bvHixRo8erYqKCi1atEhhYWH66KOP3I6pr6/X8OHDtWjRIq+3E4A1CDMAeuSyyy7Trl27tG3bNr399tuaPHmybr/9dl1xxRU6fPiw1+rk5OToqquu0vnnn69PPvlEGzdu1MyZM3XLLbfozjvvdDt2y5YtuvDCC5WUlKSJEydqwYIFmjt3rurr613H3HXXXQoPD1d2drbX2tjq0KFDXn9NAF1g9eJQAOxnzpw55qqrrmqz/W9/+5uRZJ5//nljjDGPPfaYGTVqlImMjDRJSUnm1ltvNd9//70xxpi6ujoTHR1tXnvtNbfXyM/PN5GRkaa2ttZs377dhIWFmcWLF7eptXz5ciPJfPzxx+0uSJibm2t++OEHc+qpp5rbbrvNGGPMmjVrTFhYmPnss89MS0uLeeSRR0xqaqrp16+fOfPMM93acvjwYXPjjTeaYcOGmX79+pkRI0aYJ554ot1xeOihh8xJJ51kUlJSTmRYAfQQZ2YAeM2FF16ojIwMrVy5UpLUp08fLV++XF988YX++Mc/as2aNbrrrrskSVFRUZo5c6Zyc3PdXiM3N1c/+clPFB0drb/+9a9qampqcwZGkm6++Wb1799feXl5Sk5O1q5duxQTE6MnnnhCu3bt0owZM9SvXz+tWLFCv//977V69WrdeOONWrZsmc455xzdd999ys3N1bPPPqsvv/xSixYt0vXXX6/CwkJJUktLi5KSkvTqq6/qq6++0i9/+UstW7ZMr776qls7/va3v+nrr7/W+++/rzfffNMXwwqgM1anKQD209GZGWOMmTFjhjn11FPb3ffqq6+auLg418+ffPKJCQkJMd99950xxpiqqioTFhZmPvzwQ2OMMbfccotxOp0dtuPMM880U6dOdf3sdDpNbm5um+N++ctfmj59+pizzz7bNDU1mbq6OtOvXz+zbt06t+Nuuukmk5WV1WG9efPmmWuuucb185w5c8zgwYNNY2Njh88B4HucmQHgVcYYORwOSVJBQYEuueQSnXzyyYqOjtbs2bO1d+9e1zUs48aN0+mnn64VK1ZIkv70pz9p6NChyszM7HYtT+677z61tLTonnvuUWhoqL766isdPHhQl1xyifr37+96rFixQlu2bHE977nnntM555yjhIQE9e/fX88//7y2b9/u9tpnnHGG+vbt26X2AvANwgwAr/r666+VmpqqyspK/fjHP9aoUaP0+uuv6/PPP9fTTz8tSWpqanId/9Of/tT1VVNubq5uuOEGV0AZMWKEampqtHPnzjZ1Dh06pK1btyo9Pb3TNoWFhUmSQkNDJck14+qtt95SaWmp6/HVV1/pr3/9qyTp1Vdf1aJFi3TjjTfqvffeU2lpqW644YY2F/lGRUV1a3wAeB9hBoDXrFmzRps2bdI111yj9evX6/Dhw3rsscd03nnnacSIEe2Gkuuvv17bt2/X8uXL9eWXX2rOnDmufddcc41CQ0P12GOPtXnec889p/r6emVlZXW7naeddprCw8O1fft2paWluT2Sk5MlSWvXrtX555+vefPmacyYMUpLS3M7awPAf4Ra3QAA9tTY2Kjdu3erublZ//znP/XOO+8oOztbV1xxhWbPnq1Nmzbp8OHDysnJ0ZVXXqni4mI999xzbV5n4MCBuvrqq7VkyRJdeumlSkpKcu0bOnSoHn30Ud15553q16+fZs2apbCwML3xxhtatmyZfv7zn+vcc8/tdtujo6N15513atGiRWppadGECRNUW1urdevWqX///pozZ47S0tK0YsUKvfvuu0pNTdWf/vQnffbZZ0pNTT2hcQPgA1ZftAPAfubMmeOaAh0aGmoSEhLMxRdfbF588UXT3NzsOu7xxx83J510komIiDBTpkwxK1asMJLM/v373V6vdUr3q6++2m69N954w1xwwQUmKirK9OvXz5x99tnmxRdfbHNcRxcAG2OMJLNq1SrXzy0tLebJJ580I0eONGFhYSYhIcFMmTLFFBYWGmOMOXjwoJk7d65xOp1mwIAB5tZbbzX33HOPycjIcBuHji6EBtB7HMYYY2maAhD0/vznP+v222/Xzp07uZgWQLfxNRMAyzQ0NKiiokLZ2dm6+eabCTIAeoQLgAFY5tFHH9Xo0aM1ePBgLV261OrmALApvmYCAAC2xpkZAABga4QZAABga4QZAABga4QZAABga4QZAABga4QZAABga4QZAABga4QZAABga/8fzaEy1CONIYMAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "pie_pumpkins = new_pumpkins[new_pumpkins[\"Variety\"] == \"PIE TYPE\"]\n", "pie_pumpkins.plot.scatter(\"DayOfYear\", \"Price\")" ] }, { "cell_type": "markdown", "id": "2527b9e5", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "If we now calculate the correlation between `Price` and `DayOfYear` using `corr` function, we will get something like `-0.27` - which means that training a predictive model makes sense.\n", "\n", ":::{note}\n", "Before training a linear regression model, it is important to make sure that our data is clean. Linear regression does not work well with missing values, thus it makes sense to get rid of all empty cells:\n", ":::" ] }, { "cell_type": "code", "execution_count": 10, "id": "242febdf", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Int64Index: 144 entries, 70 to 1630\n", "Data columns (total 8 columns):\n", " # Column Non-Null Count Dtype \n", "--- ------ -------------- ----- \n", " 0 Month 144 non-null int64 \n", " 1 DayOfYear 144 non-null int64 \n", " 2 Variety 144 non-null object \n", " 3 City 144 non-null object \n", " 4 Package 144 non-null object \n", " 5 Low Price 144 non-null float64\n", " 6 High Price 144 non-null float64\n", " 7 Price 144 non-null float64\n", "dtypes: float64(3), int64(2), object(3)\n", "memory usage: 10.1+ KB\n" ] } ], "source": [ "pie_pumpkins = pie_pumpkins.dropna()\n", "pie_pumpkins.info()" ] }, { "cell_type": "markdown", "id": "92709a1f", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Another approach would be to fill those empty values with mean values from the corresponding column.\n", "\n", "## Simple Linear Regression\n", "\n", "To train our Linear Regression model, we will use the **Scikit-learn** library." ] }, { "cell_type": "code", "execution_count": 11, "id": "c3d17ae6", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "from sklearn.linear_model import LinearRegression\n", "from sklearn.metrics import mean_squared_error\n", "from sklearn.model_selection import train_test_split" ] }, { "cell_type": "markdown", "id": "02c3a570", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "We start by separating input values (features) and the expected output (label) into separate NumPy arrays:" ] }, { "cell_type": "code", "execution_count": 12, "id": "8a96563a", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "X = pie_pumpkins[\"DayOfYear\"].to_numpy().reshape(-1, 1)\n", "y = pie_pumpkins[\"Price\"]" ] }, { "cell_type": "markdown", "id": "6bfab9bb", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ ":::{note}\n", "Note that we had to perform `reshape` on the input data in order for the Linear Regression package to understand it correctly. Linear Regression expects a 2D-array as an input, where each row of the array corresponds to a vector of input features. In our case, since we have only one input - we need an array with shape N×1, where N is the dataset size.\n", ":::\n", "\n", "Then, we need to split the data into train and test datasets, so that we can validate our model after training:" ] }, { "cell_type": "code", "execution_count": 13, "id": "63a5c90e", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=0)" ] }, { "cell_type": "markdown", "id": "975937b9", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Finally, training the actual Linear Regression model takes only two lines of code. We define the `LinearRegression` object, and fit it to our data using the `fit` method:" ] }, { "cell_type": "code", "execution_count": 14, "id": "cc303abf", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "data": { "text/html": [ "
LinearRegression()
In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.
" ], "text/plain": [ "LinearRegression()" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "lin_reg = LinearRegression()\n", "lin_reg.fit(X_train, y_train)" ] }, { "cell_type": "markdown", "id": "4c93426e", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "The `LinearRegression` object after `fit`-ting contains all the coefficients of the regression, which can be accessed using `.coef_` property. In our case, there is just one coefficient, which should be around `-0.017`. It means that prices seem to drop a bit with time, but not too much, around 2 cents per day. We can also access the intersection point of the regression with the Y-axis using `lin_reg`.intercept_` - it will be around `21` in our case, indicating the price at the beginning of the year.\n", "\n", "To see how accurate our model is, we can predict prices on a test dataset, and then measure how close our predictions are to the expected values. This can be done using mean square error (MSE) metrics, which is the mean of all squared differences between expected and predicted values." ] }, { "cell_type": "code", "execution_count": 15, "id": "cd436fa3", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Mean error: 2.77 (17.2%)\n" ] } ], "source": [ "pred = lin_reg.predict(X_test)\n", "\n", "mse = np.sqrt(mean_squared_error(y_test, pred))\n", "print(f\"Mean error: {mse:3.3} ({mse/np.mean(pred)*100:3.3}%)\")" ] }, { "cell_type": "markdown", "id": "43014d87", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Our error seems to be around 2 points, which is ~17%. Not too good. Another indicator of model quality is the **coefficient of determination**, which can be obtained like this:" ] }, { "cell_type": "code", "execution_count": 16, "id": "2803d5e1", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Model determination: 0.04460606335028361\n" ] } ], "source": [ "score = lin_reg.score(X_train, y_train)\n", "print(\"Model determination: \", score)" ] }, { "cell_type": "markdown", "id": "d64ea140", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "If the value is 0, it means that the model does not take input data into account, and acts as the *worst linear predictor*, which is simply a mean value of the result. The value of 1 means that we can perfectly predict all expected outputs. In our case, the coefficient is around 0.06, which is quite low.\n", "\n", "We can also plot the test data together with the regression line to better see how regression works in our case:" ] }, { "cell_type": "code", "execution_count": 17, "id": "337d38c5", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiQAAAGdCAYAAAAi3mhQAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAs8klEQVR4nO3df3RU1b338c8kIQliMhBsSGISyNUWDdFUWlBsyg9LBWuD1PbxBwvE2h+KBq9yH4q46sXc9jHSe+tVa0GrFdctVXufR6CgbRRLSLRKUUIuRCygjoCQGFtkkgAJIbOfP2hGhvycyWT2mcz7tVaWK2f2zPmyd+L55Jyzz3YZY4wAAAAsirNdAAAAAIEEAABYRyABAADWEUgAAIB1BBIAAGAdgQQAAFhHIAEAANYRSAAAgHUJtgs4k8/n06FDh5SSkiKXy2W7HAAA0AfGGDU1NSkrK0txccGf73BcIDl06JBycnJslwEAAEJw4MABZWdnB/0+xwWSlJQUSaf+QampqZarAQAAfdHY2KicnBz/cTxYjgskHZdpUlNTCSQAAESZUG+34KZWAABgHYEEAABYRyABAADWEUgAAIB1BBIAAGAdgQQAAFhHIAEAANYRSAAAgHWOezAaYEO7z2ir57AamlqUnpKsiXlpio9jLSUAiJSgzpCUlZVpwoQJSklJUXp6umbPnq3du3f7X29ra9OSJUt00UUXadiwYcrKytJNN92kQ4cOhb1wIFzKa+tUtHyTbnxyi/75+Rrd+OQWFS3fpPLaOtulAUDMCCqQVFZW6o477tCWLVu0ceNGnTx5UldeeaWOHj0qSTp27Jiqq6t13333qbq6WmvWrNGePXs0a9asASke6K/y2jotWF2tOm9LwPZ6b4sWrK4mlABAhLiMMSbUN3/yySdKT09XZWWlJk+e3GWbt956SxMnTtS+ffuUm5vb62c2NjbK7XbL6/Wylg0GVLvPqGj5pk5hpINLUoY7Wa8vuYLLNwDQi/4ev/t1U6vX65UkpaWl9djG5XJp+PDhXb7e2tqqxsbGgC8gErZ6DncbRiTJSKrztmir53DkigKAGBVyIDHGaNGiRSoqKlJBQUGXbVpaWnTPPfdozpw53aalsrIyud1u/1dOTk6oJQFBaWjqPoyE0g4AELqQA0lJSYl27Nih5557rsvX29radMMNN8jn82nFihXdfs7SpUvl9Xr9XwcOHAi1JCAo6SnJYW0HAAhdSNN+Fy5cqPXr16uqqkrZ2dmdXm9ra9N1110nj8ejTZs29XgtKSkpSUlJSaGUAfTLxLw0ZbqTVe9tUVc3UnXcQzIxr/tLkgCA8AjqDIkxRiUlJVqzZo02bdqkvLy8Tm06wsjevXv16quvauTIkWErFgin+DiXlhXnSzoVPk7X8f2y4nxuaAWACAgqkNxxxx1avXq1nn32WaWkpKi+vl719fU6fvy4JOnkyZP6zne+o7ffflu//e1v1d7e7m9z4sSJAfkHAP0xsyBTK+eOV4Y78LJMhjtZK+eO18yCTEuVAUBsCWrar8vV9V+Kq1at0s0336wPP/ywy7MmklRRUaGpU6f2ug+m/cIGntQKAP3T3+N3UPeQ9JZdxowZ02sbwIni41yadB6XFwHAFhbXAwAA1hFIAACAdQQSAABgHYEEAABYRyABAADWEUgAAIB1BBIAAGAdgQQAAFhHIAEAANYRSAAAgHUEEgAAYB2BBAAAWEcgAQAA1hFIAACAdQQSAABgHYEEAABYRyABAADWEUgAAIB1BBIAAGAdgQQAAFhHIAEAANYRSAAAgHUEEgAAYB2BBAAAWEcgAQAA1hFIAACAdQQSAABgHYEEAABYRyABAADWEUgAAIB1BBIAAGAdgQQAAFhHIAEAANYRSAAAgHUEEgAAYF2C7QIABKfdZ7TVc1gNTS1KT0nWxLw0xce5YrYOOBM/HwgWgQSIIuW1dSrdsEt13hb/tkx3spYV52tmQWbM1QFn4ucDoeCSDRAlymvrtGB1dcD/5CWp3tuiBaurVV5bF1N1wJn4+UCoCCRAFGj3GZVu2CXTxWsd20o37FK7r6sWg68OOBM/H+gPAgkQBbZ6Dnf6i/N0RlKdt0VbPYdjog44Ez8f6A8CCRAFGpq6/598KO2ivQ44Ez8f6A8CCRAF0lOSw9ou2uuAM/Hzgf4gkABRYGJemjLdyepu0qRLp2YxTMxLi4k64Ez8fKA/CCRAFIiPc2lZcb4kdfqffcf3y4rzB/w5D06pA87Ezwf6g0ACRImZBZlaOXe8MtyBp7sz3MlaOXd8xJ7v4JQ64Ez8fCBULmOMo+ZfNTY2yu12y+v1KjU11XY5gOM45QmYTqkDzsTPR+zp7/GbJ7UCUSY+zqVJ5420XYZj6oAz8fOBYHHJBgAAWEcgAQAA1nHJJsK4rgoAQGdBnSEpKyvThAkTlJKSovT0dM2ePVu7d+8OaGOM0f3336+srCwNHTpUU6dO1TvvvBPWoqNVeW2dipZv0o1PbtE/P1+jG5/coqLlm1hsCgAQ84IKJJWVlbrjjju0ZcsWbdy4USdPntSVV16po0eP+tv87Gc/00MPPaTHHntMb731ljIyMvT1r39dTU1NYS8+mrACJgAA3evXtN9PPvlE6enpqqys1OTJk2WMUVZWlu666y4tWbJEktTa2qpRo0Zp+fLluvXWW3v9zME47bfdZ1S0fFO3i065dGqO/utLruDyDQAgKvX3+N2vm1q9Xq8kKS3t1GOAPR6P6uvrdeWVV/rbJCUlacqUKXrjjTe6/IzW1lY1NjYGfA02rIAJAEDPQg4kxhgtWrRIRUVFKigokCTV19dLkkaNGhXQdtSoUf7XzlRWVia32+3/ysnJCbUkx2IFTAAAehZyICkpKdGOHTv03HPPdXrN5Qq87GCM6bStw9KlS+X1ev1fBw4cCLUkx2IFTAAAehbStN+FCxdq/fr1qqqqUnZ2tn97RkaGpFNnSjIzP1uvoKGhodNZkw5JSUlKSkoKpYyo0bECZr23RV3dsNNxDwkrYAIAYlVQZ0iMMSopKdGaNWu0adMm5eXlBbyel5enjIwMbdy40b/txIkTqqys1OWXXx6eiqMQK2ACANCzoALJHXfcodWrV+vZZ59VSkqK6uvrVV9fr+PHj0s6danmrrvu0gMPPKC1a9eqtrZWN998s8466yzNmTNnQP4B0YIVMAEA6F5Q0367uw9k1apVuvnmmyWdOotSWlqqJ554Qp9++qkuvfRS/fKXv/Tf+NqbwTjt93Q8qRUAMBj19/jdr+eQDITBHkgAABiMrD6HBAAAIBwIJAAAwDoCCQAAsI5AAgAArCOQAAAA6wgkAADAOgIJAACwjkACAACsI5AAAADrCCQAAMA6AgkAALCOQAIAAKwjkAAAAOsIJAAAwDoCCQAAsI5AAgAArCOQAAAA6wgkAADAugTbBQBAf7T7jLZ6DquhqUXpKcmamJem+DiX7bIABIlAAiBqldfWqXTDLtV5W/zbMt3JWlacr5kFmRYrAxAsLtkAiErltXVasLo6IIxIUr23RQtWV6u8ts5SZQBCQSABEHXafUalG3bJdPFax7bSDbvU7uuqBQAnIpAAiDpbPYc7nRk5nZFU523RVs/hyBUFoF8IJACiTkNT92EklHYA7COQAIg66SnJYW0HwD4CCYCoMzEvTZnuZHU3udelU7NtJualRbIsAP1AIAEQdeLjXFpWnC9JnUJJx/fLivN5HgkQRQgkAKLSzIJMrZw7XhnuwMsyGe5krZw7nueQAFGGB6MBiFozCzL19fwMntQKDAIEEgBRLT7OpUnnjbRdBoB+4pINAACwjkACAACsI5AAAADrCCQAAMA6AgkAALCOQAIAAKwjkAAAAOsIJAAAwDoCCQAAsI5AAgAArCOQAAAA6wgkAADAOgIJAACwLmZW+233GZYoPw39AQBwkpgIJOW1dSrdsEt13hb/tkx3spYV52tmQabFyuygPwAATjPoL9mU19ZpwerqgIOvJNV7W7RgdbXKa+ssVWYH/QEAcKJBHUjafUalG3bJdPFax7bSDbvU7uuqxeBDfwAAnGpQB5KtnsOdzgSczkiq87Zoq+dw5IqyiP4AADjVoA4kDU3dH3xDaRft6A8AgFMFHUiqqqpUXFysrKwsuVwurVu3LuD15uZmlZSUKDs7W0OHDtWFF16olStXhqveoKSnJIe1XbSjPwAAThV0IDl69KgKCwv12GOPdfn63XffrfLycq1evVrvvvuu7r77bi1cuFC///3v+11ssCbmpSnTnazuJrO6dGp2ycS8tEiWZQ39AQBwqqADyVVXXaWf/vSnuvbaa7t8/c0339T8+fM1depUjRkzRj/84Q9VWFiot99+u9/FBis+zqVlxfmS1Okg3PH9suL8mHn+Bv0BAHCqsN9DUlRUpPXr1+vgwYMyxqiiokJ79uzRjBkzwr2rPplZkKmVc8crwx14GSLDnayVc8fH3HM36A8AgBOF/cFojz76qH7wgx8oOztbCQkJiouL01NPPaWioqIu27e2tqq1tdX/fWNjY7hL0syCTH09P4Mnk/4D/QEAcJoBCSRbtmzR+vXrNXr0aFVVVen2229XZmampk+f3ql9WVmZSktLw11GJ/FxLk06b+SA7yda0B8AACdxGWNCfgqWy+XS2rVrNXv2bEnS8ePH5Xa7tXbtWl199dX+dt///vf10Ucfqby8vNNndHWGJCcnR16vV6mpqaGWBgAAIqixsVFutzvk43dYz5C0tbWpra1NcXGBt6bEx8fL5/N1+Z6kpCQlJSWFswwAABBlgg4kzc3Neu+99/zfezwe1dTUKC0tTbm5uZoyZYoWL16soUOHavTo0aqsrNR//dd/6aGHHgpr4QAAYPAI+pLN5s2bNW3atE7b58+fr2eeeUb19fVaunSpXnnlFR0+fFijR4/WD3/4Q919991yuXq/abK/p3wAAEDk9ff43a97SAYCgQQAgOjT3+P3oF7LBgAARAcCCQAAsI5AAgAArCOQAAAA6wgkAADAOgIJAACwjkACAACsI5AAAADrCCQAAMA6AgkAALCOQAIAAKwjkAAAAOsIJAAAwDoCCQAAsI5AAgAArCOQAAAA6wgkAADAOgIJAACwjkACAACsI5AAAADrCCQAAMA6AgkAALCOQAIAAKwjkAAAAOsIJAAAwDoCCQAAsI5AAgAArEuwXUAklf3xXT1R+YH/+/SUJF2SO1yX5I7QJTnDdXH2cA1NjLdYIQAAsSlmAokxRi9sOxiwraGpVS+/87Fefufjbt/nHjrkVGjJGaFLcoerMGe43EOHDHS5AADEFJcxxtgu4nSNjY1yu93yer1KTU0N62fXe1v0s/K/avuBI/L87WjYPjd5SJw/sFySe+q/55ydFLbPBwDA6fp7/I6pQNITn8/o/U+atX3/EW0/8Kmq9x3R7o+bwvb58XEuXZIzXJfkDlfBuW4lD4nXjHEZYft8AABsIpBEkDFG+/5+TNsPfHoquOw/op0HvQO2v9JZ4zTvstGKi3MN2D4AAAgHAokD1XmP/yOw/CO4HDiidt/AdHPJtPO18GvnKymBm3EBAPYQSKLUlg/+rruer1F9Y8uA7+uGCTm69+oLlZrMzbgAgIFBIBnE/lrfqMX/d8eAXhbqkJqcoA0LizR65LAB3xcAYPAhkMS4g0eOa+manara80lE9re+5Cu6OHt4RPYFAIgeBBL06tOjJ/STF3dpzfaDvTcOg6dv/rKuuGBURPYFAHAGAgnCorn1pBas3qbX9v4tIvt74FsXac6luRHZFwBg4BFIoky7z2ir57AamlqUnpKsiXlpio+Sab1t7T7dt65Wz791ICL7u33qeVo8Y6xcrujoHwCIZQSSKFJeW6fSDbtU5/1sZk2mO1nLivM1syDTYmXhY4zRI3/aq4df3RuR/Q0dEq/t//p1JQ9h2jMA2EQgiRLltXVasLpaZ3Z2x9/+K+eOHzShpC9+99Z+LXlhZ8T2t+3H0zWSx/kDwIAhkESBdp9R0fJNAWdGTueSlOFO1utLroiayzeRsHl3g25e9VbE9vfqoik6P/3siO0PAAaT/h6/Y2a1X5u2eg53G0YkyUiq87Zoq+ewJp03MnKFOdzUsen68MGre2xT8dcGffeZ8ISW6Q9V9trm2e9fqsvPPycs+wMAfIZAEgENTX17Gmtf2+Ez0y7oPbTsOtSobzz6Wlj2N+epv/Ta5mffvljXTcgJy/4AIFYQSCIgPSU5rO0QnPys1F5Dy8eNLbr0gT+FZX8/emGHfvTCjh7b3D71PP1o5gVh2R8ADAbcQxIBHfeQ1HtbOt3UKnEPSbQ42npSF5e+MmALJZ7pigvS9ev5X+407Tmap44DGLy4qTVKdMyykRQQSmJ1ls1g1e4zKv7F69pV1xixfWakJun+WeP4+QFgFYEkisTCc0jQN//8/Hb9vuZQxPb3P/96pdxnsdozgIFDIIkynG5HX/3iT3v18417Ira/zf97qsacw2rPAEJDIAFiwJvv/103PrnFdhn6v7dN0oQxabbLAOBAEX8OSVVVlf793/9d27ZtU11dndauXavZs2cHtHn33Xe1ZMkSVVZWyufzady4cfrv//5v5eaymBoQir5OCX/khi/qmi+eG7Btywd/1w2/Ck+Y+V+Pv9lrm5/MLtC8y0aHZX8AYkfQgeTo0aMqLCzUd7/7XX3729/u9Pr777+voqIife9731NpaancbrfeffddJSczpRUIVX+mjl/2TyN7nfb8wSfNuuLnvT8Yri/uW1er+9bV9thm7mW5+unsi8KyPwCDQ78u2bhcrk5nSG644QYNGTJEv/nNb0L6TC7ZoCuxfu+NE6aOf3r0hC75ycYB+eyuFOYM1+/v+ErE9gegfxz16Hifz6eXXnpJP/rRjzRjxgxt375deXl5Wrp0aafLOh1aW1vV2trq/76xMXLTJREdmJ0kxce5tKw4X7f9Y+r4mYykZcX5AxrSRgxL9J9p6W5M7v3GBVr4XE1Y9vc/B45ozD0v9dru/Qe+EVPhFBiswhpIGhoa1NzcrAcffFA//elPtXz5cpWXl+vaa69VRUWFpkyZ0uk9ZWVlKi0tDWcZGES6WyW53tuiBaureX6LBT2NyZ3P1ejxXsbEGKOi5RU6eOR4WOo5794/9NqmtnSGzk7iwdSAk4X1ks2hQ4d07rnn6sYbb9Szzz7rbzdr1iwNGzZMzz33XKfP6OoMSU5ODpdswCrJp3FKX0Syjtt+s03l79T36zOC8fqSacoecVbE9gcMNo66ZHPOOecoISFB+fn5AdsvvPBCvf76612+JykpSUlJSeEsA4MEqyR/xil9Eck6Hp/3pV7bPPzqHj386t5+7adD0fKKXtusuf1yjc8dEZb9AQgU1kCSmJioCRMmaPfu3QHb9+zZo9GjmQaI4LBK8mec0hdOqaPDXdO/oLumf6HHNhv+55AWPrc9LPu7dsUbvbZ5bM4l+ubFWWHZHxBLgg4kzc3Neu+99/zfezwe1dTUKC0tTbm5uVq8eLGuv/56TZ48WdOmTVN5ebk2bNigzZs3h7NuxABWSf6MU/rCKXUEo7gwS8WFPQeEmgNHNPuXfw7L/kqe3a6SZ3sOQPd9M1/fK8oLy/6AwSLoe0g2b96sadOmddo+f/58PfPMM5Kkp59+WmVlZfroo480duxYlZaW6pprrunT5zPtFx2cMNXVKZzSF06pw4ZDR47r8gc3RWx/N07MUdm1F0dsf0B/8eh4DGqskvwZp/SFU+pwomMnTir/X1+O2P4Kzk3VhpIiuVyDK/whOhFIMOjxHJLPOKUvnFJHNPL5jK565DXt/rgpIvubfuEorZw7XkPi4yKyP8QuAgliQqw/qfV0TukLp9QxWN3x22q9tLMuIvsaOypF/2/BJKUkD4nI/jA4EUgAIEY9+qe9emjjnojsa1hivP70L1OV4XbODctwFgIJAKBbb314uE+rNIdL+V1f1QUZ/L87FhFIAAD9Es7VnnsSH+fSLV8Zo/mXj+GpuIMQgQQxgfsVPkNfOJNTxmWg6vhbc6uueezPYVuDqCfXfzlHtxTlaWxGyoDvC+FDIMGgx4yOz9AXzuSUcbFdx7ETJzX3qb+oev+RAd/XNy/O1C1FeTzK30EIJBjUultZNhafeUFfOJNTxsUpdfTG5zOq2vuJnv7zh6ra88mA7mvyFz6nW74yRpM//znFcRZxwBFIMGg5ZYVbJ6AvnMkp4+KUOsKpev+nevp1j17cMbBTn7+YM1y3FOXpqoIMntXST45a7RcIJ6escOsE9IUzOWVcnFJHOI3PHaHxc0bosTndt9nzcZOeeeNDPfuX/SHvp+bAEd3Zy+KL//S5YbrlK3n69vhsDU2MD3lf6BmBBI7ltJVlbaIvnMkp4+KUOiLtC6NS9MC3LtID37qo2zYHjxzXb97cp6f/7NGJk76Q9vPBJ0f143W1+vG62m7bnHN2km4pGqM5E3M1/KzEkPYT6wgkcKxoXFl2oNAXzuSUcXFKHU507vChuueqC3TPVRd02+bvza169i/79fSfPfr0WFtI+/lbc6t+Vr5bPyvf3W2bsxLj9dB1hSr6/Od0dhKH3zPRI3CsiXlpynQn97qy7MS8tEiXFnH0hTM5ZVycUke0Gnl2khZ+7fNa+LXPd9umufWkXtj2kZ7+s0f7/n4spP0cO9Gu2/6xMOWZzk5K0CW5w3VJznBdkjtCX8wZrhHDYutMCze1wtFYWfYz9IUzOWVcnFJHLGs92a4/7KzT069/qJ0HvWH//CHxLl2SM+JUcMk9FVxGpTrnrBezbDDo2X62gpPQF87klHFxSh3onjFGnr8d1fb9R7T9wKfavv+I3jnUGNZ9FGa7dUnuP4JLzgjlpA2VyzXws6sIJIgJTnkKphPQF87klHFxSh0InTFGh7wt2r7/VGDZvv9TbT9wROE8Wv9kdoHmXTY6fB8oAgkAADGpobFF2w8cCQgtwcwk+vDBq8NaD88hAQAgBqWnJmvGuAzNGJfRbZsjx06opiO0HDgVXJpaTuraS86NYKV9QyABAGCQGn5WoqaOTdfUsem2S+kVz8kFAADWEUgAAIB1BBIAAGAdgQQAAFhHIAEAANYRSAAAgHUEEgAAYB2BBAAAWEcgAQAA1hFIAACAdQQSAABgHYEEAABYRyABAADWsdovIKndZ7TVc1gNTS1KT0nWxLw0xce5bJcFADGDQIKYV15bp9INu1TnbfFvy3Qna1lxvmYWZFqsDABiB5dsENPKa+u0YHV1QBiRpHpvixasrlZ5bZ2lygAgthBIELPafUalG3bJdPFax7bSDbvU7uuqBQAgnAgkiFlbPYc7nRk5nZFU523RVs/hyBUFADGKQIKY1dDUfRgJpR0AIHQEEsSs9JTksLYDAISOQIKYNTEvTZnuZHU3udelU7NtJualRbIsAIhJBBLErPg4l5YV50tSp1DS8f2y4nyeRwIAEUAgQUybWZCplXPHK8MdeFkmw52slXPH8xwSAIgQHoyGmDezIFNfz8/gSa0AYBGBBNCpyzeTzhtpuwwAiFlcsgEAANYRSAAAgHVcsgGiDCsTO5NTxoU60BMnjwuBBIgirEzsTE4ZF+pAT5w+Li5jjKNWDmtsbJTb7ZbX61VqaqrtcgDH6FiZ+Mxf2I6/bZimbIdTxoU60JNIjEt/j9/cQwJEAVYmdianjAt1oCfRMi5BB5KqqioVFxcrKytLLpdL69at67btrbfeKpfLpYcffrgfJQJgZWJncsq4UAd6Ei3jEnQgOXr0qAoLC/XYY4/12G7dunX6y1/+oqysrJCLA3AKKxM7k1PGhTrQk2gZl6Bvar3qqqt01VVX9djm4MGDKikp0csvv6yrr7465OIAnMLKxM7klHGhDvQkWsYl7PeQ+Hw+zZs3T4sXL9a4ceN6bd/a2qrGxsaALwCBWJnYmZwyLtSBnkTLuIQ9kCxfvlwJCQm68847+9S+rKxMbrfb/5WTkxPukoCox8rEzuSUcaEO9CRaxiWsgWTbtm165JFH9Mwzz8jl6ts/bOnSpfJ6vf6vAwcOhLMkYNBgZWJncsq4UAd6Eg3j0q/nkLhcLq1du1azZ8+WJD388MNatGiR4uI+yznt7e2Ki4tTTk6OPvzww14/k+eQAD1z8pMWY5lTxoU60JOBHJf+Hr/D+qTWefPmafr06QHbZsyYoXnz5um73/1uOHcFxCxWJnYmp4wLdaAnTh6XoANJc3Oz3nvvPf/3Ho9HNTU1SktLU25urkaODPyHDhkyRBkZGRo7dmz/qwUAAINS0IHk7bff1rRp0/zfL1q0SJI0f/58PfPMM2ErDAAAxI6gA8nUqVMVzG0nfblvBED04R4BAOHEar8Agub0VUMBRB8W1wMQlI5VQ89cG6Pe26IFq6tVXltnqTIA0YxAAqDPomXVUADRh0ACoM+iZdVQANGHQAKgz6Jl1VAA0YdAAqDPomXVUADRh1k2APqsY9XQem9Ll/eRuHRqbQzbq4bCPqaFB6I/ekcgAdBnHauGLlhdLZcUEEqctGoo7GJaeCD6o2+4ZAMgKNGwaijsYVp4IPqj7/q12u9AYLVfIDpwChpnavcZFS3f1O1MrI5Leq8vuSImflZirT8ctdovgNjh5FVDYUcw08Jj4WeH/ggOl2wAAGHBtPBA9EdwCCQAgLBgWngg+iM4BBIAQFh0TAvv7m4Il07NLomVaeH0R3AIJACAsOiYFi6p00E4FqeF0x/BIZAAAMKGaeGB6I++Y9ovACDsmBYeKBb6g2m/AADHYVp4IPqjd1yyAQAA1hFIAACAdQQSAABgHYEEAABYRyABAADWEUgAAIB1BBIAAGAdgQQAAFhHIAEAANYRSAAAgHUEEgAAYB2BBAAAWEcgAQAA1hFIAACAdQQSAABgHYEEAABYRyABAADWEUgAAIB1BBIAAGAdgQQAAFhHIAEAANYRSAAAgHUEEgAAYB2BBAAAWEcgAQAA1hFIAACAdQQSAABgHYEEAABYRyABAADWEUgAAIB1BBIAAGBd0IGkqqpKxcXFysrKksvl0rp16/yvtbW1acmSJbrooos0bNgwZWVl6aabbtKhQ4fCWTMAABhkgg4kR48eVWFhoR577LFOrx07dkzV1dW67777VF1drTVr1mjPnj2aNWtWWIoFAACDk8sYY0J+s8ultWvXavbs2d22eeuttzRx4kTt27dPubm5vX5mY2Oj3G63vF6vUlNTQy0NAABEUH+P3wkDUFMAr9crl8ul4cOHd/l6a2urWltb/d83NjYOdEkAAMBhBvSm1paWFt1zzz2aM2dOt2mprKxMbrfb/5WTkzOQJQEAAAcasEDS1tamG264QT6fTytWrOi23dKlS+X1ev1fBw4cGKiSAACAQw3IJZu2tjZdd9118ng82rRpU4/XkpKSkpSUlDQQZQAAgCgR9kDSEUb27t2riooKjRw5Mty7AAAAg0zQgaS5uVnvvfee/3uPx6OamhqlpaUpKytL3/nOd1RdXa0XX3xR7e3tqq+vlySlpaUpMTExfJUDAIBBI+hpv5s3b9a0adM6bZ8/f77uv/9+5eXldfm+iooKTZ06tdfPZ9ovAADRJ+LTfqdOnaqeMkw/HmsCAABiFGvZAAAA6wgkAADAOgIJAACwjkACAACsI5AAAADrCCQAAMA6AgkAALCOQAIAAKwjkAAAAOsIJAAAwDoCCQAAsI5AAgAArCOQAAAA6wgkAADAOgIJAACwjkACAACsI5AAAADrCCQAAMA6AgkAALCOQAIAAKwjkAAAAOsSbBcAAIPBiZM+/ebND7Xv8DGNTjtL8yaNUWJC5P/ma/cZbfUcVkNTi9JTkjUxL03xca6I1wEEi0ACAP1U9oddevI1j3zms23/5w/v6gdfzdPSb+RHrI7y2jqVbtilOm+Lf1umO1nLivM1syAzYnUAoeCSDQD0Q9kfdumJqsAwIkk+Iz1R5VHZH3ZFpI7y2jotWF0dEEYkqd7bogWrq1VeWxeROoBQEUgAIEQnTvr05GueHts8+ZpHJ076BrSOdp9R6YZdMl281rGtdMMutZ+ZmgAHIZAAQIh+8+aHnc6MnMlnTrUbSFs9hzudGTmdkVTnbdFWz+EBrQPoDwIJAIRo3+FjYW0Xqoam7sNIKO0AGwgkABCi0WlnhbVdqNJTksPaDrCBQAIAIZo3aYx6m1Eb5zrVbiBNzEtTpjtZ3ZXi0qnZNhPz0ga0DqA/CCQAEKLEhDj94Kt5Pbb5wVfzBvx5JPFxLi0rPjW9+MxQ0vH9suJ8nkcCRyOQAEA/LP1Gvm6dnNfpTEmcS7p1cuSeQzKzIFMr545XhjvwskyGO1kr547nOSRwPJcxxlHzwBobG+V2u+X1epWammq7HADoE57UiljX3+M3T2oFgDBITIjT9776T7bLUHycS5POG2m7DCBoXLIBAADWEUgAAIB1BBIAAGAdgQQAAFhHIAEAANYRSAAAgHUEEgAAYB2BBAAAWEcgAQAA1jnuSa0dT7JvbGy0XAkAAOirjuN2qCvSOC6QNDU1SZJycnIsVwIAAILV1NQkt9sd9Psct7iez+fToUOHlJKSIpfLpcbGRuXk5OjAgQMstmcJY2AfY2AfY2AfY2BXb/1vjFFTU5OysrIUFxf8HSGOO0MSFxen7OzsTttTU1P5AbSMMbCPMbCPMbCPMbCrp/4P5cxIB25qBQAA1hFIAACAdY4PJElJSVq2bJmSkpJslxKzGAP7GAP7GAP7GAO7Brr/HXdTKwAAiD2OP0MCAAAGPwIJAACwjkACAACsI5AAAADrIh5IysrKNGHCBKWkpCg9PV2zZ8/W7t27u21/6623yuVy6eGHHw7Y3traqoULF+qcc87RsGHDNGvWLH300UcDXP3g0NcxePfddzVr1iy53W6lpKTosssu0/79+/2vMwah68sYNDc3q6SkRNnZ2Ro6dKguvPBCrVy5MqANYxC6lStX6uKLL/Y/5GnSpEn64x//6H/dGKP7779fWVlZGjp0qKZOnap33nkn4DPo//7paQza2tq0ZMkSXXTRRRo2bJiysrJ000036dChQwGfwRj0T2+/B6cb8OOxibAZM2aYVatWmdraWlNTU2Ouvvpqk5uba5qbmzu1Xbt2rSksLDRZWVnmP//zPwNeu+2228y5555rNm7caKqrq820adNMYWGhOXnyZIT+JdGrL2Pw3nvvmbS0NLN48WJTXV1t3n//ffPiiy+ajz/+2N+GMQhdX8bg+9//vjnvvPNMRUWF8Xg85oknnjDx8fFm3bp1/jaMQejWr19vXnrpJbN7926ze/duc++995ohQ4aY2tpaY4wxDz74oElJSTEvvPCC2blzp7n++utNZmamaWxs9H8G/d8/PY3BkSNHzPTp083vfvc789e//tW8+eab5tJLLzVf+tKXAj6DMeif3n4POkTieBzxQHKmhoYGI8lUVlYGbP/oo4/Mueeea2pra83o0aMDOuDIkSNmyJAh5vnnn/dvO3jwoImLizPl5eWRKn3Q6GoMrr/+ejN37txu38MYhFdXYzBu3Djzb//2bwHtxo8fb3784x8bYxiDgTBixAjz1FNPGZ/PZzIyMsyDDz7of62lpcW43W7z+OOPG2Po/4HSMQZd2bp1q5Fk9u3bZ4xhDAbKmWMQqeOx9XtIvF6vJCktLc2/zefzad68eVq8eLHGjRvX6T3btm1TW1ubrrzySv+2rKwsFRQU6I033hj4ogeZM8fA5/PppZde0he+8AXNmDFD6enpuvTSS7Vu3Tr/exiD8Orq96CoqEjr16/XwYMHZYxRRUWF9uzZoxkzZkhiDMKpvb1dzz//vI4ePapJkybJ4/Govr4+oG+TkpI0ZcoUf9/S/+F15hh0xev1yuVyafjw4ZIYg3DragwieTy2GkiMMVq0aJGKiopUUFDg3758+XIlJCTozjvv7PJ99fX1SkxM1IgRIwK2jxo1SvX19QNa82DT1Rg0NDSoublZDz74oGbOnKlXXnlF3/rWt3TttdeqsrJSEmMQTt39Hjz66KPKz89Xdna2EhMTNXPmTK1YsUJFRUWSGINw2Llzp84++2wlJSXptttu09q1a5Wfn+/vv1GjRgW0P71v6f/w6G4MztTS0qJ77rlHc+bM8S/sxhiER09jEMnjsdXVfktKSrRjxw69/vrr/m3btm3TI488ourqarlcrqA+zxgT9HtiXVdj4PP5JEnXXHON7r77bknSF7/4Rb3xxht6/PHHNWXKlG4/jzEIXldjIJ0KJFu2bNH69es1evRoVVVV6fbbb1dmZqamT5/e7ecxBn03duxY1dTU6MiRI3rhhRc0f/58f+iW1Kkf+9K39H9wuhuD00NJW1ubbrjhBvl8Pq1YsaLXz2QMgtPdGBw/fjyix2NrZ0gWLlyo9evXq6KiQtnZ2f7tr732mhoaGpSbm6uEhAQlJCRo3759+pd/+ReNGTNGkpSRkaETJ07o008/DfjMhoaGTn/RoHvdjcE555yjhISETn+lXHjhhf5ZNoxBeHQ3BsePH9e9996rhx56SMXFxbr44otVUlKi66+/Xv/xH/8hiTEIh8TERJ1//vn68pe/rLKyMhUWFuqRRx5RRkaGJHX6C+/0vqX/w6O7MejQ1tam6667Th6PRxs3bgxY9p4xCI/uxiDSx+OIBxJjjEpKSrRmzRpt2rRJeXl5Aa/PmzdPO3bsUE1Njf8rKytLixcv1ssvvyxJ+tKXvqQhQ4Zo48aN/vfV1dWptrZWl19+eUT/PdGotzFITEzUhAkTOk1D3bNnj0aPHi2JMeiv3sagra1NbW1tiosL/BWNj4/3n8FiDMLPGKPW1lbl5eUpIyMjoG9PnDihyspKf9/S/wOjYwykz8LI3r179eqrr2rkyJEBbRmDgdExBhE/Hgd7921/LViwwLjdbrN582ZTV1fn/zp27Fi37znzrl5jTk0zys7ONq+++qqprq42V1xxBVO9+qgvY7BmzRozZMgQ86tf/crs3bvX/OIXvzDx8fHmtdde87dhDELXlzGYMmWKGTdunKmoqDAffPCBWbVqlUlOTjYrVqzwt2EMQrd06VJTVVVlPB6P2bFjh7n33ntNXFyceeWVV4wxp6b9ut1us2bNGrNz505z4403djntl/4PXU9j0NbWZmbNmmWys7NNTU1NwO9Ja2ur/zMYg/7p7ffgTAN5PI54IJHU5deqVau6fU9XHXD8+HFTUlJi0tLSzNChQ803v/lNs3///oEtfpDo6xj8+te/Nueff75JTk42hYWFAc+/MIYx6I++jEFdXZ25+eabTVZWlklOTjZjx441P//5z43P5/O3YQxCd8stt5jRo0ebxMRE87nPfc587WtfC/ifsM/nM8uWLTMZGRkmKSnJTJ482ezcuTPgM+j//ulpDDweT7e/JxUVFf7PYAz6p7ffgzMN5PHYZYwxoZzSAQAACBfrzyEBAAAgkAAAAOsIJAAAwDoCCQAAsI5AAgAArCOQAAAA6wgkAADAOgIJAACwjkACAACsI5AAAADrCCQAAMA6AgkAALDu/wMPmoSHGBASlgAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.scatter(X_test, y_test)\n", "plt.plot(X_test, pred)" ] }, { "cell_type": "markdown", "id": "9f93331c", "metadata": {}, "source": [ "## Polynomial Regression" ] }, { "cell_type": "markdown", "id": "7d7f9463", "metadata": {}, "source": [ "

\n", "\n", "A demo of logistic-regression. [source]\n", "

" ] }, { "cell_type": "markdown", "id": "c2831e2f", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Another type of Linear Regression is Polynomial Regression. While sometimes there's a linear relationship between variables - the bigger the pumpkin in volume, the higher the price - sometimes these relationships can't be plotted as a plane or straight line.\n", "\n", ":::{seealso}\n", "Here are [some more examples](https://online.stat.psu.edu/stat501/lesson/9/9.8) of data that could use Polynomial Regression\n", ":::\n", "\n", "Take another look at the relationship between Date and Price. Does this scatterplot seem like it should necessarily be analyzed by a straight line? Can't prices fluctuate? In this case, we can try polynomial regression.\n", "\n", ":::{note}\n", "Polynomials are mathematical expressions that might consist of one or more variables and coefficients.\n", ":::\n", "\n", "Polynomial regression creates a curved line to better fit nonlinear data. In our case, if we include a squared `DayOfYear` variable in input data, we should be able to fit our data with a parabolic curve, which will have a minimum at a certain point within the year.\n", "\n", "Scikit-learn includes a helpful [pipeline API](https://scikit-learn.org/stable/modules/generated/sklearn.pipeline.make_pipeline.html?highlight=pipeline#sklearn.pipeline.make_pipeline) to combine different steps of data processing together. A **pipeline** is a chain of **estimators**. In our case, we will create a pipeline that first adds polynomial features to our model, and then trains the regression:" ] }, { "cell_type": "code", "execution_count": 18, "id": "bf0c99b8", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "data": { "text/html": [ "
Pipeline(steps=[('polynomialfeatures', PolynomialFeatures()),\n",
       "                ('linearregression', LinearRegression())])
In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.
" ], "text/plain": [ "Pipeline(steps=[('polynomialfeatures', PolynomialFeatures()),\n", " ('linearregression', LinearRegression())])" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from sklearn.preprocessing import PolynomialFeatures\n", "from sklearn.pipeline import make_pipeline\n", "\n", "pipeline = make_pipeline(PolynomialFeatures(2), LinearRegression())\n", "\n", "pipeline.fit(X_train, y_train)" ] }, { "cell_type": "markdown", "id": "5181be18", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Using `PolynomialFeatures(2)` means that we will include all second-degree polynomials from the input data. In our case, it will just mean $DayOfYear^2$, but given two input variables X and Y, this will add $X^2$, XY and $Y^2$. We may also use higher-degree polynomials if we want.\n", "\n", "Pipelines can be used in the same manner as the original `LinearRegression` object, i.e. we can `fit` the pipeline, and then use `predict` to get the prediction results. Here is the graph showing test data, and the approximation curve:\n", "\n", ":::{figure} https://static-1300131294.cos.ap-shanghai.myqcloud.com/images/ml-regression/poly-results.png\n", "---\n", "name: 'Polynomial regression'\n", "width: 70%\n", "---\n", "Polynomial regression{cite}`Polynomial_regression`\n", ":::\n", "\n", "Using Polynomial Regression, we can get slightly lower MSE and higher determination, but not significantly. We need to take into account other features!\n", "\n", ":::{seealso}\n", "We can see that the minimal pumpkin prices are observed somewhere around Halloween. How can we explain this? \n", ":::\n", "\n", "🎃 Congratulations, we just created a model that can help predict the price of pie pumpkins. We can probably repeat the same procedure for all pumpkin types, but that would be tedious. Let's learn now how to take pumpkin variety into account in our model!\n", "\n", "## Categorical Features\n", "\n", "In the ideal world, we want to be able to predict prices for different pumpkin varieties using the same model. However, the `Variety` column is somewhat different from columns like `Month`, because it contains non-numeric values. Such columns are called **categorical**.\n", "\n", "Here we can see how the average price depends on variety:\n", "\n", ":::{figure} https://static-1300131294.cos.ap-shanghai.myqcloud.com/images/ml-regression/price-by-variety.png\n", "---\n", "name: 'Average price by variety'\n", "width: 70%\n", "---\n", "Average price by variety{cite}`Average_price_by_variety`\n", ":::\n", "\n", "To take variety into account, we first need to convert it to numeric form or **encode**** it. There are several ways we can do it:\n", "\n", "* Simple **numeric encoding** will build a table of different varieties, and then replace the variety name by an index in that table. This is not the best idea for linear regression, because linear regression takes the actual numeric value of the index, and adds it to the result, multiplying by some coefficient. In our case, the relationship between the index number and the price is clearly non-linear, even if we make sure that indices are ordered in some specific way.\n", "* **One-hot encoding** will replace the `Variety` column by 4 different columns, one for each variety. Each column will contain `1` if the corresponding row is of a given variety and `0`` otherwise. This means that there will be four coefficients in linear regression, one for each pumpkin variety, responsible for the \"starting price\" (or rather \"additional price\") for that particular variety.\n", "\n", "The code below shows how we can one-hot encode a variety:" ] }, { "cell_type": "code", "execution_count": 19, "id": "246bfdd0", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
FAIRYTALEMINIATUREMIXED HEIRLOOM VARIETIESPIE TYPE
700001
710001
720001
730001
740001
...............
17380100
17390100
17400100
17410100
17420100
\n", "

415 rows × 4 columns

\n", "
" ], "text/plain": [ " FAIRYTALE MINIATURE MIXED HEIRLOOM VARIETIES PIE TYPE\n", "70 0 0 0 1\n", "71 0 0 0 1\n", "72 0 0 0 1\n", "73 0 0 0 1\n", "74 0 0 0 1\n", "... ... ... ... ...\n", "1738 0 1 0 0\n", "1739 0 1 0 0\n", "1740 0 1 0 0\n", "1741 0 1 0 0\n", "1742 0 1 0 0\n", "\n", "[415 rows x 4 columns]" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "pd.get_dummies(new_pumpkins[\"Variety\"])" ] }, { "cell_type": "markdown", "id": "968fa4ea", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "To train linear regression using one-hot encoded variety as input, we just need to initialize `X` and `y` data correctly:" ] }, { "cell_type": "code", "execution_count": 20, "id": "f0bc3720", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "X = pd.get_dummies(new_pumpkins[\"Variety\"])\n", "y = new_pumpkins[\"Price\"]" ] }, { "cell_type": "markdown", "id": "cc63ccc2", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "The rest of the code is the same as what we used above to train Linear Regression. If we try it, we will see that the mean squared error is about the same, but we get the much higher coefficient of determination (~77%). To get even more accurate predictions, we can take more categorical features into account, as well as numeric features, such as `Month` or `DayOfYear`. To get one large array of features, we can use `join`:" ] }, { "cell_type": "code", "execution_count": 21, "id": "6fd52a7f", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [], "source": [ "X = (\n", " pd.get_dummies(new_pumpkins[\"Variety\"])\n", " .join(new_pumpkins[\"Month\"])\n", " .join(pd.get_dummies(new_pumpkins[\"City\"]))\n", " .join(pd.get_dummies(new_pumpkins[\"Package\"]))\n", ")\n", "y = new_pumpkins[\"Price\"]" ] }, { "cell_type": "markdown", "id": "e113e574", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "Here we also take into account `City` and `Package` type, which gives us MSE 2.84 (10%), and determination 0.94!\n", "\n", "## Putting it all together\n", "\n", "To make the best model, we can use combined (one-hot encoded categorical + numeric) data from the above example together with Polynomial Regression. Here is the complete code for our convenience:" ] }, { "cell_type": "code", "execution_count": 22, "id": "fffe46b9", "metadata": { "attributes": { "classes": [ "code-cell" ], "id": "" }, "slideshow": { "slide_type": "" }, "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Mean error: 2.23 (8.25%)\n", "Model determination: 0.9652869939759445\n" ] } ], "source": [ "# set up training data\n", "X = (\n", " pd.get_dummies(new_pumpkins[\"Variety\"])\n", " .join(new_pumpkins[\"Month\"])\n", " .join(pd.get_dummies(new_pumpkins[\"City\"]))\n", " .join(pd.get_dummies(new_pumpkins[\"Package\"]))\n", ")\n", "y = new_pumpkins[\"Price\"]\n", "\n", "# make train-test split\n", "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=0)\n", "\n", "# setup and train the pipeline\n", "pipeline = make_pipeline(PolynomialFeatures(2), LinearRegression())\n", "pipeline.fit(X_train, y_train)\n", "\n", "# predict results for test data\n", "pred = pipeline.predict(X_test)\n", "\n", "# calculate MSE and determination\n", "mse = np.sqrt(mean_squared_error(y_test, pred))\n", "print(f\"Mean error: {mse:3.3} ({mse/np.mean(pred)*100:3.3}%)\")\n", "\n", "score = pipeline.score(X_train, y_train)\n", "print(\"Model determination: \", score)" ] }, { "cell_type": "markdown", "id": "28aa766b", "metadata": { "slideshow": { "slide_type": "" }, "tags": [] }, "source": [ "This should give us the best determination coefficient of almost 97%, and MSE=2.23 (~8% prediction error).\n", "\n", "| Model | MSE | Determination |\n", "|-------|-----|---------------|\n", "| `DayOfYear` Linear | 2.77 (17.2%) | 0.07 |\n", "| `DayOfYear` Polynomial | 2.73 (17.0%) | 0.08 |\n", "| `Variety` Linear | 5.24 (19.7%) | 0.77 |\n", "| All features Linear | 2.84 (10.5%) | 0.94 |\n", "| All features Polynomial | 2.23 (8.25%) | 0.97 |\n", "\n", "🏆 Well done! We created four Regression models in one section and improved the model quality to 97%. In the final section on Regression, we will learn about Logistic Regression to determine categories.\n", "\n", "## Self study\n", "\n", "In this section, we learned about Linear Regression. There are other important types of Regression. Read about Stepwise, Ridge, Lasso and Elasticnet techniques. A good course to study to learn more is the [Stanford Statistical Learning course](https://online.stanford.edu/courses/sohs-ystatslearning-statistical-learning)\n", "\n", "## Your turn! 🚀\n", "\n", "Test several different variables in this notebook to see how correlation corresponds to model accuracy.\n", "\n", "Assignment - [Create a regression model](../../assignments/ml-fundamentals/create-a-regression-model.md)\n", "\n", "## Acknowledgments\n", "\n", "Thanks to Microsoft for creating the open-source course [ML-For-Beginners](https://github.com/microsoft/ML-For-Beginners). It inspires the majority of the content in this chapter.\n", "\n", "---\n", "\n", ":::{bibliography}\n", ":filter: docname in docnames\n", ":::" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.16" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": true, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false } }, "nbformat": 4, "nbformat_minor": 5 }