diff --git a/notebooks/statsmodels_TP.ipynb b/notebooks/statsmodels_TP.ipynb
index d876b1a113039f10883c8cec8b9226493fe0527f..256c2a7dcdc33957d374695e65f5a77deee3b38d 100644
--- a/notebooks/statsmodels_TP.ipynb
+++ b/notebooks/statsmodels_TP.ipynb
@@ -2,8 +2,8 @@
"cells": [
{
"cell_type": "code",
- "execution_count": null,
- "id": "d5065981",
+ "execution_count": 1,
+ "id": "4e16caf7",
"metadata": {},
"outputs": [],
"source": [
@@ -22,7 +22,7 @@
},
{
"cell_type": "markdown",
- "id": "4ddd8902",
+ "id": "358dce7a",
"metadata": {},
"source": [
"# Multi-way ANOVA"
@@ -30,7 +30,7 @@
},
{
"cell_type": "markdown",
- "id": "9c4680b9",
+ "id": "a1face9f",
"metadata": {
"heading_collapsed": true
},
@@ -42,7 +42,7 @@
},
{
"cell_type": "markdown",
- "id": "d2d4fa47",
+ "id": "965cc75c",
"metadata": {},
"source": [
"## A"
@@ -51,14 +51,14 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "ae476b31",
+ "id": "7eb30cc9",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "b32ced14",
+ "id": "154c1460",
"metadata": {
"heading_collapsed": true
},
@@ -69,12 +69,12 @@
"\n",
"Treat `Pclass` and `Sex` as factors.\n",
"\n",
- "Print an ANOVA table."
+ "Print an ANOVA table for different types of sum of squares.."
]
},
{
"cell_type": "markdown",
- "id": "24492fad",
+ "id": "3ad89dec",
"metadata": {},
"source": [
"## A"
@@ -83,30 +83,56 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "0ad3c464",
+ "id": "39f58924",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "ded672e8",
- "metadata": {
- "heading_collapsed": true
- },
+ "id": "ecf3bcd9",
+ "metadata": {},
"source": [
"## Q\n",
"\n",
- "Let us ignore the not-normal residuals and play with post-hoc tests instead.\n",
+ "Because we have a large sample, we will ignore the not-normal residuals and play with post-hoc tests instead.\n",
"\n",
- "Split the ANOVA for levels of `Pclass` and `Sex`, perform all pairwise comparisons if it make sense, and correct for multiple comparisons.\n",
+ "First proceed considering type-3 sums of squares."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "63677a0d",
+ "metadata": {},
+ "source": [
+ "## A"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "09794922",
+ "metadata": {},
+ "outputs": [],
+ "source": []
+ },
+ {
+ "cell_type": "markdown",
+ "id": "d30cdb9c",
+ "metadata": {},
+ "source": [
+ "## Q\n",
+ "\n",
+ "Let us suppose we want to use type-1 sums of squares instead.\n",
"\n",
- "We are not interested in the significance of the slope of `Age` for the different levels of the factors."
+ "Proceed again to performing with `Sex` first, `Pclass` second, and `Age` last.\n",
+ "\n",
+ "In the post-hoc comparisons, we will disregard the effect of the slop of `Age`."
]
},
{
"cell_type": "markdown",
- "id": "00ace9f5",
+ "id": "7a417d76",
"metadata": {},
"source": [
"## A"
@@ -115,14 +141,14 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "9f6645bb",
+ "id": "034acc94",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "9a863196",
+ "id": "b7b20012",
"metadata": {},
"source": [
"# Linear model with multiple variables"
@@ -130,21 +156,19 @@
},
{
"cell_type": "markdown",
- "id": "ef97cb7f",
+ "id": "144b0584",
"metadata": {
"heading_collapsed": true
},
"source": [
"## Q\n",
"\n",
- "Load the `mi.csv` file and plot the variables `Temperature`, `HeartRate` and `PhysicalActivity`.\n",
- "\n",
- "We will try to «explain» `Temperature` from `HeartRate` and `PhysicalActivity`."
+ "Load the `mi.csv` file and plot the variables `Temperature` vs `HeartRate` and `PhysicalActivity`."
]
},
{
"cell_type": "markdown",
- "id": "c66b1c3b",
+ "id": "35949307",
"metadata": {},
"source": [
"## A"
@@ -153,29 +177,25 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "4bdeea0c",
+ "id": "47cf88a3",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "358d3903",
- "metadata": {
- "heading_collapsed": true
- },
+ "id": "62449fb6",
+ "metadata": {},
"source": [
"## Q"
]
},
{
"cell_type": "markdown",
- "id": "f3003272",
- "metadata": {
- "hidden": true
- },
+ "id": "c71d2a4c",
+ "metadata": {},
"source": [
- "The `PhysicalActivity` variable exhibit a long-tail distribution. This is usually undesirable for an explanatory variable, because we cannot densely sample a large part of its domain of possible values, and therefore a model based on the data cannot be reliable.\n",
+ "The `PhysicalActivity` variable is very asymmetric. This is usually undesirable for an explanatory variable, because we cannot densely sample a large part of its domain of possible values, and therefore a model based on the data cannot be reliable.\n",
"\n",
"We will proceed to transforming `PhysicalActivity` using a simple natural logarithm. `log` is undefined at $0$ and tends to the infinite near $0$, which renders its straightforward application to `PhysicalActivity` inappropriate. Therefore we will also add $1$ to the `PhysicalActivity` measurements prior to applying `log`.\n",
"\n",
@@ -184,7 +204,7 @@
},
{
"cell_type": "markdown",
- "id": "34dff28d",
+ "id": "21e8879e",
"metadata": {},
"source": [
"## A"
@@ -193,14 +213,14 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "16e7787e",
+ "id": "0b8cac58",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "7887928e",
+ "id": "2c6d5225",
"metadata": {
"heading_collapsed": true
},
@@ -214,7 +234,7 @@
},
{
"cell_type": "markdown",
- "id": "b9ed6879",
+ "id": "c47489e6",
"metadata": {},
"source": [
"## A"
@@ -223,14 +243,14 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "1044c8c6",
+ "id": "bf32c7e6",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "49408adc",
+ "id": "1c7871a2",
"metadata": {
"heading_collapsed": true
},
@@ -244,10 +264,8 @@
},
{
"cell_type": "markdown",
- "id": "1b977a53",
- "metadata": {
- "heading_collapsed": true
- },
+ "id": "12d7e9a5",
+ "metadata": {},
"source": [
"## A (with nested Q&A)"
]
@@ -255,32 +273,20 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "b891fa3e",
- "metadata": {
- "hidden": true
- },
- "outputs": [],
- "source": []
- },
- {
- "cell_type": "code",
- "execution_count": null,
- "id": "c469b3d5",
- "metadata": {
- "hidden": true
- },
+ "id": "82e9f8d6",
+ "metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "475eb53d",
- "metadata": {
- "hidden": true
- },
+ "id": "66ddaa3e",
+ "metadata": {},
"source": [
"### Q\n",
"\n",
+ "The interaction term is not significant, but most importantly, the increase in log-likelihood is very small; the interaction term does not help to better fit the model to the data.\n",
+ "\n",
"To get a better intuition about the log-likelihood, plot it (with a dot plot) for different models, with one variable, with two variables, with and without interaction.\n",
"\n",
"Feel free to introduce one or two extra explanatory variables such as `BMI`."
@@ -288,10 +294,8 @@
},
{
"cell_type": "markdown",
- "id": "f7a51e03",
- "metadata": {
- "hidden": true
- },
+ "id": "05844c56",
+ "metadata": {},
"source": [
"### A"
]
@@ -299,16 +303,14 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "3a4295ad",
- "metadata": {
- "hidden": true
- },
+ "id": "2588d8d6",
+ "metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "2474a218",
+ "id": "665a7a5c",
"metadata": {},
"source": [
"# White test for homoscedasticity"
@@ -316,7 +318,7 @@
},
{
"cell_type": "markdown",
- "id": "1bf04c2f",
+ "id": "6d0bdd7b",
"metadata": {},
"source": [
"To keep things simple, let us use the `'Heart + PhysicalActivity'` or `'Heart + logPhysicalActivity'`.\n",
@@ -328,7 +330,7 @@
},
{
"cell_type": "markdown",
- "id": "77c0d822",
+ "id": "a6f37b2e",
"metadata": {},
"source": [
"## A"
@@ -337,14 +339,14 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "4d7f9d1f",
+ "id": "77e350be",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "20fdd05b",
+ "id": "789bd3f7",
"metadata": {
"heading_collapsed": true
},
@@ -354,12 +356,12 @@
"We will further inspect the residuals for heteroscedasticity, using the [White test](https://itfeature.com/heteroscedasticity/white-test-for-heteroskedasticity).\n",
"\n",
"`statsmodels` features an implementation of this test, but the [documentation](https://www.statsmodels.org/stable/generated/statsmodels.stats.diagnostic.het_white.html) is scarce on details.\n",
- "Try to apply the `het_white` function, but do not feel ashamed if you fail."
+ "Try to apply the `het_white` function (if you can 'x))."
]
},
{
"cell_type": "markdown",
- "id": "34e7a050",
+ "id": "0a822074",
"metadata": {},
"source": [
"## A"
@@ -368,14 +370,14 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "5de56d54",
+ "id": "85d73982",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "98ca812e",
+ "id": "e3ccb464",
"metadata": {
"heading_collapsed": true
},
@@ -393,7 +395,7 @@
},
{
"cell_type": "markdown",
- "id": "3ebf1176",
+ "id": "7ce19023",
"metadata": {},
"source": [
"## A"
@@ -402,15 +404,17 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "46c53e4e",
+ "id": "7dec1d91",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
- "id": "46bd2a33",
- "metadata": {},
+ "id": "b588b5c1",
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## Q\n",
"\n",
@@ -434,7 +438,7 @@
},
{
"cell_type": "markdown",
- "id": "374e25eb",
+ "id": "5224316d",
"metadata": {},
"source": [
"## A"
@@ -443,7 +447,7 @@
{
"cell_type": "code",
"execution_count": null,
- "id": "08216a83",
+ "id": "1bac9cde",
"metadata": {},
"outputs": [],
"source": []
diff --git a/notebooks/statsmodels_TP_solutions.ipynb b/notebooks/statsmodels_TP_solutions.ipynb
index 83257005c05c6324f1eac1dd6533be8e5c884aab..1e094a2a7705815459e26fbf5cb01c58f8dc8e14 100644
--- a/notebooks/statsmodels_TP_solutions.ipynb
+++ b/notebooks/statsmodels_TP_solutions.ipynb
@@ -43,16 +43,20 @@
{
"cell_type": "markdown",
"id": "965cc75c",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A"
]
},
{
"cell_type": "code",
- "execution_count": 2,
+ "execution_count": 56,
"id": "7eb30cc9",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -192,7 +196,7 @@
"4 0 373450 8.0500 NaN S "
]
},
- "execution_count": 2,
+ "execution_count": 56,
"metadata": {},
"output_type": "execute_result"
}
@@ -204,9 +208,11 @@
},
{
"cell_type": "code",
- "execution_count": 3,
+ "execution_count": 57,
"id": "17f37d48",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -214,7 +220,7 @@
"(891, 12)"
]
},
- "execution_count": 3,
+ "execution_count": 57,
"metadata": {},
"output_type": "execute_result"
}
@@ -225,9 +231,11 @@
},
{
"cell_type": "code",
- "execution_count": 4,
+ "execution_count": 58,
"id": "31c6e99b",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"df['LogFare'] = np.log(1+df['Fare'])"
@@ -235,13 +243,15 @@
},
{
"cell_type": "code",
- "execution_count": 5,
+ "execution_count": 59,
"id": "e5b46f36",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
- "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAABd30lEQVR4nO29e5gc1Xnn/32rqu/dc5/RCCRZGiSQwQaDBQZbq8jYm+DYIZssm0XZONms/YM8SQzOrr12srHXIZtd88SJjZ3LwmLn5g1soo3XxDE4xjLGJGAjBOI6IHkkoREazX363nU7vz9OVXV1TXV3dU/3dPfM+TzPSNM91VWnqqve8573SowxCAQCgWDzIHV6AAKBQCBYX4TgFwgEgk2GEPwCgUCwyRCCXyAQCDYZQvALBALBJkPp9ADcjIyMsJ07d3Z6GAKBQNAzPPPMM/OMsdFGPtNVgn/nzp04evRop4chEAgEPQMRnWn0M8LUIxAIBJsMIfgFAoFgkyEEv0AgEGwyhOAXCASCTYYQ/AKBQLDJ6KqoHoGgWb746Gu4/4lTyKkGEmEZ79k7ipm0irNLeWwfjOP2AxM4uHes08MUCLoCofELep4vPvoa7jlyEgXNgCIBuZKOrz13HpMzKxiIhTCbKeLTD72ExyZnOz1UgaArEIJf0PPc/8QpSAQokgSJJNiFxlcKOogI8bCCkEy49/Gpjo5TIOgWhOAX9Dw51YBE5dd2iwnT1WoiFpIxvZRf34EJBF2KEPyCnicRliuEPFmTgHsyKGgGtg3G13dgAkGXIgS/oOf58P5dMBmgmyZMZsKW9/0xBYwx5FUdmsFw+4GJjo5TIOgWRFSPoOe5472XAkA5qieiOFE900t5bBNRPQJBBdRNPXf37dvHRJE2gUAgCA4RPcMY29fIZ9pq6iGiASI6TESTRPQKEd3QzuMJBAKBoD7tNvXcA+ARxtgtRBQGILxrAoFA0GHaJviJqB/AAQD/HgAYYyoAtV3HEwgEAkEw2mnq2QVgDsCfEdGzRHQ/ESW8GxHRbUR0lIiOzs3NtXE4AoFAIADaK/gVANcA+FPG2NUAcgA+6d2IMXYfY2wfY2zf6GhD3cMEAoFA0ATtFPzTAKYZYz+wXh8GnwgEAoFA0EHaJvgZYzMAzhLRZdZb7wHwcruOJxAIBIJgtDuq5yMA/rcV0TMF4JfbfDyBQCAQ1KGtgp8x9hyAhhILBAKBQNBeRK0egUAg2GQIwS8QCASbDCH4BQKBYJMhBL9AIBBsMkRZ5h7isclZ3Pv4lNNA/IaJITw5tVi1obh3e7/SxI3uUyDoJEHuaUF9RFnmHuGxyVl8+qGXEJIJsZCMhVwJsxkVo8kwRpIRFDQDmsFw181X4ODesVXbe//ut8/5bAlzWRVjqTCGE6v3KRB0kiD39Gak68oyC1rHvY9PISTzxuFEhHRBh0RApujfUNy7vV/Dce82maIOMIYL6RJevZDBzEoRmmGIJuVdzmOTszh031PYf/cRHLrvKTw2OdvpIbWFIPe0IBhC8PcIZ5fyiIVk57VqmJCI/2/jbiju3d77d79tiroJg/Em5bJE0E2G+YyKE7OZdp2WYI3YWvBspoiBWAizmSI+/dBLG1L4B7mnBcEQgr9H2D4YR0EznNdhWYLJ+P827obi3u29f/fbxjb7SQQQCBIRQICqmxB0J5tJCw5yTwuCIQR/j3D7gQmsFDScmM1gciYNzTChGwypqH9D8dsPTEAz+PvVGo57t4Hl7pGIwBiDafI3wjKtGs9GoBUmkk6bWTaTFhzknhYEQwj+HoIAgHHNXCJCX0zBQCyElYKGsVS0wsl1cO8Y7rr5Coylor5/99smGVUwGFcQUSQYjEGRCcOJMPZs6evMCbeRVphIusHMspm04CD3tCAYIqqnQzQalnbovqcwmykiHi5H4OZVHWOpKB647fqWjanRqIlmQka7IQSvFdfTbx9zmSLyqoG+WGhdzlVEughEVE+P0IymuB5L+kY1qiDn0Q1asR+tuJ7efaQLGhZyKnKqvm7nKrRgQTOIBK4O4HbIAUA8rCCv6rj38amqD+z2wfgq7bIdS/qDe8cCC40g59HMua4Hrbie3n3MZ0sAgKgiO47W9TjXRr4zgQAQGn9HaEbb7EbHVjMho37bdIJWXE/vPkq6CTBgNBVxtumGcxUIvAjB3wGacch145K+mZBRv206QSuup3cf8bCMkVQYqWjI2aYbzlUg8CJMPR3g9gMT+PRDLyGv6hUOuXraZrct6YOcR7Pnuh604nq692H7M7rxXAUCN0Lj7wDdqL03QzMho716rkHYTOcq6G1EOKdAIBD0MM2EcwpTj6DjdGOcf7vYTOcq6F6EqUfQUbo1zr8dbKZzFXQ3QuMXdJROxvnX075brZ13a05DqxCrmd6hrYKfiE4DyAAwAOj17FCZoo5D9z0lukG1iWYezHY/zGeX8hiIhSreW4/Yd3epA7f2fRewqpGN39+boVPnuh6043oJ2sd6mHrezRh7WxDnwxvLBWcZfHohi3uOnMSp+axYFreAZswM62Ga6FScf71yxu0od9ytOQ2t4N7Hp6AZBmZWiqKJTw/QVTZ+IgTuMCVojGYE2XrUeu9URnK9jOJ2ZBx3Y/Z1qzgxm8F8RoVuMtHEpwdot42fAfhHImIA7mWM3VdrY4nKdd/rdZgSNMbZpTxkAqbmslANE2FZwkgyXPN6rodp4uDeMdwFPslML+WxrUlzUqNN4+vV6tk+GMfphSzSBd25XhFFgm4y7L/7SFNmr1adazei6iZA5WeYCDCJiSY+XUq7Bf9+xtg5IhoD8G0immSMPe7egIhuA3AbAEQGtzjvh2XJeeBsNsqyuBOkIgpOzGYhS+RoZOeWi9gzlqz6mW4sDOeH1758aj6LH55edJrG+9mb62UU3zAxhB+eXoREvCNZUTOQUw0MxpU12bDbkX3dDU7VkEwoaIBpMhABdnrQRm3i0+u01dTDGDtn/T8L4GsArvPZ5j7G2D7G2L5QvN9ZBvfFFJgMVTtMCRrDSdRjrh/3+z70imnCr2m8REC6UN1MWC/L9smpRYwmw06LSwY+Aag66yrTY7eEiF66pQ/DiTAUmTZ8E5+NQNs0fiJKAJAYYxnr9x8HcFetz1w0EMNYKorppTx2Didx6Fq+XN9oy+JOkFUNXDwQxXxWdVZS48kIcqpR9TO9YprwmqSCmglrad9nl/IYSUYwmooCACZn0l1peuyWEFF7BTXer4g6RT1AO009WwB8jbjNTwHw14yxR2p9IBVVVnU/uqNtw2uebgyLrIdttpkYLZt27I5Ttei2wnB+eE1SYVlyJrQXzq1AIqAvquDNW/vXtM9uND12S4horygJAk7bTD2MsSnG2FXWzxWMsd9r17HWk24Ni6xHr5htmsF7bl7zlcmA5YKO8b5w0/tMRbnpsS/WXabHbgoRPbh3DA/cdj2+/4kb8cBt1wuh38WIzN0GaWZp3cxnWr1CqKaRAahImuuElrbWc/WeW8lgkCyfImM8woQI+M7kXEPjuOWaix1T466RJH7+uu4wPbrHmYooWCloANBxE0unV7WC4IjqnA2y/+4jGIiFQK7QU8YYVgoavv+JG1vymccmZ/Hxw8eRKerQTROKJCEVVfDB69/U0kzmbmjU3Y4xXPJb34QiARKVF7QmM6GbwI/++0+u2zjagd840wUNw4kwcqrRsQmpV67fRkQ0W18HmllaN/qZux+ZxFJeAwOgyBIYgMWcii9992RLzUXrkaDlx2OTszh031PYf/cR3PHgs1B1o6VjSIRlmB59xmT8/Wp06lo0it84FZnwxkoRnVTheuX6CThC8DdIM7byRj8zNZ+z4scJBIJEBIMBmsFa+mB1oh+u19+RU3Us5FSkLXNFK8bw4f27YDJAN01L0zdhMv5+Nbq1N7AX7zgzRQ3zGRV51ehoOGevXD8BRwj+Bmmmy1KrOjN5U2HW+mB1wjHo1QyjChcW89lSy8Zwx3svxZ037kYsJEM3+XW688bduOO9l1b9TDc5SWvhHedcpgQQEFGkjmravXL9BBzh3G2CZkIcG/nMruE4Ts7lQJ4sSMUzTa/1wepEP1xv+OFoKoJzSwWUdBOMsZaN4Y73XlpT0Hvp5t7AbrzjLOoGJCKMJCPONp3QtHvl+gk4QvB3IZ9835vxscPHkS3pMKyiVylFRjQkt/TB6kTstV8NnFRUgW5yZ3ennJO9EofuHWcirCAeltHnmkz9FIJ2R9z0yvUTcERUT5diP6je0Ev3e+vRr6DVAuOLj76Ge46cdGrgmIz/1DPFtKNpynqEH7bjGO59JsMyFnIq+mKhqtE0IuJmY9NMVI8Q/F1CowJiPR7mdhzj0H1PYXJmBSsFHSbjwr8/pmDveP+qrG33OLwroGREweduuWpV05Sg4+zV6+e3z5WChtFkBNmS7qtpH7rvqVXF9uys7WrXXNA7iHDOHqWZzN71CJ9rxzFOzGaQLRoIyRKiIQkhWUK2aNSs2/7Zh1/Bcl4DMwGZCMwElvMaPvvwK844G20C0qnrpxkG7njwWey/+wgO3fdUw9E3fvvsj4UwEA9XzZgVETcCLz1v498I2YLNZPa2q0aL+3rOZUoY74tU/D3IMWp9J83UbT+1kOemIan8GWYynFrg4zgxm8FKXoPkKjk9n1GhGdUnE7/rpxsmjr2+1HS9/XrHsEMvGYAdQ/GmSjs38723orz2RnjOBGV6WvBvlD6fQR9m98OXLmjQDdOpHgmsPcrHmzFsMmB6qYDtREhF+fgWciXkSkZV4VjvO2lH3Xa/yURjJpbyWtVxep3MEgDNZGAMmFkpYD5TwscPH8fvW+akZvAKXCf0Ui6HXs5lirjjwWfRFwsFEqjNCPHbD0zgzv/zLNKFHBh4WHBfTMGn3n+5s00twd6q5yzI5CEmmPWhpwV/t5SkXStBHmbvw2eYJmYzKgBgJBnxjfL54qOv4f4nTiGnGkiEZXx4/66aDlQ7Y1iWCIosQTd4mYM3lgu4dIuChVwJsxkVo8mwIwA+dvg4RpMRZEo6tg/GsZQr1fxOLt3Sh1PzWWSK7qieEHaNlKuGeh/+sVQE51eKIFaeLEwG7BlNAFjdBEQ3+KQlg1UVVN5GKyWdz0D83AmMAUt5DXc/Mtn0vVQv9DJd0LCQU2EyVnMF4K3NM5cuQjVYRTmPWkJ8vC+MdEGvGFu6oOP56WXHR+Ke8L2TXiuesyCTh982Hz98HMOJMLKqISaCFtLTNv6NYrsMktnrtWOnCzr6owryquGbFPbFR1/D5x89gXSRO0TTRR2ff/QEvvjoa1XHMTWfA2MMmmGipHHhSeCa8EpBQ65kYDQZxmgqCiKCbjAs5zWcms85D+qJuSx0o9Js4/5Obj8wgbAiY7w/isu2pDDeH0VYkZ1ztR25z55dwoV0Ec+eXcJyXkUiIoPAhToBGIyH8Imb9gJY3QTEBCATEA3JVe333kYrNowxJ1taIn5NmsWbuJcIKxhOhJ3QSztpLapUH6fX/7NcUJFVDRimyWs/ESpKNfj5i/7fc+dB1vWIhXhYsCIT7n/iFAD/EiH2pAe05jkL4ofx+i8Mk2Epr+H0Yr6jWckbkZ4W/BslWzBIZq9fM+tMUUdYkXyden/82I9W1W5h1vvVMEwThqVNM5T/lwn4/iduRF8sVJEoNJ8tQSLAYK6uVJKEC5lSxX7d30m9c/Vz5BY0E6mwjKt3DGJrfwxX7xisMMF4JxPb5GOYDJMzaUxZk5FbUJ1dyiPizYhD2fTUKtylir9469UIK7IzwZd0E2A8ic3GK1C9wjBd0CFLhLAiY+94H/aMpdAfCzkC1M/562q45iARnJ4FfiVC3JNeK56zIM3YvRPMXMa6v8z17XrmriXVjAO+F+hpU89Gyhasl9nbqFO01OD7ABCWZejm6o5cYZk/jF6TlGpp3+7mJFv6IpheLtb8TmqdazVH7lxOwxO/6R966E0eCssSSrrBJ60q/YW9PYg1g4tGBlg1/fnEt3ukdUqEd5zxsIxERHb8J8BqgdpodzE/fxFh9YRWr2idm1Y8Z0Hu3yD313rVkup1v2E9elrjb1UNnF4gZDk/TZM3GjEt+0QzTtGqmky1XVnve01SskQwPRqrIku4dCxZ8Z3ccs3FuPfxqbZqUG7NesdQjJtBavQX9vYgVqh8qgZjIAkYiIfwyfe9uW3j/OKtVyMkyzVNfF5t2zZNVesE5qedD8b5RFCtaN2u4ThPpHPdWybj79tjXutzFuT+XXV/Eb+/3KvM9a4ltVGrjPa0xg/0RmvAVhDEKeomLBFUb21ii2qaTL1Vgldj3TkUx0JOhSxRRZ2dT71/b9MRIRMjCZyYzVZ15AYhSH9h7zaRkIzBsIR00cBYKrIuJQeClDm4/cAEPnb4OM4tF7jJA4BuMAwnFN/aRn7aeTIawo9dOoLvTM75Ovr9SoQMRConvbU+Z0HuX+/12DWSwFy2ZDnbW1fHqRbd0sqy3fS84N8sNNrM+tdv3I0/fPTEqve3pCKOJpNXddz9yKQTAWJUmSjc73sFgF9pCfffG40I+cRNe8sRJgaPWnE7coMQpL9wtW12j61vNmsQgUoAwGCtsiT0xWQMxEK+tY2aqZlzcO8YPnfLVW2tsxP0/m30/mo1rch56AWE4O8RGn2gbW3ODuc0TYbRZBhjfWXhpxsmTi/ksXM4joFYCNNLBd99yVJ1c1I9wdWoBnVw7xh+f41CKIhNulv8Q/Xi1u99fAp9sRDG+2POe3lVx2Aigkd+o/YE1Yifut0r52aLuNUbV6vj/rvlvmg3olbPJsGvXsuJCxmAgD1jKQDAi+dWfIVFKiLjhd+5qWXHXY86MUE0xfXWJv3GWK+WTzNtO5upD9SLiVPtqrfU6fuiUUSRth7G++C1uvKm30NyeiGHbQMx9MXCAICpuSwKqgHb0h+kgFozx7V7xNpJOetRZbRT1BKoQSZFv23ms0XkSoZvtm8zE22vVu8Uxec4zQh+YerpArwO0FPzWfzw9CLGUmEMJyItCSk7uHcMt0wvV2Tybk1FoLiiQ0aSEZxbLiAsEXaPJR0hvZQrOaUPgghpr7C75ZqL8eTUIqaX8khGeFy5ZrKq53rHA8egm0DJMANlHLeKVmu99RzbQcxgXtODX/a0O7u1mfpKvZoBv1kcse2g7eGcRCQT0bNE9I12H6tX8YaQZYo6JOJp9a0KKXtschaHj53DaCqCN4+nMJqKoGRwk4EdPqfIhIF4CLtGElgpaAhbWZxuIX3PkZM4vZCtmknplzl6+Ng53H5gAt//xI0YiIfRHwtVPddMQUO6ZCCvGVAk7li758jJmhnHraCZCqn1qBcaGCQxyhtK6c2e9ma3EoBzy8WKHsb1nJO9mgG/URI4O8F6aPx3AngFQN86HKsnaTRJx6YRDdVPqwN4PPhAPOzYMz/1/ssrzAaqYTrbuoX0SDLqW2SsXq0ev3MFY8irBiZn0k4iFQGQSIJEPP78f35vquZKo15donp/b1brrbXfehppNUfieF8YV37mW84+37N3FAB31pZ004nLB1Znt473RzG9VMCFTBGpqFLVOekt+JcraVB15oRa9sUU7Bz2DxXuFjaLI7YdtFXjJ6JtAN4P4P52HqfXaTRJB2hcQ62m1WVLupNQ5C374P2Md0Kyi4zlVD1wrR7vuUpE0HjlAkfor4Ix5DWj6rnaXb0KVVYJ9f5e6/rU0nrr7beeRuqXGPX2Hf146PkZZ5+5ko6vPXcekzMrvhq9/V3Y90oqGsLFA1EwhqrJVt57R5EIS3kdBc1wvt/ZjIobJoaqnns3sJkSOFtNuzX+LwD4zwBS1TYgotsA3AYAO3bsaPNwuhOv5pKKKpjLquiL+SfpAI1rqNsH4z4JNMqqBLBapZ/DsuR8FqhSZMyq1WM7jIFKYec9V61K0ph7CtBM7miudq73P3EKEgGKxMdlrxLuf+IU7njvpXX/bl+fRuO36+03iEbqDVe88jPfqtgnA584Vgo6tg2u1uhl4iUp3Nmtiizhmh2DVR2c3nunpJuQqVyfyb43npxaxB1Vz747WK8Ezl6MeqpF2wQ/EX0AwCxj7BkiOlhtO8bYfQDuA4DL3vI2dui+pzp+cde7v6uf4/XmK8cxk1arhpQ16tjyliFWDRNzWRU/f11Zq6tX+tk7IZV0XkvFXbJhS18EZ5cKOHEh41s22BvPzQiQYdfI4fuwhb7JTKdy5nCi+rnmVK4du3EXIav3d6A5s4HffhnjlVBtZ7jbsR0kNNC7T/ua2NeBa/QMM+kSVgpaU9mtfuY2ReblEfaO9znn0e02/vViI9bvaafG/y4ANxPRTwKIAugjoq8yxn6h2gfeWC5gm2c5v94Xt96X3MxNEGSftuN1hyV0nnl9peaytVEN9ckpHjljNx6x7bhurc6rCY4kuaafKxmWkEni568bcgSZX5Ex1TB5whfBt2wwUKmlXfbbDzsZqTYl3YDJAN3khcS2pBREPGYY97kmwrJjprBxFyGr93d7TI0mGHn3a/cvIKDCsd2I+cG7T7tshXvsXo2+0bhz773jXckBwknqplejnmrRNsHPGPtNAL8JAJbG/7FaQp9vV305v17U+5KbuQnasc9GNdSzS3kMJyKOMAdWa3V+q4jhRASKVJksZE8U9oTmHsNijsfouzuD1ToXv9o8EhEu25LEwx89UPU47nP98P5duOfISegm90GYVn0fuwhZvb/bNGo28O5Xt9TyMU9ZjEbuYe8+bXnfX8Ps1+i4vfdOX0zBbEZFKlr9GJuZjRg22lXVOSVXdiLQmYtbz8nXjBOwHfts1LEVJPSt0fA4vzEkI3KFvbneuXzipr0YjIeqNlkJcq53vPdS3HnjbsRCMnSTH+/OG3c79vt6f28W734BYMxTFqPRe9i7z0REwc+8bSv2jve3zIHpvZ47h5O488bd2DWSFE5SHzZi2GhXZe72bb+M/dh//rLz2puFtx4OlnrZgM1kC7Zjn41id7ZyV2BMRhR8ztXQJEiWbb1r3mzmaC+lyFdDZJJuTLo9s7mZzN2u0vgZQ9Xa5O1IsPGjXhvEIG0S12OfzeCu8gi2uvy+VxP0JnAFuebNnIu7Rr03pLRZOtFFab2+R8H6shHDRrtK498ycTlL3vo5RyO97k0DAElOaGE8LK+yHYckwmAiEngVEGTVUE8DbUZDrZdA9BsPHsNDz884537zleP4/K3XNHgFq9OOlUo1Oq3Bt0pDa6Z+UivOvd11mzrFRguJ7BZ6vkhbZOsetuM/3AOJeDKPyYCheAgXDcTwykwaEhEu6o85zarTBRXTywXsHE4EesA7tWSrd1w7EcgOtbSdj62wQ9s0WuWx2c90A60wuXi/M3eNnJFkpG33jve489kS5rKqU8uo28wMQel2c0kv0/OmHoCnnqs6c+KWlwsaiAhRhTs/7aQhALiQKSEkSYHbpPnVTlF1A3c8+GxbTQL1ara4E4Ekkqz/+fu1aMSc0YyDar2cWq02y7Si9oxfk3OJeNmKdrbkW4+6TZ1gs7Q07BW6TvCbrDLu254ARlMRgPEsQ7f9dEsDlQi9AsGv5EA7/Ab1BFFOrYwxB1YnGHlp1OfRDt9ErbEFFeTt8N20YsKqV64CaC7qrN61addx18paJ+deLQTX7TRrsQkcx09EMQA7GGOvNnWkwMfh/3vPJxUNYSRlOMlE2wbjTuKJm1oPuDdxxa/kgLcdYStskfWSrRJhGTlVB2MGGOPXgAhIhKt/PUFi/2uVRw7akq/RpKZGE9zakRzTiuJd7UhyCnJtvKU1DJNBt+6JyZl01VIb7aQVmau93tJwPf0TpslgMAbDtH6s5vS6yZy/6QaDaW3TDIEEPxH9FIDPAQgD2EVEbwNwF2Ps5qaOWgvPeUgEJ6kkJMv44q1Xrgo/rPWAu7+wVETBilXcKhaSfUsOeNsRtiKDuJ4ges/eUXztufPlS8D4dbCrMvpRL6nE72FtNIsUaDw56N7Hp6DqBhaylTWBqgnydiTHNNvmz007kpyCTHLe0hqa9TzwaqX+pTbaTSsmZ79nIF3QEJLIKW/Rrc7eVkx8hsmgmyZME1yoG9WF+3r4XYNq/J8BcB2AxwCAMfYcEe2q9YFmkCUCWc5NXpRLhkTk21QaqP+Ae7+wgmaAAIQkvk+/kgMX0rXLCjdDvXHOpFUMxBSki7pz7n1RBTNpteo+62lQnUozf+1CGumiDgnEC4gZDAs5FbqRbuo8qlFPA1tr8S7vd7ZzOIlD1w41tGLyEmSS85bW8JaZ8Cu10W5aMTl7r2cyoqCoGTi9mIdhMsxnS/jY4eMVeSXdQrVn6X9+70e4YfdwWXiblYLc/brbCCr4NcbYClVm1rb8bGQi7BpJwC95CACen15e9bB7B+PeJl3QIEuoqDOeiipOo2rfFYNpIhGS8NIbK44QHk6EoLqqSDaz7KsliM5agsQbPVPrwaq3iuhUmrldWlmyJBYRX7qqrpLLtVZhQTTpVhXN8n6P431hfGdyzrcOPgBcuW1gTVFWQSY5b2mNyZm0E+nVqQJq2wfjOL2QXVXnqdF6/e5n4KbPfw8FzYRMXEFgJrCc1/DZh1/pmOC3hbRtQrFNK6cXcuiLKtAMEzwFhkEi4PRCDueWCh0Z61oJKvhfIqKfByAT0R7wci3/3OrBXDQQw1gqGrhF38cOHwcB6IuFMBAL4fQC38YOuXtjqQATgCIRFGm19umniS/nSljI686YTAbMZTVs6+e23WaFTq3Johmt16+i54f371rTPoOM2xtT7n0NAGCA6aq9Y5oMBdXA/ruPIBmWsZBTne/MuwrbZu3z3sen8NtffxGpCDetuDOHm1nN+J3H4WPnnO9xcmYFT07pkAlQZHLq4A/GFVw8EF8Xkx/QvgJqjSor7u3BGC6kS5AlqqjXf+jaoaqfqXeMUwt5bs5yKQjMZDi10NoJrZqJxX7P/bdqJpYtqSgWcqUK53RRMzHeF2vpWBvlh1OLePDpswiN7nxro58NGtXzEQBXACgB+GsAKwA+2ujB6jGXKTnaLwCEZMLMShGvXshgNlMCwCrC2rIlHSsFzdnmQppvY4fc2ampXAgRJIlgMoalgu5EJzw/vQygrNVlSnx1QbCcrNb7iwU+GTQTllYvcqWZ6Bm/VoqHj52ruU93/9wgkRnecXtbL55e4K9PzZdbMZZ0E6moAkXibQHt6ydbAvb0Yh5Lec3pGBUPK+iLhTCYiOD7n7gRtx+YwOFj5zCbKUIm4MRsFifncpAJznU7MZtpKELE7/r/8WM/gqobzve4Yn2/BuMrRHuBstLCMMogGaDe7y0VVWAyOKWw/e6NehE3jUZOebefz6qAVTzPNjeNJsN4cmqx6WM0i2kyaIaJomYgV9KRLmpYzqtYyJYwmyni/EoB00t5vL6Qx6n5HM5YWvn5lQJm00Us5EpYzqvIFnXkVR2qbkI3zZp29Vuv3Q7d5H4dBv6/bjLceu32lp5bI/xwahH3HDmBhVwJYKZe/xOV1NX4iUgG8A+MsXcD+C/NDDIoumk6N82p+Rx3cpIrwscK57TRDL4k0wyjbHdiQNHahiyV07SWZ3wpB0hUfRVhmySY8w8nb2ms7Whm3YwzstF9eldQQTRY7zG8rRfdse2jKd6KcTAewlJew7bBGGIhGSdnsyAibHH1iJWIT/K2b8V9/dzHnJrL8hLPDJjPqpgYTToPa0EzAq9m/K6VbprWuPk2bjOsN5zYjqYZSYZbZmLJqzpeemMFHzt8HHvGUs737f3evKWw6/mx/L7XRldI3u0Nqx+zIhEmrAvmNTc1egx3VVaAgTGunV8ymsBKXqvQyt2aeicSTq+bGMKd2IMHnz6LmXQB430x3HrtdlzXhg5lJmPWpGYgW9KRK+nIq3ySy6k6siX++7demkFeNVaVXQlKXcHPGDOIyCSifsbYSpPHCYRE5GhXpkvYu9FdT6jp+p1cE4RhMkzOpMEYc5xjtv1OJiAaklclx9g9ZGsxEAthPlPCueUiAHIyiIM0s65nb2/UGdnoPr39c6s9mO7luneS88aU+8WYjyQj0A3TMdkxAINxBfPZEt5YKTirAM0wHIHqthe7z8up7Q8+mU/NZa0evajrF6h1HgAQkaUKJcKN+14CuIKhGwamlwrYM1a9Y1m9pj2piMKbpkiEdJEraQXVwOmFbIWw9rsXqjlygwjcRv093u3DsgTNMCu+Zz/fhP0ZWzhHFAmvL+aQKWqrnJ4f2r8L//2br1QUDeyLhfAf3rWLa7FdxnUTQ3UFvWGyVcI6W9KRs37PqzpyJVuIV76XtQR7vmS03nnqQ1AbfxbAC0T0bQA5+03GWNsDC+wZzX0x7HC6irnBJwzUAFaVfQB4yObkTBq6wSeGvEsIycSX+97jhmRqqJm1m3bEMDfqcAvy8Hu1R+8k57U3V7M/2/WUGHik1mJOgyJLkCXimpu1rT1pzKyUsJBVcclvfRMAkC1q2D6UQFiWoNsON5PHLhO4XVgzTLyxXIBqMCTCMq7YmsIdDz6LnGogokjQdAMGK2vyZxYLeNNQebLuj4ewmNMce7vtQHUG7sK+FwxWmVTnrXg6ny3hIw8cw7bBODIl3RH0/ZY/4+Rs1lmNukkXdIz3Vw95rTW5+H2vumHi2OtLTphk0mruEvT+896vI8kwzi4WwAC8cn4FskRIRBR89D17MJcpwTAZxlK8Xla2qEMzTIRkCcmIjIsHEpjLcEFu26TPpwvY2hfDT191EZ49u1KhRb86k8bv/sPLKGgGYiEZP/f2bfjgO3f6jtNvn81o4nY3ObdAtn/Pl3RkXe+Vte/V2/l9t63AntgTERkJ6/+puRx0kyEsS3ijiX0GFfx/Z/2sO36z3yszGSTCMiRY2hlgedvLmIwLtpDMe5KuFDREZBlFzQBkgkyADuYIeVni2xmM/85YuWyE7FpPeVvfBTHLBHHqNep8u2FiCE9NLTjnrBkG8qqBQ9f69y0OMvnc+/gUNKMcg2/3cz2zWKkZEky8cn4FkhWuCWbglfMrUCQJkZDEnXTzOeim6dSpNzwPBaEsaBkA1WCIKPx4ywUd6XMrIGsStrdXTR5VEZLtEgZ88sgWdTx5asnZd75KxvOZxTxkiRCSeUnqXzt4iWNCSUYUqLqJom763nNE3CE2mylro599+BUs5VQn21wzGIoATs5mcNl4H07OZlHSTSxk1ZpaXE418PpiHtNLBccm7l4lTC/loVkToHdySRc0ZIsaNINHrkngq+KQLDmmn3RBc44fC8mOn+CX37kTuZIOgzH86XdP4q+eOoO8aiAsS1AkYCgZQUTh2zPw54Cx8k+uZCBT5Cuvqy7ux3Nnl3nzGAI0w8RCzsQH3toPoGyTViRCX1TBQq6ER16+gDtv3OMI6r/659P4i6fOQCJAlngntr946gwAOMLfLegTYQWLuRISEQWxkITz6QJ+/x9fxU9duRUXD8aRU/UK4cx/N5z33UK8XSGX0ZCERERB0hbcEaVCiPOQcsX6ka3tyu8nIwrC3v6enuvZDIEEP2PsL5rae4OYlg3Pbj3n911IBLx5PIWCZqCgGTAMfpMTlW37EZlwqSv0zS4q9r4vPI4Ts1nLc4vVKwbrDQkMey/iN6xtXnBrtfWaWXtpNN8giP39b55+fZUwYdb7fiGHQSafE7MZrOQ1SBJBlgiq5q/BmAxO6ClD+XdY9n4wIKRIUGQJulndnLJ3vA8vvbFifxQSSZCJwYAJEzzSRwafFBgAYvw9O2SUMSCsSCi6juE103ixzX2EyvBMu1BeSOaRKyXdCksFnLaPhudcTs7lHMWBUL6fNJNfk1KVScT3elg/3ki1V2fS0MzKyLSiZjqTS76kYzGv8cgoGbBiEzAUVWAwrhHGwzJCsoRUNITzK2XteveWJC6ki6sErm6aKOlAQjOg6ibyqonhRAhDibK5rKAZePDps47QfvbsCobiIeRUw9H4E2EZz55dwQcBPPj0WSgSOU55ewL66lNncPFgDDlVxwNPn3WuhqMUMIa/+sHrOLWQxxvLBZxayDnf7wWTT8LLhUrf5p8/eSbgVa+ORLwJjq1h28I6ab0Xd2nfCdd7tuC2PyM3KZjr4fY7gKSGOykGzdzdA+B/ALgcvH8uAIAx1tJC44okOVr02cV8Rey3jS104mEFI8kw5rMqSILjNAQDxvvLYVZurTZT0nHxQBTzWRWq1e3JfmDth6QvFsJivmwCSEUVzGZK0Cwt19s8PCi1bPjNhCeeS/Ob3p1awVj5fb/j13Mgq7rJY5gNtmoFxTOdDatyKMMVW/sxNZcFA3/IJ0a5iemFc1yQ2zHPbtxC2WTlbd3n4Z4o7Lj1F63topbQ0MyyRk8NurciisS1zlio4vo+ObWIVERGuqg72bKAvZrk18NkwO6R8grJ0RL97JGel/UmJICbEWfSRYABW6wuXvbca5jMiZSy39eMspJkjw/gmnm2pGMwHgbAJ8dMUcd9v+hfwPFvnpm2hL7kOh8T2ZKBv//Ifhz6X0+hL1oWFYwxhGXCueU8Ti/kkC8ZOL2YQ1SRkIpKMF1ZqJMX0vjU11/Ey+dXQETW3ywlD8D0chEf/MoPK8Zjei6UbjI89tpcrYu3isF4qEJwV2jRFYK70oRifyYakuDJW+o6bL/D1z9y+oVGPxt0pvgzAP8VwOcBvBvAL6MNBd4mRhNOud9LfuubUCSraBtb9UwB4P1gs1boZs4wEA3JThy2X1q9be6whZRbm7ffy6s6hhMRDCYimF7KYzAeRq6kQzMZmAnf5uFrJajzzW0OqiZEagmXeg5kk5nwmWtX7dsWMH7O3fK+GhsbebZxK0qSZfKxcwPcFLXqhex8j2ON13t9X7uQRk41EJL46lG30+fBlQJZIgxEQvjk+97sfEaRrdUHc51Ak4z1hRENydCtCVP1mMYYKhal/D3GQxtliV+3iZEkzi7loRsmNNd34o45tx2QjqNR1R3TmHvSJXAT1G997QVkizoWc6o1ubCK7/Y//PnRuuf2TycXXGcRDPfllCTC+94yjiOTs4iFJEhW0tdiXnUmwu2DcUgSQdUNjCSj+MN/e1XgY21Gggr+GGPsO0REjLEzAD5DRM8A+HS7BpawHFIRy8Ria5tugTCfLUE1GLYNRrHDU//Dz/4etP7Kp95/edVoGKB28/BmCGJ/95qDzi0VHN8GuewM8XBlfHsjSCRBIv7w+0225Wbo/LWfc9dt8vDi/Zv7NRfsZdPIcCLsbBeSJciMa7z2Sq2KH7aOVm2Zh2Rp1fX1ZhyHZAJ3QxPG+6K+K6SJ4ThOzOZgegaiEKDqRuU4a4wrGZIQDykwTcaP79o2LJOz8jVZWTArEpBTde5nMRlkIizkVOs132Z6Ke/4rZYLGn7ynu8HdkDaQ3jKFatfC3vYElGF3fnyrX1403AcKwUNR88sQZEIEUWyatIAH3zHm/COS4aQiCj4u6Nn8Zc/eN3xE9irmF98xw588J07Mb1YqEikkq08H1kmhBRCUeOKSyfj63uFoIK/REQSgBNE9OsAzgFoa3nAD+/fhXuOnIRums6NAMZr2NhCeimvIRGWMbNSrIhssUsyeGmm/sp6lD4IYn/3moPGUhFcsByNbi35V9bQ5i+sSJBVnuhGBKgat7UDXOh4v4NUVMFcVnWSiwqaUVPyhxVybOeKRAjZk7pmwASgm3zikonbuO19pqI8B6E/FkIsJOO1C5kKM2CtycaLaTL0JUKrrm9YkZAv6TBgOkKHGNAXV/D1X9/vxJSfWy7wCokmwy+/awJ3f2sSOVV3hHaICFv67AStMGYzaoUg866oFAJkRcJctsS1dGvCtR3q7vBl9+JGN2FFXNnfD89Kd5N3+Wgyxeo5PtWuXywk4f1XbkUirGAhq+L49DIyRQ3DiQhuumIc104MWTZvGWFFwtOnlmrGutuO2Zl0AdsG4qv+/ovv2gUiwt88M+0b1XPrtdtxz5ETKGgGopbm3xcLYTAWQqaotzW+fqMRVPDfCSAOHk78u+Dmnl9q16AA4I73XopT81k89PwMNIMvtfduieONlZIT1SNLViy+5YzUTYb5jArNyFTdr9fc8djkLJ6cWqwqONajnGwQ+7t3Ahrri4IxhrmsCkki33aOjbJnLFURIhoLywAYijqDbvIy0e/ZO4qZtOokF90wwevbzGb4dxK2HLp25IdddA8oO0yBcqKKXYK6L6Lg+c/8BIDV7Qttf4r9HhEhpjDYPj0GICQBOuMaJwBEFQJjPOHPjrqRCEjFFFzUH8cvXL8Dl46nML2Uh2mC10oCQ7ZUdk72xxRs7Y/j9Hy2Ihqk/LuBfW8awvGzy8iUNEQUGVv6IlBkCSZ41EuFOcgHncHJGq4k2FQmEyEWLoe9xkIyJkYS2D4cRyLMo0Lq2bh/4cs/gKabWMprjtVqMB5CWJHwqwd3BxoHUD/WPUgs/AffubNq+KZfItWvHdy9roKeiFxZ/eSYHvl9vPpvBACe197tQPD9G/+o5xju/cIVWNEgNQU/Ef0jY+zHGWNPE9FvMsb+B7h9vy1kijoO3fdURbLLzuG404LuxFweY6kwdgzxRhtT8zlIABSSrPECJrFV9tFqVIumuWV62ak/kwzLSDdYRGwt1JqA/OL2r58YDhxdVA975THer1Rtj2cLZQZgOa/ixGwGo6mIY2qbmstBkuDYyt225ohL4y9HbfD/veV5//eH38EdzZZj9XuTs04CEDNNFDxfsa3cEnEhm1MZRhIKtvRFUNBMaAbDT191EbYNxZArGXj1QgbPvr6MrCXMNcPESoGn5NsROQXNxFx2BT/7p08Gun4FzcSyda/UI6zwqBd3REg8XI4acYf+JV1ORztaJBlpnQNya18MC7kSRpKVUTvDiUiNT60fbmH4zt0jeNeeEev96sKWqIqg9XlfourC1fvZjUI9jd9dEP7fgEf2tI03lgvYZtX6ODmbhW4yJMIKKFwly5YBJvjy3R01EZaDfUF+0TRzVh2XbYMxp4gYA7cLB43bBxqLyw8Szumt016tUNZaaDTs1PsdxcMKIiEJhmFCkQklK4kKsMMVCQTma+fPqjqGkxGcW87jk3/3PG77F5dg79YUciUdT59ewt8+c9bJt1BrzOvuOX8+p2M+V9amv/LPp1tynWwk4iuF8b6IFe7nHyESj6wW6CG55bERTVNpQuHRW4bJ8O/escMJla6nzUpUKUirabFuDRY++7OF60YUtt1EPcHfdAALEUUBPA4gYh3nMGPsv9b6jGqY+NGckxiMsEyYz5bQFwtBNeyWi+Us25DEH3RFJkcLDskSVIMFau7gZ7/PFHXoZmVpAwAYiIfx8EcPACgXxaqVpt/qLlRPTvGqo3ZnJrvEdCvqsjNml6MFjr2+iBfPLSOnGljOq3j69AKuftMgGGP4o++eRFHVMV/iY9BNO6mp6IRaDsZ4H4HhqIIhOYSzS9wOLRGvE+S+oZIRmUeUMG4G+tF8zsrlAO76h5fXeFarkS2TWLUwv4Q7mcaTKXnXQy9hernoVJS0I1wUAj77r69sekxVtU6yhKlHyEpUDmD1E5L1hK37tb2vXSMJjPdH19S4RtBb1BP8E0T0EPi9ZP/uUKcDVwnAjYyxLBGFADxBRA8zxp4KOjieuMNVOAlAyeQDsRt8GAyQJGC8P+qYg+ayKsZS4UAC189+X9JNJ5LIpl5nq7UWxTq7lEe/5TAF+GwbUSScXcyhoHLzw5nFHAbj4YokGsYYzizkrFA75spgZo6T0BZQDOV4dHsbHiprm14Y/uz7p/DXT591nJE51cCfPDaFY68v463bBvDCuWWoeqUuYIKbOE7OZir6Jc948gkMezAusqXKUMxa2ZMELnAlKkfgAFiV6GdHlDDGAAK+/EvXImkJ94gS3DTiFbwzmRJkCZAkyRHCMBmmlwvoi4UgubeX4AhzW6DbArdCM25Tck8zrLVxjaC3qCf4f9r1++ca2THjEiVrvQxZPw2vIOzyCfYHFTvtEvwB2tofdQqC5VUDY6mw08SinsD1i6bhxaLKQhjg4ZsX9cdQ1Az8yWM/giwBEUWGyfj/hqnjj757Em/d1g8GVDRucJ/H6YUszlrRGm6BPJKI4NxyrsKxaNc5Ob/CGz2MJaM4NZ9Fuqg7Drg+q/fqfKZY6Xz0qydS5f28UyyqXPvI8ASMPzm1WFGC149a8f9+hCRCPCIjU9TLIanggjsWlrAlFcOnf+pyJCIKPvP1l7CYV3kYHwELWdWJYFEkwK7QwE1B5eiZiZEk3rZ9wHGg2UJXorITWKqmJXsmCHJt736PiCps4wJBL1BT8DPGvreWnVslnZ8BsBvAHzPGfuCzzW0AbgMAuW91j9ntAzEs5flDPpIIIa8ajqljtC8Cw2T44393DRhj+MCXnkC/R+AqEuHMQg4zK0UwlDNSGWO4ZCyJX3/3JXjgh2d5Knt/DD9/7QgeefkCVgoaoiEJRc2EbjL8zNUX443lAs4scqGuu44RkgnTS3mnYqRf44aCZmBLKlYxNsYYirqJy7YkcHya1zkB8VVHUTOxc5jw5/90GjlVx4VMESuukDwGYKWo4/lzy/jxL3y/iW+nMSJKuZqlO+rAUqyRiikYTkTwE5dvwb5dQ040SSwk46tPncHfHp1GXjMQD8vYf8kw5rIazq/weitvrPDa+1x7Z8gUDfzi9Vuwb+cQiICPvGc3fufvX4bBGGKKjIF4CLppQtVNqAZDPCyBmSYYCLppQpEkRMM8Q/cDX3qiJf1c3WWE3THme0YTa7+4AsE6Q0HqWxPRC1itra8AOArgvzHGFlZ/quLzAwC+BuAjjLEXq20X2bqHbf2lL5Q/B+A7/+nHAAD/8f8c9xWmIVlCXzSE8+kCskUd8bCEwfjq6IRGMvnc8cbe2GB7HFFFcsoX5FUDqWgIv/yunbzO+rk0vjM562iZdt+Ai/pjkCRUaOftascZC8kuB6Pbrq2Uo0c8RaN++2sv+taXIQCXbknh9cUc/zvjmrXtZ9g+lMBXfvlaJ6PSth/bHZtqmVcO3fcUTs1nV/kudo0kK6KVvCGe1cofTy/lkXB1+aoWndQoj03O4uOHjzs+ILt0x+93YY9YweaCiJ5hjPnX46hC0Dj+h8GrHP+19fpW8Lj+GQB/DuCnan2YMbZMRN8FcBOAqoLfi93zFOCRB3d/axIX0kWnfndY4c5czTDRF1VgGCYWclzrHoiHHW3dm8mn6qZPk4NyfexXZzI4t1xAtqShqJn4yj+dwl8+dRq5koGlvMrNE56xXsiU8F/+X+1TO7WQq/l3G7cp4vKL+pAIK3hyasGyIQOwTBJ2BdEv/Nu3VRSNsotDuSewiCLjA1duxTt3j3Dbs8R9Jfx3LqC/N3kBf//ChVXjSUV4gs6Wvihet1rjERFMcJPKr797N/pcDettIfzahTQ0gyGsSBXNRvxq5Y+OlvMB0wW1oqywX29lL97eA5rJWtpo/uDeMfz+LVc17ABtpj9zN9Cr4xYEI6jgfy9j7BrX6xeI6Bhj7Boi+gW/DxDRKHiT9mUiigH4lwDuDjqwgaiM9755vOI9buu3nJFgyJcMJCM8VpybIXg6+FJeR6ZkICLLGE6G8eV/OoUvffek08VGa8Ag3UhsdkSRPIkyMuKeCJLKyBH++5e+cxLpooZ4WHZszd6Vyk996QmUdKNcSAvcGZoIy3jXnhHIRPjnk/P4wqNnML3M8w8WLa13JBlBuqjhi0dOYiQZqfoAP3t2pco14EZ02xbPbKcwW10izXZ+q7rh22zklunlil633pr/6YKGs4sFgIAL6SLmsyXc8cAxhEOyU9e+ntP+7FIeMlXWYmqme9YXH31tVU/jRnImWtUUfr3p1XELghNU8MtEdB1j7IcAQETXArBtLtVywbcC+AvLzi8B+BvG2DdqHUQiHm5HBBR1hk899CKGEmEwxuvyVJhFLIdeumQgXVpdqEs3efnalWJ9oU2AozGnCxp0kzlhhZJl0uiPhfFv9m3zmE4qhbjSYGw2N4UQ/r8Du/AH334NmskQs/wKDMCvHrwEO4bikCXC7QcmcM+RkzCtrmI8goZw27+YwEgygscmZ3H3t15FSCYMxcM4OZeFbjAkoyHY1UzrddyatkoAOJYZyx+imVzQX8gUQUR401DMaZvo3acd0bSQ1SGBryZ008SFdAlEwD1HTmI0GUZ/jDvgvY1t3ljOwwSgWGYjZgJpzUBIM7DVqrpaT4NPhmWrTy85EWDnlovY3YA93i7TLBF3IBc0A/ccOQkAgbOjm6m6CnRe22523ILeIajg/zCArxBRElxOpgF8iIgSqJLUxRh7HsDVjQzGZKyiwxEAq4F6bWTLdGGbRyKKhOsnhivqYiddKeruJgeJiIxoSHY++7N/8k/IlnRuq7YiRDTDREnT8TNXX1x1DES8OJUkkTMeRZKc320h/8SJOdz//VM4MZupMIPcum97Rc2gGyaG8Gf/fBp3/cPLzsN/85XjFSUsrnvTAJ6cWsTf3H0E6YKGRER2BGq93rbAas1ueolHEPm5fVYKGhjj3czmMiWcWy5Ytd4lTC8VHLPMaxfSSEYUJ0qIjLJ5JiLzzN2FnIqIIqMvFlrV2EZnXNDaCU5EAIzKxCxgdZcpt3B09weAaxJrJBno/idOWUJfstoG8tIP9xw5iSu3DQQSgEG6Y/n5Kty+hPlMCR8/fHxdfQnrUZ9K0FmCNmJ5GsBbiajfeu22CfxNKwdk1xYHeIIPAPzKwd3426NnkSvZbfK4ME0XVKwUDYz3RSsicNxdfRrFW6WRiAtRzWToi4W4Filz4S5bNYJky05ej8cmZ3HXN17xNYOcWy44zke/pbbdoMMuYbGQK+GHZ5YxmgxjJBnBzEoRBdVARJGRioYC9Un1anYEfxt6RKaKRjb2eZd0Azmr4bNtllF1E0u58irL3h8BACNIxN+xE/OAysY2l/32w74Ntd3vpAsazi0XoUj+pghv34WwLGG8L4JsqXqhMi851YBi9XlwR2IZJsOH/vIoLh1L4hM37a0pjL15IvXGDQB3PzKJpbwGWSIosgTGgKW8hrsfmWyr4HevMtIFDbphOi00gdbXpxJ0lkC2CSLqJ6I/BPAdAN8hoj+wJ4FWQuCaXthqllEyTFwymsL73jKOXzlwCX8QAIQUXvgqpCj4sT0jWMyp+NEcT2S66fItVYU+WcI6JEuIWZEtfbEQhhJhjKYi2NofQzQkOVKmrDgSoiEZI8kIBhNh9EV5k4doiHc2CpqIYwvaTJGbQRSJfzZd0BGSCfc+PlWxXTysOGaabElHpqg776ULvIRFxupHEFF41pDd33QkGYHJ+GqIZzzrq2oM2bXbp+ayvDl9lXGPpniUlCOQrTh/WwtngGOWMa0yGt59MXC/zEgibGXqmr7jmhhJONFSPMmMd8wKycRbAFomJ4Cbiezr475+2wfjUKweC3vH+zAxmoQiSw0JrkSY52nYNerd5yMTcGo+h08/9JLTKtGP2w9MQDNYxbhNy1756oUMryqrG864AfD6U05WLTmJYVPzwQIDqmFnm++/+wgO3fdUxbhtRWPWKpeSiMiYy6qYyxSr3jtrPaagswQ19XwFPBrn56zXHwRvzvKzrR5QudkGg26Ua2uvqszXH8M12wfwyEszGEnyJhZF3cQ/vjKLfTuHcOCyUce80ohWfumWvtXF0BKhqk3MG8FeQtv9bAE+uXgbg/gttQ2r2fhLb6w4vg4J5UbxknXtiroBxhgUmTAQD2E0GalaYygZljF5IYt6TC8XsfOT/wAAiCk8e9XbfCVInfeL+mNWhzMVugm8+EYaALCtP4J7H5/Cb3/9RaQiCsIyoWD5OXgGLBBVJMxlSgjLvFyCBIbTC2XTQ1whpzjf7Qcm8Kt/fayi9248LK/qmua1pY/38UqjXNsn3yAA25xnmAwhmfCxv30OSwXdiTS7+cpxfP5WHgfhrX2kG6ZjPrR9Dws5FbqRrnndTKukRi3zUC2fgF9T+I8dPo7PWeYj78rPToDMlYyG6lN5r61wEHcvQQX/JYyxf+16/TtE9FzLByNL1gPHBeOukQTef9VWR4D/7Nu34d9cu92JDT9031OIhWXnhg1bjaG/+oPX8f6rLmpqDNUqVLaiGqe99A/LEnSjnAjkbQziV0rCaYzCyiYZW9TKVu0YZjBEw7LzsLobyvjRjM22oAN9Ed4W0d06MQipqILXZlZW2eunV0qYy6nYM5bCQq6Eom7y7msgmFbDjqFEGIosQTN4B6mSZx95nUHO89XO15+bXtVwPa8a+Ppz01ULzr1yfgVPTumQJS7cq/UKloic72w+XcSKK7DAMBm+9tx5AMcqhL99zCs/8y3eKtFlRjRNVtFbYNdwHCfnciCr8CDvBMbLTvsJUD8B+/HDxzGcCCOrGtg+GMfZxRyW85rj7GYmsJzX8NmHX8HBvWO+isZwIgJF0pyOeI0iHMTdTVDBXyCi/YyxJwCAiN4FoNDqwUgEXDKadITtB67citv+8pmqmkwQ59kNE0NOieUgERJBauM3iz2ppKIKFnIqX/YzrGoM4ldKwm7I4WuOsUwvkkTYMRjDI7/xY4HGk61V5tIDWbMNA4+kCpL452WloMEnAAsAL9Jmm7BkiZyuXrp1jWZWilBknj2sm/7HzlizwdePn/f9+9ePn8fnb+W/ewWT7XPhHcYkMMYHajt4NavHgGaaCEkSRlMRZ8Xh9hkzBjz0/IxzHDdhRUJBNZxVLbO+t7BStrh+8n1vrtDOTcYgAbh4MO4bneU9D8NkWMpryJR07B5NYjZTxPRyEYpVZ8geLzMZTlnjb7bnRK2VhnAQdzdB4w9/BcAfE9FpIjoN4I8A3N7qwdjN1sdSUdxyzcU4fOycY3e0NR23nXD7YNyJMQfKzjMC145OL2Rxz5GTODWfrboPPw7uHcMDt12P73/iRjxw2/Ut01AO7h3DXTdfgV0jSfRH+YqiP87NSO6s0oN7x3DLNRdjLsObzsxlSiDwaBcnicuFYZl2Lh6IroqKaimuA68EzG1wE0R7dPfxVQ3eDctkJkoGN3UFiZitlg3tfv/sUr4iC9zbH8DdX/iigSjCTvE3YGt/BHINs2G1YnN7xlIYSYWdpumKRBhJhbFnLOVsc3DvGD53y1W4evsgxvuiCMkStrvCZwGsMgu6z2MuU+I9ik3mTBTA6qgoN15fRBCbvtcv4H22vM8mIBzE3UTQqJ7jAK4ioj7rdZqIPgrg+VYOxt1s/dB9T9VdKno1Y6/Tz+0AHU1Fu2K5GaQK4mOTszh87FxFg5NMSfd9eCXiZheAx9SPuSIx6tl+vZUt3TgaqQ+yxKN8Lvmtb9asqNkM3j6+tqnDrmZp2vauNeLVcu1r4S5fbPcXTkVDSI2HMJcpIq/y3s9jqShOzed8r1+1SSGoGdGbhTybKVb8vZZZ0O5J7O6DbNdZcq803HWGmlnl1jPlBGknKugcQU09ALjAd738jwC+0MrBuDtw2an8brxLxYN7x3DL9LKTXWkyhpFEyNGO3NpjtX20giAJN40k5fg9VGEJvmaSkISKRvH2gxXEufaOnYN48tSS7xi8Qt/9+uYreUZ1PCQhU81242E4Hqq7DWMMfTEFsxkVqaiCiCLh3HLRaa1oJ9WFJfg2Y7lh16B1LAUL+dWhm8Px8u3uFUx9UQXLBd1ZYdimLXeP57Ai47M/e6Vz/X7jwWP42nPnV10r+/p4aUbA1hOg3r/LxFuQuiuGDsZDWMipIHBTqCJJGIyH8Imb9laMrRFlqJ4pp50mU8HaaUjwewieDRMQdweu+SxP5SciR5B7l4pezfjkbBZLeR3xsIa+WGiV9ui3j7USRMA2GuHg91BpPoKOAOiMfCMvAjnXSEI8JFU05ZZQ6TTemgrjfEb1jVrZNhjH5EzGt6ib+72YQviDn3sbAC5I0z6Nv2WJn8fO4SQOXTvkJLPtGUvizEIeJcNAWCKMpiJIRUN4bWalYiK8YdcgHrj9nQCAP/i5q/ErXz2Koqt3QFQh/MHPlfMJvYLpzVv7K6J6vP2F/QQXvw7H8NDzM77Xx49GBWw9Aer9+66RBOayJSgyVUxYH3n37ooEwbUK4SB+AVHjv3tZi+BveV1JonLHqy1Wjf0zC3kr0oJXQ3SH5HmFmzf93609+mnFrSCIgA2yjTeBJlfSoOrMmbhMK5on6rLnFjXeIu/sUgFvrBQx3heucK5569X4ZdmOpiIViU6JsIyVoo7RVKRqaKC9KlvIltyJsU7o5UBMgWbCqXFzUX8EH/rLozXNQhf1R31j/wfiYbzvLeNOfZ9YiEdupWJhTCQjyJR0x4nv7or2vreMu4Q4r7PjFUKtEEyfv/UaX0duK6k3zuenl/HSGyvOuVabsNbaqc3N7Qcm8LHDx3FuueBMesmIsipktlHWWq6iFeUuOl0yYz2o12w9A38BTwBirR6Mu8kF4ASr8FR7Wj0Qr2bsTf/3ao/tWG4GiV6ot413RZAv6VjMa5CsxCW3qaqkG+XQThd2KOHfP/9Np9+AY7NmgGaUs2xnVgqYz5SgmTzLNqRIvJetYSKXNRCusjLxjvPcUsER+M44GLBU0LFrJIEdIRmnF3KYvFA7+YjAE6WGExGcXsjih6cXnYzk2UwRh4+dwy3XXOx8jwmrmJ1qmBiIhXBqnn9mLBV27UPFaDKMHUPcyXj42LnApRZ6Cb+aQg89P4M7b9yNO94bvKBcM1jpNlUL9jXKWmP/W5E7sFnyD+o1YknV+nurMV3G0vlsCbJEiFkZmMDqgmB+y013+r9NKzUdL0GWvPW2uffxKWiGgYUsTxozTOaEW5lWzDgxA6rp73B1O2K9oY6rbPWAUwrAmU/sz1rx5JrB8OqFjFMbv1rooH1cIiCq8JVIUePOz2qOTzchmUe3EIB0QcdIMlrVIf/NF85jMBEBA/DGShFxV/5Gpsg/U2sfc5ki7njwWfTFQhtKi3PXFAL4iks3Tdz/xKnAxeSa4d7Hp9AXC2G8v6z/rTVwwm9lPJ8N/r21Indgs+QfNFZOss0wBiekrKSbACuXCwBWa9LNhKG1miBjqLfNidkMZtMl5FVuirLLHhCRU3IgpMhVRlAfryZW0syK+jOaaaLoasBil2CwM0tPXOA+fW/ooF30jIfaM6eVJFDuzFULmXiUjmGVcAD8HfK6YeLEXNYJHcypOhZyKtJWSKn3M97X6YKGhZyKnKo3FNbbCyUHcqoBbxCRRGhvWC9W3wvA2gMnvPvMFDXMZ1TkVSPQ99aKMbXjvLqRrhL8Fw3EMJaKYqXAa9OPpMIV8ct+zqO7br7C+cxYKrqmLkvNEGQM9bbJq4bTs9Zt7XKvgAqWmSYiS4gqjX1tUc+NzFAZxhmSJV6jyEVJN6GZZkVmqTc2O6pIvAG5FTeuuCWQuzJmFez+t0D5XG1/htshfyFTQkiSnDpF9upiPlvy/Yz3tb1dVJF9a/v4US9OvZ00MuHYNYXcmIy/307aEafv3edcpgQQD0cN8r21YkybJf9gLc7dlpOKKo6Jxn7w6sUBd0PkQJAx1NpG1a1GJ84/HMZQURZBlspVQ2UqNzivl0TrvZFXUeXzjJUjfAC+cvn44eM4t1SAbvJ4ccMEr58Dnkjm/mw9CprhzA0mAyZn0lbIIcNwotIhnwzLTp0iIoAYUGK80FsqqmAuq6IvpqwKCbVXj4TVq8cTs5kKh3DDUVFovSOwURvzh/fvwj1HTkI3TScXwWT8/XbSjjh97z6LugHJ08y+lvbdijFtlvyDrtL4J2cyjobTDdr8eiFLktNsvJaSbJhAUTe4Hb2Fx1dkqhDaNkRYlSnLAIDKZh4nqsfjca53LhX7c+1DIkLMqjdkZy3HQxJWinrZfGRNSHYI6K6RJO68cTd2Dicdp/6dN+7GrpFk1dXjQq6ETFGvqtEHWfLbtfOffX0JMysFPPv6Ej5++PiaVgV+lVlrabl3vPdS3HnjbsRCMnSTj5E7dttn3wfas9r27jMRVjCcCDvlu4Ha2ncrxrRZ5E5XafyKRKs0nI12wW3cmqLk0t692KYfbymBVgbT2s5zd+VPZh/L8rvY41WsSqcGGAwrszXkcsCfW8pj0aonLxGg6izQUCUCLh6Igwg4u5gHs95TddOp/OndT0k38cZKkZePYCZeOp9BTuUVJW+YGKq5elzMaRiMh6pq9EGc9n618xeyKn71r49hKBFeVfHTHWrZSP2pejbmK7cN4IqL+p19Aqi6kmmUWiuadjyf7n0GXfVX+3wrxrBRoWaKbbWLgR172b/4+P1O6YFG+pt2O+4HKBmWsWD1w42FZEzNZSuSqHqZfitBa73uqojVm8Ft57YTyLYNRAEiXwHMGMNQIlyRw6AQkFa5eSGiSDBNEyC+f2ZlDe8dTzkNWOymMXYfZN0woVkDectFfTi3nMeSTwaxBDgd2lJRpaK7ll2iwT3hzGeLyJUM38gWt2nIbtBzIV2CTACskuTJiOKUYG4E775twbueGrD93FQLx25nzP1a9r2euQBE9AxjbF9Dn+lGwc8Yw0qh+ZKw3Yb3AbL74W4b5MW3Gi1vLChjC6RqhGQeX64ZDFv6IhhJRlDQDLy+mIdpMiiyxHsiaKZjPnNPJu56RrIEjCYjCCsy7rr5Ctz+1WeguqKhbAjAWy7ur1hBeWsfRUN2SC3DnrEkHv7oAQD+gnw2ozp5DV7h650oTlzIoGj5NCKKdQzGsHs0Ebhqq43fJNRNSlk7J6a17Hu9J8xmBH9X2fhtNoIX3R2ZcceDz0LVDcduy3u48kYiQui3Fm9bXdnW2AEs5dRyxUrGTVV2lqAt9Hn2cfmxcJreWA7lTLHcLS0ZlqpmN7o/64dp8l7OuskwOZNx/ALeyqyzGRV9URmjKf9uY15fRMlVyY+s5kMSwSnB3AhrDW1sdzhso/6Q9dp3O8fVKrpO8HciFr/VeEMBvXHnzAxm9xY0jncB626w7m54Ys8Qdg6DF78QSZ3x+PjXF/M49vpS9X4GAbzaqmFWHMN2LLvrT715PAUCkC0ayBTLZbDdwtcv/BDg81lRM1DSjaYrqK4ltHE9wmHbGXO/ln33Qi5AVwl+w2QbwovunfG9cefuR8mroQqqQ7T6epmsvm/E9xpbM0RIqp8X4f24naJQ0s2KaCy3ps+sXsGew/miSHA0Qu+94+2lDFQKX29yoDuVws7XMBgwmgzXPEc/1pIguR5abztj7tey717IBWib4Cei7UT0XSJ6mYheIqI7633msvFUSxufdArvjD+ailQ0GLeFgB2aKAgGY6sFaEmvfQEZmFMfPySXG8+DyHKA8h/J2b5yMvEvrU8Y74860VgRRUY0xH/s7V+ZySASkhGSqu2jjETkaITee8eOYbd7KXuFrzf8MBKSKyYhAs/5cIeyBmUtoY3rofW2M3N/LfvuhooC9WhnOKcO4D8xxo4RUQrAM0T0bcbYy208ZlfgDQVMRUMYSRlO82ouiBjCcmWlzY0wB/RHFd42cA0noxA3q9hs64/g2l1DTvljiYCITOiPh7GQU32botvY9efjMRlb+6JOCWseJy47UT2xsOxE9egmzxGIWQ5Yd7tHiXhHrlQ0hOFECHNZzUme0i2H8GBcwcUDXOtLFzQMJ8LIqQZmMyX0RWUs5DTHcSwT9y+4NUL3vdMXC6GkG8ir1Rufu8MP9999BHISFdFKI8kwsqXV0UVBaDa0MUg47FojX9pZ838t++6FXgTrFtVDRF8H8EeMsW9X22bfvn3s6NGj6zKedlLPq++uqOjOtrQTby7/1CPIu7JagfLCQCJX43WLkEw1hV8tJPCV1muzWd7YPMBnIoqEV//b+2pGfdgZkLOZIkyrDSBj3KzhF7r6M2/bWrOOvRv7uDMrRegGQ8mo3J+tYRMRLuqP+j549SJW6kVi2duXNAPpou7UxumLytg2mPDdp31M3WB4Y6UACQQGnrQ21se1aQBrigjplkices9AN4SKbhS6NqqHiHYCuBrAD9bjeJ2m3hK5Xrblr/zYBBfwKP8A/MuyzQqOcANvvbhzOM7jwwlO8hQAnrWaDOPN46sLrRKAr/z7a/HwRw/gy7+4D3vHU85kJMH/5iAAv3bwEgC1l7T2Ndg5FIfBeAXOrf0RbB2IIR6SnPHJEjUk9IGyGcHuyat4bCm8rAPhzht3V+2bXG857v0Odw7FMRgPQZaoYvv/9q/eiuc/8xP40X//SWztj+LigUo7rtu8YR9TkfmERBI3Fe0aSTj3x1ozR7vFzFDvPHoh8mUj03aNn4iSAL4H4PcYY3/n8/fbANwGADt27Hj7mTNn2jqeXuGLj77mtJRMhLmteKVoQCYCkd2LliEkAZeN960yK9RKdqm1/PRu4018+vD+XRXlAJrZ51qXvV6NX5IIql4udCdLFKhsQaPjqrd9EG271deiFefVCfbffQQDsVBF1NVGy99ZL7ougYuIQgC+AeBbjLE/rLf9RjH1tIPHJmfxscPHuf3c6ngUlggXDcSqCvqNim0mUHUDCzmVv8mAkVQYIVnumLlAmC+C0y0mqY1AM4K/bc5d4lP5lwG8EkToC2pzcO8YPnfLVV2vya0HbueZbqShGgxhRcLO4WRHr0kvOPW6hc1SBbNbaZvGT0T7AXwfwAsoJ0b+FmPsm9U+IzR+QTU2Qx/UzUYvmKR6ga4z9TSKEPwCP4QJRSCoTtdG9QgEa0FEgAgErUUIfkHX0wu1TwSCXkIIfkHX0wu1TwSCXkIIfkHX0y1JSQLBRqGrWi8K1sZGjXwRYZICQWsRgn+D4I58cdc/t3sX9zqboQ+qQLBeCFPPBkFEvggEgqAIwb9BEJEvAoEgKELwbxBE5ItAIAiKEPwbBBH5IhAIgiIE/wZhrXXcBQLB5kFE9WwgROSLQCAIgtD4BQKBYJMhBL9AIBBsMoTgFwgEgk2GEPwCgUCwyRCCXyAQCDYZQvALBALBJkMIfoFAINhkCMEvEAgEmwwh+AUCgWCTIQS/QCAQbDKE4BcIBIJNRtsEPxF9hYhmiejFdh1DIBAIBI3TTo3/zwHc1Mb9CwQCgaAJ2ib4GWOPA1hs1/4FAoFA0Bwdt/ET0W1EdJSIjs7NzXV6OAKBQLDh6bjgZ4zdxxjbxxjbNzo62unhCAQCwYan44JfIBAIBOuLEPwCgUCwyWhnOOcDAJ4EcBkRTRPRh9p1LIFAIBAEp209dxljh9q1b4FAIBA0jzD1CAQCwSZDCH6BQCDYZAjBLxAIBJsMIfgFAoFgkyEEv0AgEGwyhOAXCASCTYYQ/AKBQLDJEIJfIBAINhlC8AsEAsEmQwh+gUAg2GQIwS8QCASbDCH4BQKBYJMhBL9AIBBsMoTgFwgEgk2GEPwCgUCwyRCCXyAQCDYZQvALBALBJkMIfoFAINhkCMEvEAgEmwwh+AUCgWCTIQS/QCAQbDKE4BcIBIJNhhD8AoFAsMlQ2rlzIroJwD0AZAD3M8Y+287jbTYem5zFvY9P4exSHtsH47j9wAQO7h3ruWO0gyDj9m5zw8QQnpxarPq6W869V78TQfdAjLH27JhIBvAagH8JYBrA0wAOMcZervaZffv2saNHj7ZlPBuNxyZn8emHXkJIJsRCMgqaAc1guOvmK1omBNbjGO0gyLi92yzkSpjNqBhNhjGSjKx63S3n3qvfiaB9ENEzjLF9jXymnaae6wCcZIxNMcZUAA8C+Ok2Hm9Tce/jUwjJhHhYARH/PyQT7n18qqeO0Q6CjNu7TbqgQyIgU9R9X3fLuffqdyLoLtop+C8GcNb1etp6rwIiuo2IjhLR0bm5uTYOZ2NxdimPWEiueC8WkjG9lO+pY7SDIOP2bqMaJiTi//u99ttHJ+jV70TQXXTcucsYu48xto8xtm90dLTTw+kZtg/GUdCMivcKmoFtg/GeOkY7CDJu7zZhWYLJ+P9+r/320Ql69TsRdBftFPznAGx3vd5mvSdoAbcfmIBmMORVHYzx/zWD4fYDEz11jHYQZNzebfpiCkwGpKKK7+tuOfde/U4E3UU7nbsKuHP3PeAC/2kAP88Ye6naZ4RztzHs6I7ppTy2tTmqp53HaAdBxu3dxo7iqfa6W869V78TQXtoxrnbNsEPAET0kwC+AB7O+RXG2O/V2l4IfoFAIGiMZgR/W+P4GWPfBPDNdh5DIBAIBI3RceeuQCAQCNYXIfgFAoFgkyEEv0AgEGwyhOAXCASCTUZbo3oahYjmAJzp9DhqMAJgvtODCIAYZ2sR42wtYpyt5TLGWKqRD7Q1qqdRGGNdnbpLREcbDZvqBGKcrUWMs7WIcbYWImo4Bl6YegQCgWCTIQS/QCAQbDKE4G+M+zo9gICIcbYWMc7WIsbZWhoeZ1c5dwUCgUDQfoTGLxAIBJsMIfgFAoFgkyEEfxWI6CtENEtEL7reGyKibxPRCev/wQ6PcTsRfZeIXiail4jozm4cpzWmKBH9kIiOW2P9Hev9XUT0AyI6SUT/h4jCXTBWmYieJaJvdOsYAYCIThPRC0T0nB3S16Xf/QARHSaiSSJ6hYhu6LZxEtFl1nW0f9JE9NFuG6c11t+wnqEXiegB69lq6B4Vgr86fw7gJs97nwTwHcbYHgDfsV53Eh3Af2KMXQ7gegC/RkSXo/vGCQAlADcyxq4C8DYANxHR9QDuBvB5xthuAEsAPtS5ITrcCeAV1+tuHKPNuxljb3PFm3fjd38PgEcYY3sBXAV+bbtqnIyxV63r+DYAbweQB/A1dNk4iehiAHcA2McYewt4yftb0eg9yhgTP1V+AOwE8KLr9asAtlq/bwXwaqfH6Bnv1wH8yx4YZxzAMQDvAM+MVKz3bwDwrQ6PbRv4A34jgG8AoG4bo2uspwGMeN7rqu8eQD+AU7ACSbp1nJ6x/TiAf+rGcaLcy3wIPAH3GwB+otF7VGj8jbGFMXbe+n0GwJZODsYNEe0EcDWAH6BLx2mZUJ4DMAvg2wB+BGCZMaZbm0yD39id5AsA/jMAu8v6MLpvjDYMwD8S0TNEdJv1Xrd997sAzAH4M8t8dj8RJdB943RzK4AHrN+7apyMsXMAPgfgdQDnAawAeAYN3qNC8DcJ41NrV8TCElESwP8F8FHGWNr9t24aJ2PMYHwpvQ3AdQD2dnZElRDRBwDMMsae6fRYArKfMXYNgPeBm/kOuP/YJd+9AuAaAH/KGLsaQA4ec0mXjBMAYNnGbwbwt96/dcM4LR/DT4NPqBcBSGC1SbouQvA3xgUi2goA1v+zHR4PiCgELvT/N2Ps76y3u26cbhhjywC+C74kHbD6MwN8QjjXqXEBeBeAm4noNIAHwc0996C7xuhgaX9gjM2C26OvQ/d999MAphljP7BeHwafCLptnDbvA3CMMXbBet1t43wvgFOMsTnGmAbg78Dv24buUSH4G+MhAL9k/f5L4Db1jkFEBODLAF5hjP2h609dNU4AIKJRIhqwfo+B+yJeAZ8AbrE26+hYGWO/yRjbxhjbCb7cP8IY+3foojHaEFGCiFL27+B26RfRZd89Y2wGwFkiusx66z0AXkaXjdPFIZTNPED3jfN1ANcTUdx6/u3r2dg92mlHSrf+gH/55wFo4FrLh8Dtvd8BcALAowCGOjzG/eBLz+cBPGf9/GS3jdMa65UAnrXG+iKAT1vvTwD4IYCT4MvrSKfHao3rIIBvdOsYrTEdt35eAvBfrPe78bt/G4Cj1nf//wAMduk4EwAWAPS73uvGcf4OgEnrOforAJFG71FRskEgEAg2GcLUIxAIBJsMIfgFAoFgkyEEv0AgEGwyhOAXCASCTYYQ/AKBQLDJEIJfIABARP+KiBgRdVU2sUDQDoTgFwg4hwA8Yf0vEGxohOAXbHqsWkf7wZP0brXek4joT6wa8t8mom8S0S3W395ORN+ziqN9y07pFwh6BSH4BQJe9OoRxthrABaI6O0Afha8LPflAD4IXlfIro30JQC3MMbeDuArAH6vE4MWCJpFqb+JQLDhOQRejA3gxdkOgT8bf8sYMwHMENF3rb9fBuAtAL7NS6VABi/tIRD0DELwCzY1RDQEXoXzrUTEwAU5A6926fsRAC8xxm5YpyEKBC1HmHoEm51bAPwVY+xNjLGdjLHt4B2jFgH8a8vWvwW8aBvAOzKNEpFj+iGiKzoxcIGgWYTgF2x2DmG1dv9/AYyDV2V9GcBXwVtFrjDGVPDJ4m4iOg5eEfWd6zZagaAFiOqcAkEViCjJGMsS0TB4ydt3MV5fXiDoaYSNXyCozjes5jFhAL8rhL5goyA0foFAINhkCBu/QCAQbDKE4BcIBIJNhhD8AoFAsMkQgl8gEAg2GULwCwQCwSbj/wfn/9YPutUutgAAAABJRU5ErkJggg==\n",
+ "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAABeU0lEQVR4nO29e5gc1Xnn/32rqu/dc5/RCF0sjQXIYHOzwMZWsIy9CcQbnM2yWZSN42TtQJ7NGuxde+1kY6/DJrvmFyc29mazsNi5eQNJtHHMOoBjLGOMAwYhEOYyIDESaCRGc5++d93O749TVV1dU91d3dM93T1zPs8jzXRPddWp6qr3vOe9EmMMAoFAINg8SJ0egEAgEAjWFyH4BQKBYJMhBL9AIBBsMoTgFwgEgk2GEPwCgUCwyVA6PQA3IyMjbNeuXZ0ehkAgEPQMTz/99DxjbLSRz3SV4N+1axeOHDnS6WEIBAJBz0BErzX6GWHqEQgEgk2GEPwCgUCwyRCCXyAQCDYZQvALBALBJkMIfoFAINhkdFVUj0DQLF95+BXc89hJ5FQDibCM9+0dxUxaxemlPHYMxnHLNRM4sHes08MUCLoCofELep6vPPwK7jx8AgXNgCIBuZKObz77BiZnVjAQC2E2U8Tn7n8Bj0zOdnqoAkFXIAS/oOe557GTkAhQJAkSSbALja8UdBAR4mEFIZlw16NTHR2nQNAtCMEv6HlyqgGJyq/tFhOmq9VELCRjeim/vgMTCLoUIfgFPU8iLFcIebImAfdkUNAMbB+Mr+/ABIIuRQh+Qc/z0f27YTJAN02YzIQt7/tjChhjyKs6NIPhlmsmOjpOgaBbEFE9gp7n1vdfAADlqJ6I4kT1TC/lsV1E9QgEFVA39dzdt28fE0XaBAKBIDhE9DRjbF8jn2mrqYeIBojoEBFNEtFLRHR1O48nEAgEgvq029RzJ4CHGGM3ElEYgPCuCQQCQYdpm+Anon4A1wD4VQBgjKkA1HYdTyAQCATBaKepZzeAOQB/SkTPENE9RJTwbkRENxPRESI6Mjc318bhCAQCgQBor+BXAFwB4E8YY5cDyAH4jHcjxtjdjLF9jLF9o6MNdQ8TCAQCQRO0U/BPA5hmjP3Yen0IfCIQCAQCQQdpm+BnjM0AOE1EF1pvvQ/Ai+06nkAgEAiC0e6ono8B+D9WRM8UgF9r8/EEAoFAUIe2Cn7G2LMAGkosEAgEAkF7EbV6BAKBYJMhBL9AIBBsMoTgFwgEgk2GEPwCgUCwyRBlmXuIRyZncdejU04D8asnhvD41GLVhuLe7f1KEze6T4GgkwS5pwX1EWWZe4RHJmfxuftfQEgmxEIyFnIlzGZUjCbDGElGUNAMaAbD7TdcjAN7x1Zt7/273z7nsyXMZVWMpcIYTqzep0DQSYLc05uRrivLLGgddz06hZDMG4cTEdIFHRIBmaJ/Q3Hv9n4Nx73bZIo6wBjOpUt4+VwGMytFaIYhmpR3OY9MzuLg3U9g/x2HcfDuJ/DI5Gynh9QWgtzTgmAIwd8jnF7KIxaSndeqYUIi/tPG3VDcu733737bFHUTBuNNymWJoJsM8xkVx2cz7TotwRqxteDZTBEDsRBmM0V87v4XNqTwD3JPC4IhBH+PsGMwjoJmOK/DsgST8Z827obi3u29f/fbxjb7SQQQCBIRQICqmxB0J5tJCw5yTwuCIQR/j3DLNRNYKWg4PpvB5EwammFCNxhSUf+G4rdcMwHN4O9Xazju3QaWu0ciAmMMpsnfCMu0ajwbgVaYSDptZtlMWnCQe1oQDCH4ewgCAMY1c4kIfTEFA7EQVgoaxlLRCifXgb1juP2GizGWivr+3W+bZFTBYFxBRJFgMAZFJgwnwjh/S19nTriNtMJE0g1mls2kBQe5pwXBEFE9HaLRsLSDdz+B2UwR8XA5Ajev6hhLRXHvze9s2ZgajZpoJmS0G0LwWnE9/fYxlykirxroi4XW5VxFpItARPX0CM1oiuuxpG9UowpyHt2gFfvRiuvp3Ue6oGEhpyKn6ut2rkILFjSDSODqAG6HHADEwwryqo67Hp2q+sDuGIyv0i7bsaQ/sHcssNAIch7NnOt60Irr6d3HfLYEAIgqsuNoXY9zbeQ7EwgAofF3hGa0zW50bDUTMuq3TSdoxfX07qOkmwADRlMRZ5tuOFeBwIsQ/B2gGYdcNy7pmwkZ9dumE7Tienr3EQ/LGEmFkYqGnG264VwFAi/C1NMBbrlmAp+7/wXkVb3CIVdP2+y2JX2Q82j2XNeDVlxP9z5sf0Y3nqtA4EZo/B2gG7X3ZmgmZLRXzzUIm+lcBb2NCOcUCASCHqaZcE5h6hF0nG6M828Xm+lcBd2LMPUIOkq3xvm3g810roLuRmj8go7SyTj/etp3q7Xzbs1paBViNdM7tFXwE9EpABkABgC9nh0qU9Rx8O4nRDeoNtHMg9nuh/n0Uh4DsVDFe+sR++4udeDWvm8HVjWy8ft7M3TqXNeDdlwvQftYD1PPexljlwVxPpxdLjjL4FMLWdx5+AROzmfFsrgFNGNmWA/TRKfi/OuVM25HueNuzWloBXc9OgXNMDCzUhRNfHqArrLxEyFwhylBYzQjyNaj1nunMpLrZRS3I+O4G7OvW8Xx2QzmMyp0k4kmPj1Au238DMA/EhEDcBdj7O5aG0tUrvter8OUoDFOL+UhEzA1l4VqmAjLEkaS4ZrXcz1MEwf2juF28ElmeimP7U2akxptGl+vVs+OwThOLWSRLujO9YooEnSTYf8dh5sye7XqXLsRVTcBKj/DRIBJTDTx6VLaLfj3M8bOENEYgO8S0SRj7FH3BkR0M4CbASAyuMV5PyxLzgNns1GWxZ0gFVFwfDYLWSJHIzuzXMT5Y8mqn+nGwnB+eO3LJ+ezePLUotM03s/eXC+j+OqJITx5ahES8Y5kRc1ATjUwGFfWZMNuR/Z1NzhVQzKhoAGmyUAE2OlBG7WJT6/TVlMPY+yM9XMWwDcBXOWzzd2MsX2MsX2heL+zDO6LKTAZqnaYEjSGk6jHXP/c7/vQK6YJv6bxEgHpQnUzYb0s28enFjGaDDstLhn4BKDqrKtMj90SInrBlj4MJ8JQZNrwTXw2Am3T+IkoAUBijGWs338awO21PnPeQAxjqSiml/LYNZzEwSv5cn2jLYs7QVY1sG0givms6qykxpMR5FSj6md6xTThNUkFNRPW0r5PL+UxkoxgNBUFAEzOpLvS9NgtIaL2Cmq8XxF1inqAdpp6tgD4JnGbnwLgrxhjD9X6QCqqrOp+dGvbhtc83RgWWQ/bbDMxWjbt2B2natFtheH88JqkwrLkTGg/ObMCiYC+qIK3bO1f0z670fTYLSGivaIkCDhtM/UwxqYYY5da/y5mjP1+u461nnRrWGQ9esVs0wzec/Oar0wGLBd0jPeFm95nKspNj32x7jI9dlOI6IG9Y7j35nfih5++Fvfe/E4h9LsYkbnbIM0srZv5TKtXCNU0MgAVSXOd0NLWeq7ecysZDJLlU2SMR5gQAd+bnGtoHDdesc0xNe4eSeKXruoO06N7nKmIgpWCBgAdN7F0elUrCI6oztkg++84jIFYCOQKPWWMYaWg4YefvrYln3lkchafOnQMmaIO3TShSBJSUQUfeuebWprJ3A2Nutsxhjf/9gNQJECi8oLWZCZ0E3j1v/3suo2jHfiNM13QMJwII6caHZuQeuX6bUREs/V1oJmldaOfueOhSSzlNTAAiiyBAVjMqfjq90+01Fy0HglafjwyOYuDdz+B/Xccxq33PQNVN1o6hkRYhunRZ0zG369Gp65Fo/iNU5EJZ1eK6KQK1yvXT8ARgr9BmrGVN/qZqfmcFT9OIBAkIhgM0AzW0gerE/1wvf6OnKpjIacibZkrWjGGj+7fDZMBumlamr4Jk/H3q9GtvYG9eMeZKWqYz6jIq0ZHwzl75foJOELwN0gzXZZa1ZnJmwqz1gerE45Br2YYVbiwmM+WWjaGW99/AW67dg9iIRm6ya/Tbdfuwa3vv6DqZ7rJSVoL7zjnMiWAgIgidVTT7pXrJ+AI524TNBPi2Mhndg/HcWIuB/JkQSqeaXqtD1Yn+uF6ww9HUxGcWSqgpJtgjLVsDLe+/4Kagt5LN/cGduMdZ1E3IBFhJBlxtumEpt0r10/AEYK/C/nM9W/BJw8dQ7akw7CKXqUUGdGQ3NIHqxOx1341cFJRBbrJnd2dck72Shy6d5yJsIJ4WEafazL1UwjaHXHTK9dPwBFRPV2K/aB6Qy/d761Hv4JWC4yvPPwK7jx8wqmBYzL+r54pph1NU9Yj/LAdx3DvMxmWsZBT0RcLVY2mERE3G5tmonqE4O8SGhUQ6/Ewt+MYB+9+ApMzK1gp6DAZF/79MQV7x/tXZW27x+FdASUjCr5446WrmqYEHWevXj+/fa4UNIwmI8iWdF9N++DdT6wqtmdnbVe75oLeQYRz9ijNZPauR/hcO45xfDaDbNFASJYQDUkIyRKyRaNm3fYvPPgSlvMamAnIRGAmsJzX8IUHX3LG2WgTkE5dP80wcOt9z2D/HYdx8O4nGo6+8dtnfyyEgXi4asasiLgReOl5G/9GyBZsJrO3XTVa3NdzLlPCeF+k4u9BjlHrO2mmbvvJhTw3DUnlzzCT4eQCH8fx2QxW8hokV8np+YwKzag+mfhdP90wcfT1pabr7dc7hh16yQDsHIo3Vdq5me+9FeW1N8JzJijT04J/o/T5DPowux++dEGDbphO9Uhg7VE+3oxhkwHTSwXsIEIqyse3kCshVzKqCsd630k76rb7TSYaM7GU16qO0+tklgBoJgNjwMxKAfOZEj516Bj+wDInNYNX4Dqhl3I59HIuU8St9z2DvlgokEBtRojfcs0EbvvrZ5Au5MDAw4L7Ygo++4GLnG1qCfZWPWdBJg8xwawPPS34u6Uk7VoJ8jB7Hz7DNDGbUQEAI8mIb5TPVx5+Bfc8dhI51UAiLOOj+3fXdKDaGcOyRFBkCbrByxycXS7ggi0KFnIlzGZUjCbDjgD45KFjGE1GkCnp2DEYx1KuVPM7uWBLH07OZ5EpuqN6Qtg9Uq4a6n34x1IRvLFSBLHyZGEy4PzRBIDVTUB0g09aMlhVQeVttFLS+QzEz53AGLCU13DHQ5NN30v1Qi/TBQ0LORUmYzVXAN7aPHPpIlSDVZTzqCXEx/vCSBf0irGlCzqem152fCTuCd876bXiOQsyefht86lDxzCcCCOrGmIiaCE9bePfKLbLIJm9Xjt2uqCjP6ogrxq+SWFfefgVfOnh40gXuUM0XdTxpYeP4ysPv1J1HFPzOTDGoBkmShoXngSuCa8UNORKBkaTYYymoiAi6AbDcl7Dyfmc86Aen8tCNyrNNu7v5JZrJhBWZIz3R3HhlhTG+6MIK7JzrrYj95nTSziXLuKZ00tYzqtIRGQQuFAnAIPxED593V4Aq5uAmABkAqIhuar93ttoxYYx5mRLS8SvSbN4E/cSYQXDibATemknrUWV6uP0+n+WCyqyqgHDNHntJ0JFqQY/f9HfP/sGyLoesRAPC1Zkwj2PnQTgXyLEnvSA1jxnQfwwXv+FYTIs5TWcWsx3NCt5I9LTgn+jZAsGyez1a2adKeoIK5KvU++PH3l1Ve0WZr1fDcM0YVjaNEP5p0zADz99LfpioYpEoflsCRIBBnN1pZIknMuUKvbr/k7qnaufI7egmUiFZVy+cxBb+2O4fOdghQnGO5nYJh/DZJicSWPKmozcgur0Uh4Rb0YcyqanVuEuVfyVmy5HWJGdCb6kmwDjSWw2XoHqFYbpgg5ZIoQVGXvH+3D+WAr9sZAjQP2cv66Gaw4SwelZ4FcixD3pteI5C9KM3TvBzGWs+8tc365n7lpSzTjge4GeNvVspGzBepm9jTpFSw2+DwBhWYZuru7IFZb5w+g1SamW9u1uTrKlL4Lp5WLN76TWuVZz5M7lNDz2W/6hh97kobAsoaQbfNKq0l/Y24NYM7hoZIBV059PfHtGWqdEeMcZD8tIRGTHfwKsFqiNdhfz8xcRVk9o9YrWuWnFcxbk/g1yf61XLale9xvWo6c1/lbVwOkFQpbz0zR5oxHTsk804xStqslU25X1vtckJUsE06OxKrKEC8aSFd/JjVdsw12PTrVVg3Jr1juHYtwMUqO/sLcHsULlUzUYA0nAQDyEz1z/lraN8ys3XY6QLNc08Xm1bds0Va0TmJ92PhjnE0G1onW7h+M8kc51b5mMv2+Pea3PWZD7d9X9Rfz+cq8y17uW1EatMtrTGj/QG60BW0EQp6ibsERQvbWJLappMvVWCV6NdddQHAs5FbJEFXV2PvuBvU1HhEyMJHB8NlvVkRuEIP2FvdtEQjIGwxLSRQNjqci6lBwIUubglmsm8MlDx3BmucBNHgB0g2E4ofjWNvLTzpPREN5zwQi+Nznn6+j3KxEyEKmc9Nb6nAW5f73XY/dIAnPZkuVsb10dp1p0SyvLdtPzgn+z0Ggz639/7R780cPHV72/JRVxNJm8quOOhyadCBCjykThft8rAPxKS7j/3mhEyKev21uOMDF41IrbkRuEIP2Fq22zZ2x9s1mDCFQCAAZrlSWhLyZjIBbyrW3UTM2cA3vH8MUbL21rnZ2g92+j91eraUXOQy8gBH+P0OgDbWtzdjinaTKMJsMY6ysLP90wcWohj13DcQzEQpheKvjuS5aqm5PqCa5GNagDe8fwB2sUQkFs0t3iH6oXt37Xo1Poi4Uw3h9z3surOgYTETz0idoTVCN+6navnJst4lZvXK2O+++W+6LdiFo9mwS/ei3Hz2UAAs4fSwEAnj+z4issUhEZP/nd61p23PWoExNEU1xvbdJvjPVq+TTTtrOZ+kC9mDjVrnpLnb4vGkUUaethvA9eqytv+j0kpxZy2D4QQ18sDACYmsuioBqwLf1BCqg1c1y7R6ydlLMeVUY7RS2BGmRS9NtmPltErmT4Zvs2M9H2avVOUXyO04zgF6aeLsDrAD05n8WTpxYxlgpjOBFpSUjZgb1juHF6uSKTd2sqAsUVHTKSjODMcgFhibBnLOkI6aVcySl9EERIe4XdjVdsw+NTi5heyiMZ4XHlmsmqnuut9x6FbgIlwwyUcdwqWq311nNsBzGDeU0PftnT7uzWZuor9WoG/GZxxLaDtodzEpFMRM8Q0bfbfaxexRtClinqkIin1bcqpOyRyVkcOnoGo6kI3jKewmgqgpLBTQZ2+JwiEwbiIeweSWCloCFsZXG6hfSdh0/g1EK2aialX+booaNncMs1E/jhp6/FQDyM/lio6rlmChrSJQN5zYAiccfanYdP1Mw4bgXNVEitR73QwCCJUd5QSm/2tDe7lQCcWS5W9DCu55zs1Qz4jZLA2QnWQ+O/DcBLAPrW4Vg9SaNJOjaNaKh+Wh3A48EH4mHHnvnZD1xUYTZQDdPZ1i2kR5JR3yJj9Wr1+J0rGENeNTA5k3YSqQiARBIk4vHn/+sHUzVXGvXqEtX7e7Nab6391tNIqzkSx/vCuOTz33H2+b69owC4s7akm05cPrA6u3W8P4rppQLOZYpIRZWqzklvwb9cSYOqMyfUsi+mYNewf6hwt7BZHLHtoK0aPxFtB/ABAPe08zi9TqNJOkDjGmo1rS5b0p2EIm/ZB+9nvBOSXWQsp+qBa/V4z1UigsYrFzhCfxWMIa8ZVc/V7upVqLJKqPf3WtenltZbb7/1NFK/xKi37+zH/c/NOPvMlXR889k3MDmz4qvR29+Ffa+koiFsG4iCMVRNtvLeO4pEWMrrKGiG8/3OZlRcPTFU9dy7gc2UwNlq2q3xfxnAfwKQqrYBEd0M4GYA2LlzZ5uH0514NZdUVMFcVkVfzD9JB2hcQ90xGPdJoFFWJYDVKv0cliXns0CVImNWrR7bYQxUCjvvuWpVksbcU4BmckdztXO957GTkAhQJD4ue5Vwz2Mncev7L6j7d/v6NBq/XW+/QTRSb7jiJZ//TsU+GfjEsVLQsX1wtUYvEy9J4c5uVWQJV+wcrOrg9N47Jd2ETOX6TPa98fjUIm6tevbdwXolcPZi1FMt2ib4ieifA5hljD1NRAeqbccYuxvA3QBw4VsvYwfvfqLjF3e9+7v6OV5vuGQcM2m1akhZo44tbxli1TAxl1XxS1eVtbp6pZ+9E1JJ57VU3CUbtvRFcHqpgOPnMr5lg73x3IwAGXaNHL4PW+ibzHQqZw4nqp9rTuXasRt3EbJ6fweaMxv47ZcxXgnVdoa7HdtBQgO9+7SviX0duEbPMJMuYaWgNZXd6mduU2ReHmHveJ9zHt1u418vNmL9nnZq/O8GcAMR/SyAKIA+IvoGY+yXq33g7HIB2z3L+fW+uPW+5GZugiD7tB2vOy2h8/TrKzWXrY1qqI9P8cgZu/GIbcd1a3VeTXAkyTX9XMmwhEwSv3TVkCPI/IqMqYbJE74IvmWDgUot7cLfedDJSLUp6QZMBugmLyS2JaUg4jHDuM81EZYdM4WNuwhZvb/bY2o0wci7X7t/AQEVju1GzA/efdplK9xj92r0jcade+8d70oOEE5SN70a9VSLtgl+xthvAfgtALA0/k/WEvp8u+rL+fWi3pfczE3Qjn02qqGeXspjOBFxhDmwWqvzW0UMJyJQpMpkIXuisCc09xgWczxG390ZrNa5+NXmkYhw4ZYkHvz4NVWP4z7Xj+7fjTsPn4Buch+EadX3sYuQ1fu7TaNmA+9+dUstH/OUxWjkHvbu05b3/TXMfo2O23vv9MUUzGZUpKLVj7GZ2Yhho11VnVNyZScCnbm49Zx8zTgB27HPRh1bQULfGg2P8xtDMiJX2Jvrncunr9uLwXioapOVIOd66/svwG3X7kEsJEM3+fFuu3aPY7+v9/dm8e4XAMY8ZTEavYe9+0xEFPyLy7Zi73h/yxyY3uu5aziJ267dg90jSeEk9WEjho12VeZu344L2Xv+09ec194svPVwsNTLBmwmW7Ad+2wUu7OVuwJjMqLgi66GJkGybOtd82YzR3spRb4aIpN0Y9Ltmc3NZO52lcbPGKrWJm9Hgo0f9dogBmmTuB77bAZ3lUew1eX3vZqgN4EryDVv5lzcNeq9IaXN0okuSuv1PQrWl40YNtpVGv+WiYtY8qYvOhrpVW8aAEhyQgvjYXmV7TgkEQYTkcCrgCCrhnoaaDMaar0Eok/cdxT3PzfjnPsNl4zjSzdd0eAVrE47VirV6LQG3yoNrZn6Sa0493bXbeoUGy0kslvo+SJtka3ns53/9k5IxJN5TAYMxUM4byCGl2bSkIhwXn/MaVadLqiYXi5g13Ai0APeqSVbvePaiUB2qKXtfGyFHdqm0SqPzX6mG2iFycX7nblr5IwkI227d7zHnc+WMJdVnVpG3WZmCEq3m0t6mZ439QA89VzVmRO3vFzQQESIKtz5aScNAcC5TAkhSQrcJs2vdoqqG7j1vmfaahKoV7PFnQgkkWT95O/XohFzRjMOqvVyarXaLNOK2jN+Tc4l4mUr2tmSbz3qNnWCzdLSsFfoOsFvssq4b3sCGE1FAMazDN320y0NVCL0CgS/kgPt8BvUE0Q5tTLGHFidYOSlUZ9HO3wTtcYWVJC3w3fTigmrXrkKoLmos3rXpl3HXStrnZx7tRBct2CaDJphoqgZyKs60kUNy3kVCy5FuBECx/ETUQzATsbYy00dKfBx+E+vBSoVDWEkZTjJRNsH407iiZtaD7g3ccWv5IC3HWErbJH1kq0SYRk5VQdjBhjj14AISISrfz1BYv9rlUcO2pKv0aSmRhPc2pEc04riXe1IcgpybbylNQyTQbfuicmZdNVSG+2kFZmrvd7SsNX+CdNk0E0GkzEYJoNhNaDXTf7TsN43TcBgvEF9Kwkk+Ino5wB8EUAYwG4iugzA7YyxG1o6GmBVmqdEcJJKQrKMr9x0yarww1oPuPsLS0UUrFjFrWIh2bfkgLcdYSsyiOsJovftHcU3n32jfAkYvw52VUY/6iWV+D2sjWaRAo0nB9316BRU3cBCtrImUDVB3o7kmGbb/LlpR5JTkEnOW1pDs54HXq3Uv9RGu2nF5Oz3DKQLGkISOeUtutXZ28jEZ5gMumlN2GZZmHtfd9q3GlTj/zyAqwA8AgCMsWeJaHetDzSDLBHIcm7yolwyJCLfptJA/Qfc+4UVNAMEICTxffqVHDiXrl1WuBnqjXMmrWIgpiBd1J1z74sqmEmrVfdZT4PqVJr5K+fSSBd1SCBeQMxgWMip0I10U+dRjXoa2FqLd3m/s13DSRy8cqihFZOXIJOct7SGt8yEX6mNdtOKydl7PZMRBUXNwKnFPAyTYT5bwicPHavIK+kGGGP4Xz94FYrELQMmAyKKDMPU8dXDJ7B3ax/Xzg3WFs28XQQV/BpjbIUqM2tbfoYyEXaPJOCXPAQAz00vr3rYvYNxb5MuaJAlVNQZT0UVp1G174rBNJEISXjh7IojhIcTIaiuKpLNLPtqCaLTliDxRs/UerDqrSI6lWZul1aWLIlFxJe1qqvkcq1VWBBNulVFs7zf43hfGN+bnPOtgw8Al2wfWFOUVZBJzltaY3Im7UR6daqA2o7BOE4tZFfVeWq0Xr/7GbjuSz9AQTMhE1cQmAks5zV84cGX1kXwM1apiRtGWVM3GINulM0wJxdy6Isq0Fxm5ZBMOLOcR17V2z7WdhBU8L9ARL8EQCai88HLtfxTqwdz3kAMY6lo4BZ9nzx0DASgLxbCQCyEUwt8Gzvk7uxSASYARSIo0mrt008TX86VsJAvf5kmA+ayGrb3c9tus0Kn1mTRjNbrV9Hzo/t3r2mfQcbtjSn3vgYAMMB01d4xTYaCamD/HYeRDMtYyKnOd+ZdhW239nnXo1P4nW89j1SEm1bcmcPNrGb8zuPQ0TPO9zg5s4LHp3TIBCgyOXXwB+MKtg3E18XkB7SvgFqjyop7ezCGc+kSZIkq6vUfvHKo6mfqHePkQp6bs1wKAjMZTi6sbUJjzGVWYW6hzirMMIYZXG/d2hfDQq5U4ZwuaibG+2JrGutaeXJqEfc9dRqh0V1va/SzQaN6PgbgYgAlAH8FYAXAxxs9WD3mMiVH+wX4rDqzUsTL5zKYzZQAsIqwtmxJx0pBc7Y5l+bb2CF3dmoqF0IESSKYjGGpoDvRCc9NLwMoa3WZEl9dECwnq/X+YoFPBs2EpdWLXGkmesavleKho2dq7tPdPzdIZIZ33N7Wi6cW+OuT8+VWjCXdRCqqQJF4W0D7+smWgD21mMdSXnM6RsXDCvpiIQwmIvjhp6/FLddM4NDRM5jNFCETcHw2ixNzOcgE57odn800FCHid/3/+JFXoeqG8z2uWN+vwfgK0V6grLQwjDJIBqj3e0tFFZgMTilsv3ujXsRNo5FT3u3nsypgFc+zzU2jyTAen1ps+hhBYYxBN0yUdAMF1UCmqGElr2Exp2IuU8LMShFnlgs4vZjHqfkcTs7n8PpiHmeXC5hZKWI+U8JSXkWmyFuMqrrZkNAHgJuu3AHd5H4dBv5TNxluunLHms5tLTw5tYg7Dx/HQq4EMLPhZUddjZ+IZAD/wBh7L4D/3Mwgg6KbpnPTnJzPcScnuSJ8rHBOG83gM7dmGGW7EwOK1jZkqZw8RNTynpuARNVXEbZJgjn/cfKWxtqOZtbNOCMb3ad3BRVEg/Uew9t60R3bPprirRgH4yEs5TVsH4whFpJxYjYLIsIWV49Yifgkb/tW3NfPfcypuSwv8cyA+ayKidGk8/AWNCPwasbvWummaY2bb+OWBd5wYjuaZiQZbpmJJa/qeOHsCj556BjOH0s537f3e/OWwq7nx/L7XhtdIXm3N6x+zIpEmLAumNfc1MgxTJNh93Acx2ezAExnZWiYDLtH4ji7XHC0crMLbOZXTQzhNpyP+546jZl0AeN9Mdx05Q5c1YYOZSZjKGo8ejFb0pEr6cirBnIlHTlVR7bEf//OCzPIq8aqsitBqSv4GWMGEZlE1M8YW2nyOIGQiBztynQJeze66wk1Xb+Ta4IwTIbJmTQYY45zzL6RZAKiIXlVcozdQ7YWA7EQ5jMlnFkuAiAngzhIM+t69vZGnZGN7tPbP7fag+lernsnOW9MuV+M+UgyAt0wHZMdAzAYVzCfLeHsSsFZBWiG4QhUt73YfV5ObX/wyXxqLmv16EVdv0Ct8wCAiCxVKBFu3PcSwBUM3TAwvVTA+WPVO5bVa9qTiii8aYpESBe5klZQDZxayFYIa797oZojN4jAbdTf490+LEvQDLPie3bf86bJ8PpiDv2xEEyT8aY6YAhJhNcWcphZKVaEKzLG8Kvv2o07vjOJXEmHbvAyJX2xED66fwJFrXr+Sqe4amKorqA3TLZKWGdLOnLW73lVR65kC/HK97KWYM+XjNY7T30IauPPAvgJEX0XQM5+kzHW9sACe0ZzXww7nK5ibvAJAzWAVWUfAB6yOTmThm7wiSHvEkIy8eW+97ghmRpqZu2mHTHMjTrcgjz8Xu3RO8l57c3V7M92PSUGHqm1mNOgyBJkiXj0g7WtPWnMrJSwkFXx5t9+AACQLWrYMZRAWJYc26xpcocbgduFNcPE2eUCVIMhEZZx8dYUbr3vGeRUAxFFgqYbMFhZk39tsYA3DZUn6/54CIs5zbG32w5UZ+Au7HvBYJVJdd6Kp/PZEj5271FsH4wjU9IdQd9v+TNOzGad1aibdEHHeH/1kNdak4vf96obJo6+vuSESSat5i5B77/tAzHMZoqIhRUwBgwnwpheKoABeOmNFcgSn2huuWbCWpkzjCajOLOcQ7ZkQDNMhGQJyYiMbQMJxwFq26TfSBewtS+Gn7/0PDxzeqVCi355Jo3/+g8voqAZiIVk/OLbt+ND79rlO06/fTajidvd5NwC2f49X9KRdb1X1r5Xb+f33bYCe2JPRGQkrJ9TcznoJkNYlnC2iX0GFfx/Z/1bd/xmv5dmMkiEZUiwtDNwDc27RI+FZIRk3pN0paAhIstcm5AJMgE6mCPkZYlvZzD+O2PlshGyaz3lbX0XxCwTxKnXqPPt6okhPDG14JyzZhjIqwYOXunftzjI5HPXo1PQjHIMvt3P9bXFSs2QYOKlN1YgWeGaYAZeemMFiiQhEpK4k24+B900nTr1huehIJQFLQOgGgwRhR9vuaAjfWYFZE3C9vaqaYIxICTbJQz45JEt6nj85JKz73yVjOfXFvOQJUJI5iWpf/PAmx0TSjKiQNVNFHXT954j4g4x7mvifOHBl7CUU51sc81gKAI4MZvBheN9ODGbRUk3sZBVa2pxOdXA64t5TC8VHJu4e5UwvZSHZk2A3sklXdCQLWrQDB65JoGvikOy5Jh+0gUNphVqGHXdf//mHTuxkC3BMBn+96Ov4t4nTyOvGY4CNJSMIBqSUNB0PolbKyGnHaQJJ3zx8h39eO7MMm8eQ4BmmFjImfjnb+sHULZJKxKhL6pgIVfCQy+ew23Xnu8I6r/8p1P48ydeg0SALPFObH/+xGsA4Ah/t6BPhBUs5kpIRBTEQhLeSBfwB//4Mn7ukq3YNhhHTtUrhDP/3XDedwvxRu3+QYmGJCQiCpK24I4oFUKch5Qr1j/Z2q78fjKiIOzt7+m5ns0QSPAzxv68qb03iH1z2q3n/L4LiYC3jKdQ0AwUNAOGwW9yorJtPyITLnCFvtlFxa7/8qPcrkioaAnorBisNyQw7D2P37C2ecGt1dZrZu2l0XyDIPb3v3nq9VXChFnv+4UcBpl8js9msJLXIEkEWSKomr8GYzI4oacM5d9h2fvBgJAiQZEl6GZ1c8re8T68cHbF/igkkiATgwETJnikjww+KTAAxPh7dsgoY0BYkVB0HcNrpvFim/sIleGZdqG8kMwjV0q6FZYKOG0fDc+5nJjLOYoDoXw/aSa/JqUqk4jv9bD+eSPVXp5JQzMrI9OKmonj5zK4cDyFXFHDUkHnkVEyYMUmYDCiQDMYFElCNCQjJEvoi4YqtOuLzuvDSkFbJXB100TJBJI6197zqonhRAhDibK5rKAZuO+p047Qfub0CobiIeTUssafCMt45vQKPgTgvqdOQ5HIccrHQjLyqo5vPPEatg3GkFN13PvUaedqOEoBY/jLH7+OkwvcYXtyIed8v+dMPgkvFyp9m3/2+GsBr3p1JOJNcGwN2xbWSeu9uEv7TrjeswW3/Rm5ScFcD7ffASQ13EkxaObu+QD+O4CLwPvnAgAYYy0tNK5IkqNFn17MV8R+29hCJx5WMJIMYz6rgiQ4TkMwYLy/HGbl1mozJR3bBqKYz6pQrW5P9gNrMGaZS0JYzJdNAKmogtlMCZql5Xqbhwellg2/mfDEM2l+07tTKxgrv+93/HoOZFU3eeyywVatoHims2FVDmW4eGs/puayYOAP+cQoNzH95AwX5JphrhLAbqFssvK27vNwTxR23Prz1nZRS2hoZlmjpwbdWxFF4lpnLFRxfR+fWkQqIiNd1J1sWcBeTfLrYTJgz0h5heRoiX72SM/LehMSY8CW/ghm0kWAAaOpKEzGJxH7WLIVKQUAOrNWGLppZbejYoWaU3UMsTAArnVmijru+bB/Ace/eXraEvqS63xMZEsG/t/H9uPg/34CfdGyqGCMIWzFsZ9ayCFfMnBqMYeoIiEVlWC67PmT59L47Leex4tvrICIrL9ZSh6A6eUiPvT1JyvG43Xo6ibDI6/MVb94PgzGQxWCu0KLrhDclSYU+zPRkARP3lLXYfsdvvWxUz9p9LNBZ4o/BfBfAHwJwHsB/BraUOBtYjThlPt9828/AEWyiraxVc8UAN4PNmuFbuYMA9GQ7MRh+6XV2+YOW0i5tXn7vbyqYzgRwWAigumlPAbjYeRKOjSTgZnwbR6+VoI639zmoGpCpJZwqedANpkJn7l21b5tAePn3C3vq7GxkWcbt6IkWSYfOzfATaOOQLLG672+r5xLI6caCEl89ajbjkhwpUCWCAORED5z/VuczyiytfpgrhNwzoN5Xtce11hfxMoI5ZOM7rmmVdwP0AwTssT3PzGSxOmlPHTDrEg2csec2w5Ix9Go6o5pzD3pErgJ6re/+RNkizoWc6o1ubCK7/bf/tmR2icG4EcnFqqMvjruyylJhOvfOo7Dk7OIhSRIVtLXYl51JsIdg3FIEkHVDYwko/ijf31p4GNtRoIK/hhj7HtERIyx1wB8noieBvC5dg0sYTmkIpaJxdY23QJhPluCajBsH4xip6f+h5/9PWj9lc9+4KKq0TBA7ebhzRDE/u41B52xnG12UTf7mYqHK+PbG0EiCRLxh99vsi03Q+ev/Zy7bpOHF+/f3K+5YC+bRoYTYWe7kCxBZgyKRM5KrZogrCdgGePj9l5fb8ZxSCZwNzRhvC+K7YNx/Pr+3XjXnhEUVAMGY9g1FMerczmYnoEoxO/XWtfCTTIkOQX5JKlSswjL5Kx8TVYWzIrEtXrJ8sPIRFjIqdZrvs30Ut7xWy0XNPzsnT8M7IC0h/CEK1a/FvawJaIKu/NFW/vwpuE4Vgoajry2BEUiRBTJiu4BPvSON+Edbx5CIqLg746cxl/8+HXHT2CvYn7lHTvxoXftwvRioSKRSrbyfGSZEFIIRY0rLp2Mr+8Vggr+EhFJAI4T0b8HcAZAW8sDfnT/btx5+AR003RuBDBew8YW0kt5DYmwjJmVYkVki12SwUsz9VfWo/RBEPu71xw0lorgnOVodGvJv7GGNn9hRYKs8kQ3IkDVuK0d4ELH+x2kogrmsqqTXFTQjJqSP6yQYztXJELIntQ1AyYA3eQTl0zcxm3vMxXlOQj9sRBiIRmvnMtUmAGDCliAhx72WSU4PrJ/F4qaAdOKU2cl7l+whQ4xIBVT8I2PvsPRLN3f+0f2T/CQRFWHaTJIEiFEhC19UWSKOsZSYcxm1ApB5l1RKQTIioS5bIlr6daEazvU3eHL7sWNbsKKuLK/H56V7ibv8tFkitVzfKpdv1hIwgcu2YpEWMFCVsWx6WVkihqGExFcd/E4rpwYsmzeMsKKhKdOLtWMdbcdszPpArYPxFf9/VfevRtEhL95eto3quemK3fgzsPHUdAMRC3Nvy8WwmAshExRb2t8/UYjqOC/DUAcPJz4v4Kbez7crkEBwK3vvwAn57O4/7kZaFac794tcZxdKTlRPbJkxeJbzkjdZJjPqNCMTNX9es0dj0zO4vGpxaqCYz3KyQaxv3snoLG+KBhjmMuqkCTybefYKOePpSpCRGNhGQBDUWfQTV4m+n17RzGTVp3koqsneH2b2Qz/TsKWQ9eO/rCL7gFlhylQTlSxS1D3RRQ89/mfAbC6faHtT7HfIyLEFAbbp8fABajuWo3EQhIY4wl/dtSNREAyqmBLXww37duBidEkzi4XAAA7BxM4Q5XhiP0xBVv741jIliqiQcq/G9j3piEcO72MTElDRJGxpS8CRZZgAsiVjEpzkA86g5M1XEmwqUwmQixcDnuNhWRMjCSwYziORJhHhdSzcf/y134MTTexlNccq9VgPISwIuHfHdgTaBxA/Vj3ILHwH3rXrqrhm36JVL95YE9PCXoiclUFIOs9OD+9f7efD/t3yf13olWmz6DUFPxE9I+MsZ9mjD1FRL/FGPvv4Pb9tpAp6jh49xMVyS67huNOC7rjc3mMpcLYOcQbbUzN5yABUEiyxguYxCoKqtWiWjTNjdPLTv2ZZFhGusEiYmuh1gTkF7f/zonhwNFF9bBXHuP9StX2eLZQZgCW8yqOz2Ywmoo4prapuRwkCY6t3G1rjrg0/nLUBv+pEPDuL3wP2wfj+Mi7d+Prv3oldNN06pH/8JU5qLpp+VpMFDxfsT2n2A7+nGpiJKFgS18EBc2EZjB88NLzsH0ohpzljHzxjTSyljDXTBMrBZ6Sb0fkFDQTc9kV/MKfPB7o+hU0E8vWvVKPsMKjXtwRIfFwOWrEHfqXdDkd7WiRZKR1Dki7Fs1IsjJqZzgRqfGpzhBk8gAqBahEZQFrXy9bQXCEK39R8XqVEJaqvO8S3N59kkdQdwv1NH53Qfh/BR7Z0zbOLhew3ar1cWI2C91kSIQVULhKli0DTPDluztqIiwHu8B+0TRzVh2X7YMxp4gYA7cLB43bBxqLyw8Szumt016tUNZaaDTs1P0dIcQnxrDCk7RkmY/RHYdPIBCYr50/qxkYTkVwdrmAz93/PG7+qTdj79YUciUdT51awt8+fdrJt1BrzOvuOX8+p2M+V9amv/5Pp1pynWwkArYPxjHeF7HC/fwjROKR1QI9JLc8NqJpvCaUomb61qKpqo36vQ/4CtIKgesSpEClMAYqI9ZWCVnX58i9TRcJ126mnuBvOoCFiKIAHgUQsY5ziDH2X2p9RjVMvDrnJAYjLBPmsyX0xUJQDbvlYjnLNiTxB12RydGCQ7IE1WCBmjv42e8zRR26WVnaAAAG4mE8+PFrAJSLYtVK0291F6rHp3jVUbszk11ieq112e1ENR6twXD0tUU8f2YZOdXAcl7FkycXcNnOARgmw1cPn0BB1TFX0qEZPDmLJzUVEVG46WsoHsJMmlffDMmE00vcDi0RrxPkvqGSEZlHlDBeGG1qPm/lcgC3/8OLazgrf2TLJFYtzC/hTqbxZErefv8LmF4uOhUl7QgXhYAv/MtLWjbGWtqkW9B5TQWOuPOYBVZtT2WBKlmC+IOXb8NQIoyv/egkzizlsW0wjpt/it/T3aitCtZOPcE/QUT3g9879u8OdTpwlQBcyxjLElEIwGNE9CBj7Imgg+OJO1yFkwCUTD4Qu8GHwQBJAsb7o445aC6rYiwVDiRw/ez3Jd10Iols6nW2WmtRrFoOZNMqVPX6Yg5DiTCGkxHHAcis91cKmvUecwSSHfHjFey27d1+bWMyhj/94Un81VOnHWdkTjXwJz+YwjOnl/G27QN4/uwyVL1SFzDBTRyvzmUr+iXPePIJDHtALrKlylDMWtmTBC5wJSpH4ABYlehnR5QwxgACvvbhK5G0hHtECW4a4UKxrNXOZEqQJU+su8kwvVzAUCJcziegSkFsC3L++2q7rURlTbiTwvWn3zqOn37reMeOL1hf6gn+D7p+/2IjO2Y8kDlrvQxZ/xpeQdjlE+wPKnL5qSIibO2POgXB8qqBsVTYaWJRT+D6RdPIEqE/XimEC5qBbQMxGKarG09IBgP/aTId//ORV7Fv9xBMxvDaQg59sVBFLLZdsGo2XYRbBpqMYSwVwemlPLJFfVWdk1MLfAU0mozi5HwW6aLuOOD6rN6r51YKlc5Hv3oiVd7PO8WiyrWPuJAun//jU4sVJXj9qBX/70dIIsQjMjJFvRySCi64Y2EJW1IxfO7nLkIiouDz33oBi3nVCeNbyJWwkNOs7QG7QgM3BZWjZyaGE3jbtn5HiLuFsiNsXdoxT2IiXwHsaNuuv3HnGmEgHl61vUDQzdQU/IyxH6xl51ZJ56cB7AHwx4yxH/tsczOAmwFA7lvdY3bXUNypxDia4CnhtqljS3+Ea6q/diUYA973R4+gPxqq0BzDsoTXFnJYyqmOlmyXab7ovD584n3n4y9//BreWClga38M771gFA88P4OVgoqIwmv76CbDz1+2Da8t5MrdeHR3TDXh9UUu1AFgLBVd1bihoBkYS0WRLZVtzozxzMvzRxN45nVe5wTEVx1FzcSuYcKf/egUcqqOc5kiVlwheQzASlHHc2eW8dNf/mET305jRJRyNUt7+c/Pgb9OxRQMJyL4mYu2YN/uISeaJBaS8Y0nXsPfHplGXjMQD8t495uHMZdVMbNSRCKs4OwKr73PtXeGTNHAL79jDFe8aRASEf79tXtw+7dfhG7yxKuBeBi6yZ34qsEQD0tgpgkGgm6avExBmJeM+OAf/6gl/VwnRhI4PpsFMVeDGQacP5pY+8UVCNYZCtIjkoh+gtXa+gqAIwB+jzG2sPpTFZ8fAPBNAB9jjD1fbbvI1vPZ1g9/ufw5AIc/eQCMMfyHvz7mK0ztGiRvpAvIFnXEwxIG46ujExrJ5HPHG3tjg+1xRBXJMaHkVQOpaAi/9u5dvM76mTS+NzkLWMLM7htwXn8MkoQK7bxNtaEQC8kuB6Pbrq2Uo0c8RaN+55vP+9aXIQAXbEnh9cUc/zvjmrXtZ9g5lMDXf+3KCvNIPQ3a5uDdT+DkfHaV72L3SLIiWskb4lmt/PH0Uh4JV5evatFJjfLI5Cw+deiY4wOyS3f8QZf1iBVsPojoacaYfz2OKgSN438QvMrxX1mvbwKP658B8GcAfq7Whxljy0T0fQDXAagq+L28b++ok/p+05U7cMd3JnEuXXTqloQV7szVDBN9UQWGYTomgIF4uGp0gqqbPk0OyvWxX57J4MxyAdmShqJm4us/Oom/eOIUciXD6uajrxKO5zIl/Oe/r31qJxdyNf9uI7kccxed14dEWMHjUwuWaQKAZTe2Vy9f/teXVRSNsotDuSsZhhUZH3jbVly9Z9hqoUerhPMjk+fw7Z+cWzWeVIQn6Gzpi+J1qzUeEcEEN6n85nv3VDSst4XwK+fS0AyGsCJVNBvxq5U/OlrOB0wX1Iqywn69lb14ew9oJmtpo/kDe8fwBzde2lCzHPe1aKQ/czfQq+MWBCOo4H8/Y+wK1+ufENFRxtgVRPTLfh8golHwJu3LRBQD8M8A3BF0YANRGe9/S6Wzidv6rThwMORLBpIRHivOzRA8HXwpryNTMhCRZQwnebTCV79/wuliozVgkG4kNjuiSJ5EGRlxTwRJZeQI//2r3zuBdFFDPCw79mfvSuXnvvoYSrpRdi6CO0MTYRnv2jMCiYDHX13Alx4+hTPLBSTDMhYtrXckEUamqOGr3z+B0VSk6gP87OmVKteAG9FtW7ztNAZbXSLNdn6ruuHbbOTG6eWKXrfemv/pgobTiwWAgHPpIuazJdx671GEQ7JT176e0/70Uh4yVdZiaqZ71lcefmVVT+NGciZa1RR+venVcQuCE1Twy0R0FWPsSQAgoisB2DaXarngWwH8uWXnlwD8DWPs27UOIhEPtyMCijrDZ+9/HkOJMBjjdXkqzCKWQy9dMpAurS7UpZu8fO1Ksb7QJsDRmNMFzWnUzKyoobAsoT8Wxr/at91jOqkU4kqDsdkS8YzjX/+p3fjDh1+BZhUPK1rmlFveM4Gtlnno1/fvxle+fwLEmBPJwkC4+acmMJqK4JHJWfx/33kZIZkwFA/jxFwWusGQjIYgSRLiYalux61pqwSAY5mxonQ0q+76uUwRRIQ3DcUcDd+7TzuiaSGrQwIv/6CbJs6lSyAC7jx8AqPJMPpj3AHvbWxzdjkPE4BCvBAXM4G0ZiCkGdhqVV2tp8Enw7LVp5ecCLAzy0XsacAeb5dplog7kAuagTsPnwCAwNnRzVRdBTqvbTc7bkHvEFTwfxTA14koCS4n0wA+QkQJVEnqYow9B+DyRgZjMlbR4QiA1UC9NjIRJKlsHokoEt45MVxRFzvpSlF3NzlIRGREQ7Lz2V/4nz9CtqRDJp5gwpi1otB0/IvLt9Uchy3I7X/2azv2WybCj47P42s/Oonj59LQTThmkJv27aioGXT1xBD+4vHX8PsPvOQ8/B+8dGtFCYur3jSAx6cW8Td3HEa6oCERkR2BWq+3LbBas5te4uUL/Nw+KwUNjPE4/blMCWeWCwjLfEKZXio4ZplXzqWRjChOlBAZZfNMROaZuws57jjvi4VWNbbRGRe0doITEQCjMjELWN1lyi0c3f0B4JrEGgmXvOexk5bQl6y2nbz0w52HT+CS7QOBBGCQ7lh+vgq3L2E+U8KnDh1bV1/CetSnEnSWoI1YngLwNiLqt167bQJ/08oB2bXFAZ7gAwC/cWAP/vbIaeRKdps8LujTBRUrRQPjfdGKjEN3V59G8VZpJOJCVDMZklEFiiQ5E40j5C0BX0+wPDI5i9974CWouoGMtUqxzSBnlguO89FvqW036LBLWCzkSnjytWWMJsMYSUYws1JEQTUQUWSkoqG6fVKB1Zodwd+GHpGpopGNfd4l3UDOavhsm2VU3cRSrrzKsvdHAMAIEvF37MQ8oLKxzYW/8+Cqksbu/QBAuqDhzHIRiuRvivD2XQjLEsb7IhURVfXIqQYUq8+Du+yEYTJ85C+O4IKxJD593d6awtibJ1Jv3ABwx0OTWMprkCWCIvN6Q0t5DXc8NNlWwe9eZaQLGnTDdFpoAq2vTyXoLIFsE0TUT0R/BOB7AL5HRH9oTwKthMA1vbDVLKNkmHjzaArXv3Ucv3HNm/mDACCk8MJXIUXBe84fwWJOxatzOSzmVFx30ZaaQl+2qkJGQ1zrT0VDGIjzxKgt1gQCy4btrr0RDckYS0UxlAijP8411XhYQUThJp4g2qQtaDNFbgZRJAmSREgXdIRkwl2PTlVsFw8rsJvOZEs6MkXdeS9d4CUsMlY/gojC89/nrIqdI8kITMZXQzzjWV9VY8iu3T41l+XN6auMezTFo6QcgWzF+dtaOAMcs4xpldHw7ouB+2VGEmHAKp7mN66JkYQTLcXATW4y8TLJeVV3TE4ANxPZ18d9/XYMxqFYPRb2jvdhYjQJRZYaElyJsAyTlWvUu89HJuDkfA6fu/8Fp1WiH7dcMwHNYBXjNi175cvnMryqrG444wbA60+RXYyr7ICfmg8WGFANO9t8/x2HcfDuJyrGbSsas1a5lERExlxWxVymWPXeWesxBZ0lqKnn6+DROL9ovf4QeHOWX2j1gMrNNhh0o1xb268y3+U7+vHQi+cwnAwjqvAOUf/40iwuf9Mg3nPBGCQJlnDlP4O0QbtgS9/qYmiJUNUm5o1gL6HtfrYAX1F4G4P4LbUNq9n4C2dXHF+HhHKjeMm6dkXdALPKDA/EQxhNRqrWGEqGZUyey6Ie08tF7PrMPwAAYgogSdKq5itB6ryf1x+zOpyp0E3g+bNpAMD2/gjuenQKv/Ot55GKKAjLhIJmOolqEgFRRcJcpoSwzMslSGA4tVA2PcQVcorz3XLNBP7dXx2t6L0bD8uruqZ5benjfbzSKNf2yTcIQJG438IwGUIy4ZN/+yyWCroTaXbDJeP40k08DsJb+0g3TMd8aPseFnIqdCNd87qZjIcD1zIP1fIJ+DWF/+ShY/iiZT7yrvzsBMhcyWioPpX32goHcfcSVPC/mTH2L12vf5eInm35YGTJeuBMSESYGI7jwN6yAL/+bVvxgUu3OsL8w197EvGw7NywEauP571PnsYNl9W2x1ejWoXKVlTjtJf+YVmCbpQTgbyNQfxKSTiNUVjZJGOLWtmqHcMMhmhYdh5Wd0MZP5qx2RZ0oC/C2yK6WycGIRVV8MrMyip7/fRKCXM5FeePpbCQK6Gom7z7Ggim1bBjKBGGIkvQDN5BquTZR15nkPN8tfOtZ6dXNVzPqwa+9ex01YJzL72xgsendMgSF+7VegVLRM53Np8uYsUVWGCYDN989g0ARyuEv33MSz7/Hd5P2mVGNE1W0Vtg93AcJ+ZyIKvwIO8EBoQk+ApQPwH7qUPHMJwII6sa2DEYx+nFHJbzmuPsZiawnNfwhQdfwoG9Y76KxnAiAkXSnI54jSIcxN1NUMFfIKL9jLHHAICI3g2g0OrBEPGlflHnZXRvuPQ8fOKvn62qyUwvF+o6z66eGHJKLAeJkAhSG79Z7EklFVWwkFP5sp8BfYlQxeTiV0rCbsjha46xTC+SRNg5GMNDn3hPoPFka5W59EDWbMPAI6mCJP55WSlo8AnAAsBr9dsmLFkip6uXbl2jmZUiFJlnD+um/7Ez1mzwrWNv+P79W8fewJdu4r97BZMdeso7jElgjA/UdvBqVo8BzTQRkiSMpiLOisNt5WMMuP+5Gec4bsKKhIJqOKtaZn1vYaVscf3M9W+p0M5NxiAB2DYYd8xabgHqPQ/DZFjKa8iUdOwZTWI2U8T0chGKxFdq9niZyXDSGn+zPSdqrTSEg7i7CRp/+BsA/piIThHRKQD/A8AtrR5MSJKQLuoYS0Vx4xXbcOjoGcfuaGs6bjvhjsG4E2MOlJ1nBK4dnVrI4s7DJ3ByPlt1H34c2DuGe29+J3746Wtx783vbJmGcmDvGG6/4WLsHkmiP8pXFP1xbkZyZ5Ue2DuGG6/YhrkMbzozlymBwKNdnCQuF4Zl2tk2EF0VFdVSXAdeCZjb4CaI9uju46savBuWyUyUDG7qChIxWy0b2v3+6aV8RRa4tz+Au7/weQNRhJ3ib8DW/khNs2G1YnPnj6UwkgpDsUxFikQYSYVx/ljK2ebA3jF88cZLcfmOQYz3RRGSJexwhc8CWGUWdJ/HXKbEexSbzJkogNVRUW68voggNn2vX8D7bHmfTUA4iLuJoFE9xwBcSkR91us0EX0cwHOtHIy72frBu5+ou1T0asZep5/bATqainbFcrNew3OAP1SHjp6paHCSKem+D69E3OwC8Jj6MVckRj3br7eypRtHI/VBlniUz5t/+4GaFTWbwdvH1zZ1EPiKxrTtXWvEq+Xa18JdBtnuL5yKhpAaD2EuU0Re5b2fx1JRnJzP+V6/apNCUDOiNwt5NlOs+Hsts6Ddk9jdB9mus+ReabjrDDWzyq1nygnSTlTQOYKaegBwge96+R8AfLmVg3F34LJT+d14l4oH9o7hxullJ7vSZAwjiZCjHbm1x2r7aAVBEm4aScrxe6jCEnzNJCEJFY3i7QcriHPtHbsG8fjJJd8xeIW++/UNl/CM6nhIcsJS6zHsqXjqf0yGvpiC2YyKVFRBRJFwZplXMw1JcJLqwhJ8m7FcvXvQOpaChfzq0M3hePl29wqmvqiC5YLurDBs05a7x3NYkfGFX7jEuX6fuO8ovvnsG6uulX19vDQjYOsJUO/fZeItSN3dtAbjISzkVBC4KVSRJAzGQ/j0dXsrxtaIMlTPlNNOk6lg7TQk+D20vHi4uwPXfJan8hORI8i9S0WvZnxiNoulvI54WENfLLRKe/Tbx1oJImAbjXDwe6g0H0FHAHRGvpEXgZxrJCEekiqackuodBpvTYXxRkb1jVrZPhjH5EzGt6ib+72YQvjDX7wMABekaZ/G37LEz2PXcBIHrxxyktnOH0vitYU8SoaBsEQYTUWQiobwysxKxUR49e5B3HvLuwAAf/iLl+M3vnEERVfvgKhC+MNfLOcTegXTW7b2V0T1ePsL+wkufh2O4v7nZnyvjx+NCth6AtT7990jCcxlS7x5vGvC+th791QkCK5VCAfxCzR6roL1Yy2Cv+V1JYnKHa+2WDX2X1vIW5EWvBqiOyTPK9y86f9u7dFPK24FQQRskG28CTS5kgZVZ87EZVrRPFGXPbeoGTBMhtNLBZxdKWK8L1zhXPPWq/HLsh1NRSoSnRJhGStFHaOpSNXQQHtVtpAtuRNjndDLgZgCzYRT4+a8/gg+8hdHapqFzuuP+sb+D8TDuP6t4059n5gVuZWKhTGRjCBT0h0nvrsr2vVvHXcJcV5nxyuEWiGYvnTTFb6O3FZSb5zPTS/jhbMrzrlWm7DW0qnNyy3XTOCTh47hzHLBmfSSEWVVyGyjrLVcRSvKXXS6ZMZ6UK/Zegb+Ap4AxFo9GIkqFxFW0ANPjqLVA/Fqxt70f6/22I7lZpDohXrbeFcE+ZKOxbwGyUpccpuqSrpRDu10YYcS/r/nHqjoO2CbLDSjnGU7s1LAfKYEzeRZtiGF5zioholc1kC4ysrEO84zSwVH4DvjYMBSQcfukQR2hmScWshh8lzt5CMCT5QaTkRwaiGLJ08tOhnJs5kiDh09gxuv2OZ8jwmrmJ1qmBiIhXBynn9mLBV27UPFaDKMnUPcyXjo6JnApRZ6Cb+aQvc/N4Pbrt2DW98fvKBcM1jpNlUL9jXKWmP/W5E7sFnyD+o1YknV+nurcbcCnM+WIEuEmJWBCawuCOa33HSn/9u0UtPxEmTJW2+bux6dgmYYWMjypDHDZE64lWnFjBMzoJr+Dle3I9Yb6rjKVg84pQCc+cT+rBVPrhkML5/LOLXxq4UO2sclAqIKX4kUNe78rOb4dBOSeXQLAUgXdIwko1Ud8g/85A0MJiJgAM6uFCvyNzJF/pla+5jLFHHrfc+gLxbaUFqcu6YQwFdcumninsdOBi4m1wx3PTqFvlgI4/1l/W+tgRN+K+P5bPDvrRW5A5sl/6CxcpJthjE4IWUl3QRYuVwAsFqTbiYMrdUEGUO9bY7PZjCbLiGvclOUXfaAiJySAyFFrjKC+ng1sZJmVtSf0UyTVwS1XtslGOzM0uPnuE/fGzpol6mwO5rZvX6BcmeuWsjEo3QMq4QD4O+Q1w0Tx+eyTuhgTtWxkFORtkJKvZ/xvk4XNCzkVORUvaGw3l4oOZBTDXiDiCRCe8N6sfpeANYeOOHdZ6aoYT6jIq8agb63VoypHefVjXSV4D9vIIaxVBQrBV6bfiQVrohf9nMe3X7Dxc5nxlLRNXVZaoYgY6i3TV41nJ61bmuXewVUsMw0EVlCVGnsa4t6bmSGyjBOXruocp8l3YRmmhWZpd7Y7KgiQbZyC+y4dAd3Zcwq2H1vgfK52v4Mt0P+XKaEkCQ5dYrs1cV8tuT7Ge9re7uoIvvW9vGjXpx6O2lkwrFrCrkxGX+/nbQjTt+7z7lMCSAejhrke2vFmDZL/sFanLstJxVVHBON/eDViwPuhsiBIGOotY2qW41OnP84jKGiLIIslauGylRucF4vidZ7I6+iyucZK0f4AHzl8qlDx3BmqQDd5PHihglePwc8kcz92XoUNMOZG0wGTM6krZBDhuFEpUM+GZadOkVEADGgxHiht1RUwVxWRV9MWRUSaq8eCatXj8dnMxUO4YajotB6R2CjNuaP7t+NOw+fgG6aTi6Cyfj77aQdcfrefRZ1AxJRRWhqLe27FWPaLPkHXaXxT85kHA2nG7T59UKWJKfZeC0l2TCBom5wO3oLj6/IVCG0bYiwKlOWAQCVzTxOVI/H41zvXCr259qHRISYVW/IzlqOhySsFPWy+ciakOwQ0N0jSdx27R7sGk46Tv3brt2D3SPJqqvHhVwJmaJeVaMPsuS3a+c/8/oSZlYKeOb1JXzq0LE1rQr8KrPW0nJvff8FuO3aPYiFZOgmHyN37LbPvg+0Z7Xt3WcirGA4EXbKdwO1te9WjGmzyJ2u0vgViVZpOBvtgtu4NUXJpb17sU0/3lICrQymtZ3n7sqfzD6W5Xexx6tY/QcMMBhWZmvI5YA/s5THolVPXiJA1VmgoUoEbBuIgwg4vZgHs95TddOp/OndT0k3cXalyMtHMBMvvJFBTuUVJa+eGKq5elzMaRiMh6pq9EGc9n618xeyKv7dXx3FUCK8quKnO9Sy2gqhmRo3l2wfwMXn9Tv7BFB1JdMotVY07Xg+3fsMuuqv9vlWjGGjQs0U22oXAzv3sp/61D1O6YFG+pt2O+4HKBmWsWD1w42FZEzNZSuSqHqZfitBa73uqojVm8Ft57YTyLYPRAEiXwHMGMNQIlyRw6AQkFa5eSGiSDBNEyC+f2ZlDe8dTzkNWOymMXYfZN0woVkDeet5fTiznMeSTwaxBDgVZ1NRpaK7ll2iwT3hzGeLyJUM38gWt2nIbtBzLl2CTACsBkHJiOKUYG4E775twbueGrD93FQLx25nzP1a9r2euQBE9DRjbF9Dn+lGwc8Yw0qh+ZKw3Yb3AbL74W4f5MW3Gi1vLChjC6RqhGQeX64ZDFv6IhhJRlDQDLy+mIdpMquJDqBqpmM+c08m7npGsgSMJiMIKzJuv+Fi3PKNp6G6oqFsCMBbt/VXrKC8tY+iITukluH8sSQe/Pg1APwF+WxGdfIavMLXO1EcP5dB0fJpRBTrGIxhz2gicNVWG79JqJuUsnZOTGvZ93pPmM0I/q6y8dtsBC+6OzLj1vuegaobjt2W93DljUSE0G8t3kZosq2xA1jKqeWKlYybquwsQVvo8+zj8mPhNL2xHMqZYrlbWjIsVc1udH/WD9PkvZx1k2FyJuP4BbyVWWczKvqiMkZT/t3GvL6IkquSHxFvGiMRnBLMjbDW0MZ2h8M26g9Zr323c1ytousEfydi8VuNNxTQG3fOzGB2b0HjeBew7paY7oYn9gxh5zB48QuR1BmPj399MY+jry9V72cQwKutGmbFMWzHsrv+1FvGUyAA2aKBTLFcBtstfP3CDwE+nxU1AyXdaLqC6lpCG9cjHLadMfdr2Xcv5AJ0leA3TLYhvOjeGd8bd+5+lAK06hVYEK2+Xiar7xvxvcbWDBGS6udFeD9upyiUdLMiGsut6TOrV7DncL4oEhyN0HvveHspA5XC15sc6E6lsPM1DAaMJsM1z9GPtSRIrofW286Y+7XsuxdyAdom+IloBxF9n4heJKIXiOi2ep+5cDzV0sYnncI744+mIhUNxm0hYIcmCoLB2GoBWtJrX0AG5tTHD8nlxvMgshyg/J/kbF85mfiX1ieM90edaKyIIiMa4v/s7V+aySASkhGSqu2jjETkaITee8eOYbd7KXuFrzf8MBKSKyYhAs/5cIeyBmUtoY3rofW2M3N/LfvuhooC9WhnOKcO4D8yxo4SUQrA00T0XcbYi208ZlfgDQVMRUMYSRlO82ouiBjCcmWlzY0wB/RHFd42cA0noxA3q9hs74/gyt1DTvljiYCITOiPh7GQU32botvY9efjMRlb+6JOCWseJy47UT2xsOxE9egmzxGIWQ5Yd7tHiXhHrlQ0hOFECHNZzUme0i2H8GBcwbYBrvWlCxqGE2HkVAOzmRL6ojIWcprjOJaJ+xfcGqH73umLhVDSDeTV6o3P3eGH++84DDmJimilkWQY2dLq6KIgNBvaGCQcdq2RL+2s+b+WffdCL4J1i+ohom8B+B+Mse9W22bfvn3syJEj6zKedlLPq++uqOjOtrQTby767EPIu7JagfLCQCJX43WLkEw1hV8tJPCV1iuzWd7YPMBnIoqEl3/v+ppRH3YG5GymCNNqA8gYN2v4ha7+i8u21qxj78Y+7sxKEbrBUDIq92dr2ESE8/qjvg9evYiVepFY9vYlzUC6qDu1cfqiMrYPJnz3aR9TNxjOrhQggcDAk9bG+rg2DWBNESHdEolT7xnohlDRjULXRvUQ0S4AlwP48Xocr9PUWyLXy7b8jfdMcAGP8j+Af1m2WcERbuCtF3cNx3l8OMFJngLAs1aTYbxlfHWhVQLw9V+9Eg9+/Bp87Vf2Ye94ypmMJPjfHATgNw+8GUDtJa19DXYNxWEwXoFza38EWwdiiIckZ3yyRA0JfaBsRrB78ioeWwov60C47do9Vfsm11uOe7/DXUNxDMZDkCWq2P73fv5teO7zP4NX/9vPYmt/FNsGKu24bvOGfUxF5hMSSdxUtHsk4dwfa80c7RYzQ73z6IXIl41M2zV+IkoC+AGA32eM/Z3P328GcDMA7Ny58+2vvfZaW8fTK3zl4VeclpKJMLcVrxQNyEQgsnvRMoQk4MLxvlVmhVrJLrWWn95tvIlPH92/u6IcQDP7XOuy16vxSxJB1cuF7mSJApUtaHRc9bYPom23+lq04rw6wf47DmMgFqqIutpo+TvrRdclcBFRCMC3AXyHMfZH9bbfKKaedvDI5Cw+eegYt59bHY/CEuG8gVhVQb9Rsc0Eqm5gIafyNxkwkgojJMsdMxcI80VwusUktRFoRvC3zblLfCr/GoCXggh9QW0O7B3DF2+8tOs1ufXA7TzTjTRUgyGsSNg1nOzoNekFp163sFmqYHYrbdP4iWg/gB8C+AnKiZG/zRh7oNpnhMYvqMZm6IO62egFk1Qv0HWmnkYRgl/ghzChCATV6dqoHoFgLYgIEIGgtQjBL+h6eqH2iUDQSwjBL+h6eqH2iUDQSwjBL+h6uiUpSSDYKHRV60XB2tiokS8iTFIgaC1C8G8Q3JEv7vrndu/iXmcz9EEVCNYLYerZIIjIF4FAEBQh+DcIIvJFIBAERQj+DYKIfBEIBEERgn+DICJfBAJBUITg3yCstY67QCDYPIiong2EiHwRCARBEBq/QCAQbDKE4BcIBIJNhhD8AoFAsMkQgl8gEAg2GULwCwQCwSZDCH6BQCDYZAjBLxAIBJsMIfgFAoFgkyEEv0AgEGwyhOAXCASCTYYQ/AKBQLDJaJvgJ6KvE9EsET3frmMIBAKBoHHaqfH/GYDr2rh/gUAgEDRB2wQ/Y+xRAIvt2r9AIBAImqPjNn4iupmIjhDRkbm5uU4PRyAQCDY8HRf8jLG7GWP7GGP7RkdHOz0cgUAg2PB0XPALBAKBYH0Rgl8gEAg2Ge0M57wXwOMALiSiaSL6SLuOJRAIBILgtK3nLmPsYLv2LRAIBILmEaYegUAg2GQIwS8QCASbDCH4BQKBYJMhBL9AIBBsMoTgFwgEgk2GEPwCgUCwyRCCXyAQCDYZQvALBALBJkMIfoFAINhkCMEvEAgEmwwh+AUCgWCTIQS/QCAQbDKE4BcIBIJNhhD8AoFAsMkQgl8gEAg2GULwCwQCwSZDCH6BQCDYZAjBLxAIBJsMIfgFAoFgkyEEv0AgEGwyhOAXCASCTYYQ/AKBQLDJEIJfIBAINhlKO3dORNcBuBOADOAextgX2nm8zcYjk7O469EpnF7KY8dgHLdcM4EDe8d67hjtIMi4vdtcPTGEx6cWq77ulnPv1e9E0D0QY6w9OyaSAbwC4J8BmAbwFICDjLEXq31m37597MiRI20Zz0bjkclZfO7+FxCSCbGQjIJmQDMYbr/h4pYJgfU4RjsIMm7vNgu5EmYzKkaTYYwkI6ted8u59+p3ImgfRPQ0Y2xfI59pp6nnKgAnGGNTjDEVwH0APtjG420q7np0CiGZEA8rIOI/QzLhrkeneuoY7SDIuL3bpAs6JAIyRd33dbece69+J4Luop2CfxuA067X09Z7FRDRzUR0hIiOzM3NtXE4G4vTS3nEQnLFe7GQjOmlfE8dox0EGbd3G9UwIRH/6ffabx+doFe/E0F30XHnLmPsbsbYPsbYvtHR0U4Pp2fYMRhHQTMq3itoBrYPxnvqGO0gyLi924RlCSbjP/1e++2jE/TqdyLoLtop+M8A2OF6vd16T9ACbrlmAprBkFd1MMZ/agbDLddM9NQx2kGQcXu36YspMBmQiiq+r7vl3Hv1OxF0F+107irgzt33gQv8pwD8EmPshWqfEc7dxrCjO6aX8tje5qiedh6jHQQZt3cbO4qn2utuOfde/U4E7aEZ527bBD8AENHPAvgyeDjn1xljv19reyH4BQKBoDGaEfxtjeNnjD0A4IF2HkMgEAgEjdFx565AIBAI1hch+AUCgWCTIQS/QCAQbDKE4BcIBIJNRlujehqFiOYAvNbpcdRgBMB8pwcRADHO1iLG2VrEOFvLhYyxVCMfaGtUT6Mwxro6dZeIjjQaNtUJxDhbixhnaxHjbC1E1HAMvDD1CAQCwSZDCH6BQCDYZAjB3xh3d3oAARHjbC1inK1FjLO1NDzOrnLuCgQCgaD9CI1fIBAINhlC8AsEAsEmQwj+KhDR14loloied703RETfJaLj1s/BDo9xBxF9n4heJKIXiOi2bhynNaYoET1JRMessf6u9f5uIvoxEZ0gor8monAXjFUmomeI6NvdOkYAIKJTRPQTInrWDunr0u9+gIgOEdEkEb1ERFd32ziJ6ELrOtr/0kT08W4bpzXWT1jP0PNEdK/1bDV0jwrBX50/A3Cd573PAPgeY+x8AN+zXncSHcB/ZIxdBOCdAH6TiC5C940TAEoArmWMXQrgMgDXEdE7AdwB4EuMsT0AlgB8pHNDdLgNwEuu1904Rpv3MsYuc8Wbd+N3fyeAhxhjewFcCn5tu2qcjLGXret4GYC3A8gD+Ca6bJxEtA3ArQD2McbeCl7y/iY0eo8yxsS/Kv8A7ALwvOv1ywC2Wr9vBfByp8foGe+3APyzHhhnHMBRAO8Az4xUrPevBvCdDo9tO/gDfi2AbwOgbhuja6ynAIx43uuq7x5AP4CTsAJJunWcnrH9NIAfdeM4Ue5lPgSegPttAD/T6D0qNP7G2MIYe8P6fQbAlk4Oxg0R7QJwOYAfo0vHaZlQngUwC+C7AF4FsMwY061NpsFv7E7yZQD/CYDdZX0Y3TdGGwbgH4noaSK62Xqv27773QDmAPypZT67h4gS6L5xurkJwL3W7101TsbYGQBfBPA6gDcArAB4Gg3eo0LwNwnjU2tXxMISURLA/wXwccZY2v23bhonY8xgfCm9HcBVAPZ2dkSVENE/BzDLGHu602MJyH7G2BUArgc3813j/mOXfPcKgCsA/Alj7HIAOXjMJV0yTgCAZRu/AcDfev/WDeO0fAwfBJ9QzwOQwGqTdF2E4G+Mc0S0FQCsn7MdHg+IKAQu9P8PY+zvrLe7bpxuGGPLAL4PviQdsPozA3xCONOpcQF4N4AbiOgUgPvAzT13orvG6GBpf2CMzYLbo69C93330wCmGWM/tl4fAp8Ium2cNtcDOMoYO2e97rZxvh/AScbYHGNMA/B34PdtQ/eoEPyNcT+AD1u/fxjcpt4xiIgAfA3AS4yxP3L9qavGCQBENEpEA9bvMXBfxEvgE8CN1mYdHStj7LcYY9sZY7vAl/uHGWP/Bl00RhsiShBRyv4d3C79PLrsu2eMzQA4TUQXWm+9D8CL6LJxujiIspkH6L5xvg7gnUQUt55/+3o2do922pHSrf/Av/w3AGjgWstHwO293wNwHMDDAIY6PMb94EvP5wA8a/372W4bpzXWSwA8Y431eQCfs96fAPAkgBPgy+tIp8dqjesAgG936xitMR2z/r0A4D9b73fjd38ZgCPWd//3AAa7dJwJAAsA+l3vdeM4fxfApPUc/SWASKP3qCjZIBAIBJsMYeoRCASCTYYQ/AKBQLDJEIJfIBAINhlC8AsEAsEmQwh+gUAg2GQIwS8QACCinyciRkRdlU0sELQDIfgFAs5BAI9ZPwWCDY0Q/IJNj1XraD94kt5N1nsSEf1Pq4b8d4noASK60frb24noB1ZxtO/YKf0CQa8gBL9AwItePcQYewXAAhG9HcAvgJflvgjAh8DrCtm1kb4K4EbG2NsBfB3A73di0AJBsyj1NxEINjwHwYuxAbw420HwZ+NvGWMmgBki+r719wsBvBXAd3mpFMjgpT0Egp5BCH7BpoaIhsCrcL6NiBi4IGfg1S59PwLgBcbY1es0RIGg5QhTj2CzcyOAv2SMvYkxtosxtgO8Y9QigH9p2fq3gBdtA3hHplEickw/RHRxJwYuEDSLEPyCzc5BrNbu/y+AcfCqrC8C+AZ4q8gVxpgKPlncQUTHwCuivmvdRisQtABRnVMgqAIRJRljWSIaBi95+27G68sLBD2NsPELBNX5ttU8JgzgvwqhL9goCI1fIBAINhnCxi8QCASbDCH4BQKBYJMhBL9AIBBsMoTgFwgEgk2GEPwCgUCwyfj/AYVU8YPWmr6jAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
@@ -258,9 +268,11 @@
},
{
"cell_type": "code",
- "execution_count": 6,
+ "execution_count": 60,
"id": "7cf9b21a",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -292,32 +304,38 @@
"\n",
"Treat `Pclass` and `Sex` as factors.\n",
"\n",
- "Print an ANOVA table."
+ "Print an ANOVA table for different types of sum of squares.."
]
},
{
"cell_type": "markdown",
"id": "3ad89dec",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A"
]
},
{
"cell_type": "code",
- "execution_count": 7,
+ "execution_count": 63,
"id": "39f58924",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
- "model = smf.ols('LogFare ~ Age * C(Pclass) * C(Sex)', df).fit()"
+ "model = smf.ols('LogFare ~ Age * C(Sex) * C(Pclass)', df).fit()"
]
},
{
"cell_type": "code",
- "execution_count": 8,
- "id": "bc136eb5",
- "metadata": {},
+ "execution_count": 66,
+ "id": "f09c5fea",
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -340,21 +358,140 @@
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
- " <th>sum_sq</th>\n",
" <th>df</th>\n",
+ " <th>sum_sq</th>\n",
+ " <th>mean_sq</th>\n",
" <th>F</th>\n",
" <th>PR(>F)</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
+ " <th>C(Sex)</th>\n",
+ " <td>1.0</td>\n",
+ " <td>46.721986</td>\n",
+ " <td>46.721986</td>\n",
+ " <td>126.001800</td>\n",
+ " <td>5.241061e-27</td>\n",
+ " </tr>\n",
+ " <tr>\n",
" <th>C(Pclass)</th>\n",
- " <td>318.183365</td>\n",
" <td>2.0</td>\n",
- " <td>429.045083</td>\n",
- " <td>1.856872e-122</td>\n",
+ " <td>318.408489</td>\n",
+ " <td>159.204245</td>\n",
+ " <td>429.348645</td>\n",
+ " <td>1.619836e-122</td>\n",
" </tr>\n",
" <tr>\n",
+ " <th>C(Sex):C(Pclass)</th>\n",
+ " <td>2.0</td>\n",
+ " <td>5.990345</td>\n",
+ " <td>2.995172</td>\n",
+ " <td>8.077506</td>\n",
+ " <td>3.402042e-04</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age</th>\n",
+ " <td>1.0</td>\n",
+ " <td>12.274132</td>\n",
+ " <td>12.274132</td>\n",
+ " <td>33.101391</td>\n",
+ " <td>1.306440e-08</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age:C(Sex)</th>\n",
+ " <td>1.0</td>\n",
+ " <td>1.486269</td>\n",
+ " <td>1.486269</td>\n",
+ " <td>4.008232</td>\n",
+ " <td>4.566288e-02</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age:C(Pclass)</th>\n",
+ " <td>2.0</td>\n",
+ " <td>0.296834</td>\n",
+ " <td>0.148417</td>\n",
+ " <td>0.400257</td>\n",
+ " <td>6.703004e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age:C(Sex):C(Pclass)</th>\n",
+ " <td>2.0</td>\n",
+ " <td>1.335653</td>\n",
+ " <td>0.667827</td>\n",
+ " <td>1.801023</td>\n",
+ " <td>1.658921e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Residual</th>\n",
+ " <td>702.0</td>\n",
+ " <td>260.304489</td>\n",
+ " <td>0.370804</td>\n",
+ " <td>NaN</td>\n",
+ " <td>NaN</td>\n",
+ " </tr>\n",
+ " </tbody>\n",
+ "</table>\n",
+ "</div>"
+ ],
+ "text/plain": [
+ " df sum_sq mean_sq F PR(>F)\n",
+ "C(Sex) 1.0 46.721986 46.721986 126.001800 5.241061e-27\n",
+ "C(Pclass) 2.0 318.408489 159.204245 429.348645 1.619836e-122\n",
+ "C(Sex):C(Pclass) 2.0 5.990345 2.995172 8.077506 3.402042e-04\n",
+ "Age 1.0 12.274132 12.274132 33.101391 1.306440e-08\n",
+ "Age:C(Sex) 1.0 1.486269 1.486269 4.008232 4.566288e-02\n",
+ "Age:C(Pclass) 2.0 0.296834 0.148417 0.400257 6.703004e-01\n",
+ "Age:C(Sex):C(Pclass) 2.0 1.335653 0.667827 1.801023 1.658921e-01\n",
+ "Residual 702.0 260.304489 0.370804 NaN NaN"
+ ]
+ },
+ "execution_count": 66,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "sm.stats.anova_lm(model, typ=1)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 67,
+ "id": "cf4e3edc",
+ "metadata": {
+ "hidden": true
+ },
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "<div>\n",
+ "<style scoped>\n",
+ " .dataframe tbody tr th:only-of-type {\n",
+ " vertical-align: middle;\n",
+ " }\n",
+ "\n",
+ " .dataframe tbody tr th {\n",
+ " vertical-align: top;\n",
+ " }\n",
+ "\n",
+ " .dataframe thead th {\n",
+ " text-align: right;\n",
+ " }\n",
+ "</style>\n",
+ "<table border=\"1\" class=\"dataframe\">\n",
+ " <thead>\n",
+ " <tr style=\"text-align: right;\">\n",
+ " <th></th>\n",
+ " <th>sum_sq</th>\n",
+ " <th>df</th>\n",
+ " <th>F</th>\n",
+ " <th>PR(>F)</th>\n",
+ " </tr>\n",
+ " </thead>\n",
+ " <tbody>\n",
+ " <tr>\n",
" <th>C(Sex)</th>\n",
" <td>12.597516</td>\n",
" <td>1.0</td>\n",
@@ -362,7 +499,14 @@
" <td>8.516312e-09</td>\n",
" </tr>\n",
" <tr>\n",
- " <th>C(Pclass):C(Sex)</th>\n",
+ " <th>C(Pclass)</th>\n",
+ " <td>318.183365</td>\n",
+ " <td>2.0</td>\n",
+ " <td>429.045083</td>\n",
+ " <td>1.856872e-122</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>C(Sex):C(Pclass)</th>\n",
" <td>3.489752</td>\n",
" <td>2.0</td>\n",
" <td>4.705654</td>\n",
@@ -376,13 +520,6 @@
" <td>1.306440e-08</td>\n",
" </tr>\n",
" <tr>\n",
- " <th>Age:C(Pclass)</th>\n",
- " <td>0.296834</td>\n",
- " <td>2.0</td>\n",
- " <td>0.400257</td>\n",
- " <td>6.703004e-01</td>\n",
- " </tr>\n",
- " <tr>\n",
" <th>Age:C(Sex)</th>\n",
" <td>1.421046</td>\n",
" <td>1.0</td>\n",
@@ -390,7 +527,14 @@
" <td>5.066866e-02</td>\n",
" </tr>\n",
" <tr>\n",
- " <th>Age:C(Pclass):C(Sex)</th>\n",
+ " <th>Age:C(Pclass)</th>\n",
+ " <td>0.296834</td>\n",
+ " <td>2.0</td>\n",
+ " <td>0.400257</td>\n",
+ " <td>6.703004e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age:C(Sex):C(Pclass)</th>\n",
" <td>1.335653</td>\n",
" <td>2.0</td>\n",
" <td>1.801023</td>\n",
@@ -409,17 +553,17 @@
],
"text/plain": [
" sum_sq df F PR(>F)\n",
- "C(Pclass) 318.183365 2.0 429.045083 1.856872e-122\n",
"C(Sex) 12.597516 1.0 33.973508 8.516312e-09\n",
- "C(Pclass):C(Sex) 3.489752 2.0 4.705654 9.331214e-03\n",
+ "C(Pclass) 318.183365 2.0 429.045083 1.856872e-122\n",
+ "C(Sex):C(Pclass) 3.489752 2.0 4.705654 9.331214e-03\n",
"Age 12.274132 1.0 33.101391 1.306440e-08\n",
- "Age:C(Pclass) 0.296834 2.0 0.400257 6.703004e-01\n",
"Age:C(Sex) 1.421046 1.0 3.832335 5.066866e-02\n",
- "Age:C(Pclass):C(Sex) 1.335653 2.0 1.801023 1.658921e-01\n",
+ "Age:C(Pclass) 0.296834 2.0 0.400257 6.703004e-01\n",
+ "Age:C(Sex):C(Pclass) 1.335653 2.0 1.801023 1.658921e-01\n",
"Residual 260.304489 702.0 NaN NaN"
]
},
- "execution_count": 8,
+ "execution_count": 67,
"metadata": {},
"output_type": "execute_result"
}
@@ -428,10 +572,137 @@
"sm.stats.anova_lm(model, typ=2)"
]
},
+ {
+ "cell_type": "code",
+ "execution_count": 68,
+ "id": "bc136eb5",
+ "metadata": {
+ "hidden": true
+ },
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "<div>\n",
+ "<style scoped>\n",
+ " .dataframe tbody tr th:only-of-type {\n",
+ " vertical-align: middle;\n",
+ " }\n",
+ "\n",
+ " .dataframe tbody tr th {\n",
+ " vertical-align: top;\n",
+ " }\n",
+ "\n",
+ " .dataframe thead th {\n",
+ " text-align: right;\n",
+ " }\n",
+ "</style>\n",
+ "<table border=\"1\" class=\"dataframe\">\n",
+ " <thead>\n",
+ " <tr style=\"text-align: right;\">\n",
+ " <th></th>\n",
+ " <th>sum_sq</th>\n",
+ " <th>df</th>\n",
+ " <th>F</th>\n",
+ " <th>PR(>F)</th>\n",
+ " </tr>\n",
+ " </thead>\n",
+ " <tbody>\n",
+ " <tr>\n",
+ " <th>Intercept</th>\n",
+ " <td>252.996076</td>\n",
+ " <td>1.0</td>\n",
+ " <td>682.290368</td>\n",
+ " <td>1.334116e-105</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>C(Sex)</th>\n",
+ " <td>0.819888</td>\n",
+ " <td>1.0</td>\n",
+ " <td>2.211108</td>\n",
+ " <td>1.374692e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>C(Pclass)</th>\n",
+ " <td>31.940629</td>\n",
+ " <td>2.0</td>\n",
+ " <td>43.069411</td>\n",
+ " <td>2.273902e-18</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>C(Sex):C(Pclass)</th>\n",
+ " <td>1.012273</td>\n",
+ " <td>2.0</td>\n",
+ " <td>1.364970</td>\n",
+ " <td>2.560654e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age</th>\n",
+ " <td>0.776183</td>\n",
+ " <td>1.0</td>\n",
+ " <td>2.093243</td>\n",
+ " <td>1.483980e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age:C(Sex)</th>\n",
+ " <td>0.179743</td>\n",
+ " <td>1.0</td>\n",
+ " <td>0.484739</td>\n",
+ " <td>4.865140e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age:C(Pclass)</th>\n",
+ " <td>0.437814</td>\n",
+ " <td>2.0</td>\n",
+ " <td>0.590357</td>\n",
+ " <td>5.544041e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age:C(Sex):C(Pclass)</th>\n",
+ " <td>1.335653</td>\n",
+ " <td>2.0</td>\n",
+ " <td>1.801023</td>\n",
+ " <td>1.658921e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Residual</th>\n",
+ " <td>260.304489</td>\n",
+ " <td>702.0</td>\n",
+ " <td>NaN</td>\n",
+ " <td>NaN</td>\n",
+ " </tr>\n",
+ " </tbody>\n",
+ "</table>\n",
+ "</div>"
+ ],
+ "text/plain": [
+ " sum_sq df F PR(>F)\n",
+ "Intercept 252.996076 1.0 682.290368 1.334116e-105\n",
+ "C(Sex) 0.819888 1.0 2.211108 1.374692e-01\n",
+ "C(Pclass) 31.940629 2.0 43.069411 2.273902e-18\n",
+ "C(Sex):C(Pclass) 1.012273 2.0 1.364970 2.560654e-01\n",
+ "Age 0.776183 1.0 2.093243 1.483980e-01\n",
+ "Age:C(Sex) 0.179743 1.0 0.484739 4.865140e-01\n",
+ "Age:C(Pclass) 0.437814 2.0 0.590357 5.544041e-01\n",
+ "Age:C(Sex):C(Pclass) 1.335653 2.0 1.801023 1.658921e-01\n",
+ "Residual 260.304489 702.0 NaN NaN"
+ ]
+ },
+ "execution_count": 68,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "sm.stats.anova_lm(model, typ=3)"
+ ]
+ },
{
"cell_type": "markdown",
"id": "0e98e30e",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"We can also have a look at the summary tables. This reveals that the residuals are not normally distributed."
]
@@ -440,7 +711,9 @@
"cell_type": "code",
"execution_count": 9,
"id": "8d3ea24c",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -586,67 +859,308 @@
{
"cell_type": "markdown",
"id": "ecf3bcd9",
- "metadata": {
- "heading_collapsed": true
- },
+ "metadata": {},
"source": [
"## Q\n",
"\n",
- "Let us ignore the not-normal residuals and play with post-hoc tests instead.\n",
- "\n",
- "Split the ANOVA for levels of `Pclass` and `Sex`, perform all pairwise comparisons if it make sense, and correct for multiple comparisons.\n",
+ "Because we have a large sample, we will ignore the not-normal residuals and play with post-hoc tests instead.\n",
"\n",
- "We are not interested in the significance of the slope of `Age` for the different levels of the factors."
+ "First proceed considering type-3 sums of squares."
]
},
{
"cell_type": "markdown",
- "id": "7a417d76",
- "metadata": {},
+ "id": "06c5bc0c",
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A"
]
},
+ {
+ "cell_type": "markdown",
+ "id": "f203eea0",
+ "metadata": {
+ "hidden": true
+ },
+ "source": [
+ "We found a single main effect and no interaction. Therefore a single call to `t_test_pairwise` is enough."
+ ]
+ },
{
"cell_type": "code",
- "execution_count": 10,
- "id": "b81407bd",
+ "execution_count": 69,
+ "id": "248c081c",
+ "metadata": {
+ "hidden": true
+ },
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "<div>\n",
+ "<style scoped>\n",
+ " .dataframe tbody tr th:only-of-type {\n",
+ " vertical-align: middle;\n",
+ " }\n",
+ "\n",
+ " .dataframe tbody tr th {\n",
+ " vertical-align: top;\n",
+ " }\n",
+ "\n",
+ " .dataframe thead th {\n",
+ " text-align: right;\n",
+ " }\n",
+ "</style>\n",
+ "<table border=\"1\" class=\"dataframe\">\n",
+ " <thead>\n",
+ " <tr style=\"text-align: right;\">\n",
+ " <th></th>\n",
+ " <th>coef</th>\n",
+ " <th>std err</th>\n",
+ " <th>t</th>\n",
+ " <th>P>|t|</th>\n",
+ " <th>Conf. Int. Low</th>\n",
+ " <th>Conf. Int. Upp.</th>\n",
+ " <th>pvalue-hs</th>\n",
+ " <th>reject-hs</th>\n",
+ " </tr>\n",
+ " </thead>\n",
+ " <tbody>\n",
+ " <tr>\n",
+ " <th>2-1</th>\n",
+ " <td>-1.460550</td>\n",
+ " <td>0.251403</td>\n",
+ " <td>-5.809591</td>\n",
+ " <td>9.496387e-09</td>\n",
+ " <td>-1.954142</td>\n",
+ " <td>-0.966958</td>\n",
+ " <td>1.899277e-08</td>\n",
+ " <td>True</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>3-1</th>\n",
+ " <td>-2.016639</td>\n",
+ " <td>0.217384</td>\n",
+ " <td>-9.276845</td>\n",
+ " <td>2.121528e-19</td>\n",
+ " <td>-2.443439</td>\n",
+ " <td>-1.589838</td>\n",
+ " <td>0.000000e+00</td>\n",
+ " <td>True</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>3-2</th>\n",
+ " <td>-0.556089</td>\n",
+ " <td>0.211313</td>\n",
+ " <td>-2.631590</td>\n",
+ " <td>8.685267e-03</td>\n",
+ " <td>-0.970969</td>\n",
+ " <td>-0.141208</td>\n",
+ " <td>8.685267e-03</td>\n",
+ " <td>True</td>\n",
+ " </tr>\n",
+ " </tbody>\n",
+ "</table>\n",
+ "</div>"
+ ],
+ "text/plain": [
+ " coef std err t P>|t| Conf. Int. Low \\\n",
+ "2-1 -1.460550 0.251403 -5.809591 9.496387e-09 -1.954142 \n",
+ "3-1 -2.016639 0.217384 -9.276845 2.121528e-19 -2.443439 \n",
+ "3-2 -0.556089 0.211313 -2.631590 8.685267e-03 -0.970969 \n",
+ "\n",
+ " Conf. Int. Upp. pvalue-hs reject-hs \n",
+ "2-1 -0.966958 1.899277e-08 True \n",
+ "3-1 -1.589838 0.000000e+00 True \n",
+ "3-2 -0.141208 8.685267e-03 True "
+ ]
+ },
+ "execution_count": 69,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "model.t_test_pairwise('C(Pclass)').result_frame"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "dc0a30ce",
"metadata": {},
- "outputs": [],
"source": [
- "class1_model = smf.ols(data=df, formula='Fare ~ Age * C(Sex)', subset=df['Pclass']==1).fit()\n",
- "class2_model = smf.ols(data=df, formula='Fare ~ Age * C(Sex)', subset=df['Pclass']==2).fit()\n",
- "class3_model = smf.ols(data=df, formula='Fare ~ Age * C(Sex)', subset=df['Pclass']==3).fit()\n",
- "female_model = smf.ols(data=df, formula='Fare ~ Age * C(Pclass)', subset=df['Sex']=='female').fit()\n",
- "male_model = smf.ols(data=df, formula='Fare ~ Age * C(Pclass)', subset=df['Sex']=='male').fit()"
+ "## Q\n",
+ "\n",
+ "Let us suppose we want to use type-1 sums of squares instead.\n",
+ "\n",
+ "Proceed again to performing with `Sex` first, `Pclass` second, and `Age` last.\n",
+ "\n",
+ "In the post-hoc comparisons, we will disregard the effect of the slop of `Age`."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "7a417d76",
+ "metadata": {
+ "heading_collapsed": true
+ },
+ "source": [
+ "## A"
]
},
{
"cell_type": "code",
- "execution_count": 11,
- "id": "cafed8b7",
- "metadata": {},
+ "execution_count": 78,
+ "id": "9466f266",
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
+ "text/html": [
+ "<div>\n",
+ "<style scoped>\n",
+ " .dataframe tbody tr th:only-of-type {\n",
+ " vertical-align: middle;\n",
+ " }\n",
+ "\n",
+ " .dataframe tbody tr th {\n",
+ " vertical-align: top;\n",
+ " }\n",
+ "\n",
+ " .dataframe thead th {\n",
+ " text-align: right;\n",
+ " }\n",
+ "</style>\n",
+ "<table border=\"1\" class=\"dataframe\">\n",
+ " <thead>\n",
+ " <tr style=\"text-align: right;\">\n",
+ " <th></th>\n",
+ " <th>df</th>\n",
+ " <th>sum_sq</th>\n",
+ " <th>mean_sq</th>\n",
+ " <th>F</th>\n",
+ " <th>PR(>F)</th>\n",
+ " </tr>\n",
+ " </thead>\n",
+ " <tbody>\n",
+ " <tr>\n",
+ " <th>C(Sex)</th>\n",
+ " <td>1.0</td>\n",
+ " <td>46.721986</td>\n",
+ " <td>46.721986</td>\n",
+ " <td>126.001800</td>\n",
+ " <td>5.241061e-27</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>C(Pclass)</th>\n",
+ " <td>2.0</td>\n",
+ " <td>318.408489</td>\n",
+ " <td>159.204245</td>\n",
+ " <td>429.348645</td>\n",
+ " <td>1.619836e-122</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>C(Sex):C(Pclass)</th>\n",
+ " <td>2.0</td>\n",
+ " <td>5.990345</td>\n",
+ " <td>2.995172</td>\n",
+ " <td>8.077506</td>\n",
+ " <td>3.402042e-04</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Age</th>\n",
+ " <td>1.0</td>\n",
+ " <td>12.274132</td>\n",
+ " <td>12.274132</td>\n",
+ " <td>33.101391</td>\n",
+ " <td>1.306440e-08</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>C(Sex):Age</th>\n",
+ " <td>1.0</td>\n",
+ " <td>1.486269</td>\n",
+ " <td>1.486269</td>\n",
+ " <td>4.008232</td>\n",
+ " <td>4.566288e-02</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>C(Pclass):Age</th>\n",
+ " <td>2.0</td>\n",
+ " <td>0.296834</td>\n",
+ " <td>0.148417</td>\n",
+ " <td>0.400257</td>\n",
+ " <td>6.703004e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>C(Sex):C(Pclass):Age</th>\n",
+ " <td>2.0</td>\n",
+ " <td>1.335653</td>\n",
+ " <td>0.667827</td>\n",
+ " <td>1.801023</td>\n",
+ " <td>1.658921e-01</td>\n",
+ " </tr>\n",
+ " <tr>\n",
+ " <th>Residual</th>\n",
+ " <td>702.0</td>\n",
+ " <td>260.304489</td>\n",
+ " <td>0.370804</td>\n",
+ " <td>NaN</td>\n",
+ " <td>NaN</td>\n",
+ " </tr>\n",
+ " </tbody>\n",
+ "</table>\n",
+ "</div>"
+ ],
"text/plain": [
- "((186, 4), (173, 4), (355, 4))"
+ " df sum_sq mean_sq F PR(>F)\n",
+ "C(Sex) 1.0 46.721986 46.721986 126.001800 5.241061e-27\n",
+ "C(Pclass) 2.0 318.408489 159.204245 429.348645 1.619836e-122\n",
+ "C(Sex):C(Pclass) 2.0 5.990345 2.995172 8.077506 3.402042e-04\n",
+ "Age 1.0 12.274132 12.274132 33.101391 1.306440e-08\n",
+ "C(Sex):Age 1.0 1.486269 1.486269 4.008232 4.566288e-02\n",
+ "C(Pclass):Age 2.0 0.296834 0.148417 0.400257 6.703004e-01\n",
+ "C(Sex):C(Pclass):Age 2.0 1.335653 0.667827 1.801023 1.658921e-01\n",
+ "Residual 702.0 260.304489 0.370804 NaN NaN"
]
},
- "execution_count": 11,
+ "execution_count": 78,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "class1_model.model.exog.shape, class2_model.model.exog.shape, class3_model.model.exog.shape"
+ "sstype = 1\n",
+ "model = smf.ols('LogFare ~ C(Sex) * C(Pclass) * Age', df).fit()\n",
+ "sm.stats.anova_lm(model)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 88,
+ "id": "b81407bd",
+ "metadata": {
+ "hidden": true
+ },
+ "outputs": [],
+ "source": [
+ "class1_model = smf.ols(data=df, formula='Fare ~ C(Sex) * Age', subset=df['Pclass']==1).fit()\n",
+ "class2_model = smf.ols(data=df, formula='Fare ~ C(Sex) * Age', subset=df['Pclass']==2).fit()\n",
+ "class3_model = smf.ols(data=df, formula='Fare ~ C(Sex) * Age', subset=df['Pclass']==3).fit()\n",
+ "female_model = smf.ols(data=df, formula='Fare ~ C(Pclass) * Age', subset=df['Sex']=='female').fit()\n",
+ "male_model = smf.ols(data=df, formula='Fare ~ C(Pclass) * Age', subset=df['Sex']=='male').fit()"
]
},
{
"cell_type": "code",
- "execution_count": 12,
+ "execution_count": 89,
"id": "7a8fb3ff",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -669,8 +1183,9 @@
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
- " <th>sum_sq</th>\n",
" <th>df</th>\n",
+ " <th>sum_sq</th>\n",
+ " <th>mean_sq</th>\n",
" <th>F</th>\n",
" <th>PR(>F)</th>\n",
" </tr>\n",
@@ -678,29 +1193,33 @@
" <tbody>\n",
" <tr>\n",
" <th>C(Sex)</th>\n",
- " <td>4.043594e+04</td>\n",
" <td>1.0</td>\n",
- " <td>6.627354</td>\n",
- " <td>0.010838</td>\n",
+ " <td>6.251806e+04</td>\n",
+ " <td>62518.055251</td>\n",
+ " <td>10.246561</td>\n",
+ " <td>0.001616</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Age</th>\n",
- " <td>3.572113e+04</td>\n",
" <td>1.0</td>\n",
+ " <td>3.572113e+04</td>\n",
+ " <td>35721.126820</td>\n",
" <td>5.854608</td>\n",
" <td>0.016520</td>\n",
" </tr>\n",
" <tr>\n",
- " <th>Age:C(Sex)</th>\n",
- " <td>8.202655e+02</td>\n",
+ " <th>C(Sex):Age</th>\n",
" <td>1.0</td>\n",
+ " <td>8.202655e+02</td>\n",
+ " <td>820.265529</td>\n",
" <td>0.134440</td>\n",
" <td>0.714299</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Residual</th>\n",
- " <td>1.110449e+06</td>\n",
" <td>182.0</td>\n",
+ " <td>1.110449e+06</td>\n",
+ " <td>6101.369705</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" </tr>\n",
@@ -709,27 +1228,29 @@
"</div>"
],
"text/plain": [
- " sum_sq df F PR(>F)\n",
- "C(Sex) 4.043594e+04 1.0 6.627354 0.010838\n",
- "Age 3.572113e+04 1.0 5.854608 0.016520\n",
- "Age:C(Sex) 8.202655e+02 1.0 0.134440 0.714299\n",
- "Residual 1.110449e+06 182.0 NaN NaN"
+ " df sum_sq mean_sq F PR(>F)\n",
+ "C(Sex) 1.0 6.251806e+04 62518.055251 10.246561 0.001616\n",
+ "Age 1.0 3.572113e+04 35721.126820 5.854608 0.016520\n",
+ "C(Sex):Age 1.0 8.202655e+02 820.265529 0.134440 0.714299\n",
+ "Residual 182.0 1.110449e+06 6101.369705 NaN NaN"
]
},
- "execution_count": 12,
+ "execution_count": 89,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "sm.stats.anova_lm(class1_model, typ=2)"
+ "sm.stats.anova_lm(class1_model, typ=sstype)"
]
},
{
"cell_type": "code",
- "execution_count": 13,
+ "execution_count": 90,
"id": "c5285b28",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -752,8 +1273,9 @@
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
- " <th>sum_sq</th>\n",
" <th>df</th>\n",
+ " <th>sum_sq</th>\n",
+ " <th>mean_sq</th>\n",
" <th>F</th>\n",
" <th>PR(>F)</th>\n",
" </tr>\n",
@@ -761,29 +1283,33 @@
" <tbody>\n",
" <tr>\n",
" <th>C(Sex)</th>\n",
- " <td>9.143118</td>\n",
" <td>1.0</td>\n",
- " <td>0.053774</td>\n",
- " <td>0.816901</td>\n",
+ " <td>29.733469</td>\n",
+ " <td>29.733469</td>\n",
+ " <td>0.174875</td>\n",
+ " <td>0.676346</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Age</th>\n",
- " <td>1140.725330</td>\n",
" <td>1.0</td>\n",
+ " <td>1140.725330</td>\n",
+ " <td>1140.725330</td>\n",
" <td>6.709083</td>\n",
" <td>0.010430</td>\n",
" </tr>\n",
" <tr>\n",
- " <th>Age:C(Sex)</th>\n",
- " <td>7.201063</td>\n",
+ " <th>C(Sex):Age</th>\n",
" <td>1.0</td>\n",
+ " <td>7.201063</td>\n",
+ " <td>7.201063</td>\n",
" <td>0.042352</td>\n",
" <td>0.837197</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Residual</th>\n",
- " <td>28734.566043</td>\n",
" <td>169.0</td>\n",
+ " <td>28734.566043</td>\n",
+ " <td>170.027018</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" </tr>\n",
@@ -792,27 +1318,29 @@
"</div>"
],
"text/plain": [
- " sum_sq df F PR(>F)\n",
- "C(Sex) 9.143118 1.0 0.053774 0.816901\n",
- "Age 1140.725330 1.0 6.709083 0.010430\n",
- "Age:C(Sex) 7.201063 1.0 0.042352 0.837197\n",
- "Residual 28734.566043 169.0 NaN NaN"
+ " df sum_sq mean_sq F PR(>F)\n",
+ "C(Sex) 1.0 29.733469 29.733469 0.174875 0.676346\n",
+ "Age 1.0 1140.725330 1140.725330 6.709083 0.010430\n",
+ "C(Sex):Age 1.0 7.201063 7.201063 0.042352 0.837197\n",
+ "Residual 169.0 28734.566043 170.027018 NaN NaN"
]
},
- "execution_count": 13,
+ "execution_count": 90,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "sm.stats.anova_lm(class2_model, typ=2)"
+ "sm.stats.anova_lm(class2_model, typ=sstype)"
]
},
{
"cell_type": "code",
- "execution_count": 14,
+ "execution_count": 91,
"id": "ed33eab1",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -835,8 +1363,9 @@
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
- " <th>sum_sq</th>\n",
" <th>df</th>\n",
+ " <th>sum_sq</th>\n",
+ " <th>mean_sq</th>\n",
" <th>F</th>\n",
" <th>PR(>F)</th>\n",
" </tr>\n",
@@ -844,29 +1373,33 @@
" <tbody>\n",
" <tr>\n",
" <th>C(Sex)</th>\n",
- " <td>553.194455</td>\n",
" <td>1.0</td>\n",
- " <td>6.099022</td>\n",
- " <td>0.014001</td>\n",
+ " <td>1001.995627</td>\n",
+ " <td>1001.995627</td>\n",
+ " <td>11.047098</td>\n",
+ " <td>0.000982</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Age</th>\n",
- " <td>1970.785958</td>\n",
" <td>1.0</td>\n",
+ " <td>1970.785958</td>\n",
+ " <td>1970.785958</td>\n",
" <td>21.728105</td>\n",
" <td>0.000004</td>\n",
" </tr>\n",
" <tr>\n",
- " <th>Age:C(Sex)</th>\n",
- " <td>896.982415</td>\n",
+ " <th>C(Sex):Age</th>\n",
" <td>1.0</td>\n",
+ " <td>896.982415</td>\n",
+ " <td>896.982415</td>\n",
" <td>9.889317</td>\n",
" <td>0.001804</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Residual</th>\n",
- " <td>31836.456663</td>\n",
" <td>351.0</td>\n",
+ " <td>31836.456663</td>\n",
+ " <td>90.702156</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" </tr>\n",
@@ -875,27 +1408,29 @@
"</div>"
],
"text/plain": [
- " sum_sq df F PR(>F)\n",
- "C(Sex) 553.194455 1.0 6.099022 0.014001\n",
- "Age 1970.785958 1.0 21.728105 0.000004\n",
- "Age:C(Sex) 896.982415 1.0 9.889317 0.001804\n",
- "Residual 31836.456663 351.0 NaN NaN"
+ " df sum_sq mean_sq F PR(>F)\n",
+ "C(Sex) 1.0 1001.995627 1001.995627 11.047098 0.000982\n",
+ "Age 1.0 1970.785958 1970.785958 21.728105 0.000004\n",
+ "C(Sex):Age 1.0 896.982415 896.982415 9.889317 0.001804\n",
+ "Residual 351.0 31836.456663 90.702156 NaN NaN"
]
},
- "execution_count": 14,
+ "execution_count": 91,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "sm.stats.anova_lm(class3_model, typ=2)"
+ "sm.stats.anova_lm(class3_model, typ=sstype)"
]
},
{
"cell_type": "code",
- "execution_count": 15,
+ "execution_count": 92,
"id": "4652591d",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"# we could also make these additional model to treat the 'Age:Sex' interaction, but we are not interested in 'Age' alone (hence the hint about «pairwise»)\n",
@@ -905,9 +1440,11 @@
},
{
"cell_type": "code",
- "execution_count": 16,
+ "execution_count": 93,
"id": "ecfbec15",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -930,8 +1467,9 @@
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
- " <th>sum_sq</th>\n",
" <th>df</th>\n",
+ " <th>sum_sq</th>\n",
+ " <th>mean_sq</th>\n",
" <th>F</th>\n",
" <th>PR(>F)</th>\n",
" </tr>\n",
@@ -939,29 +1477,33 @@
" <tbody>\n",
" <tr>\n",
" <th>C(Pclass)</th>\n",
- " <td>436677.926114</td>\n",
" <td>2.0</td>\n",
- " <td>109.672659</td>\n",
- " <td>4.283604e-35</td>\n",
+ " <td>460882.452649</td>\n",
+ " <td>230441.226324</td>\n",
+ " <td>115.751681</td>\n",
+ " <td>1.699917e-36</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Age</th>\n",
- " <td>4564.359068</td>\n",
" <td>1.0</td>\n",
+ " <td>4564.359068</td>\n",
+ " <td>4564.359068</td>\n",
" <td>2.292698</td>\n",
" <td>1.312222e-01</td>\n",
" </tr>\n",
" <tr>\n",
- " <th>Age:C(Pclass)</th>\n",
- " <td>5386.550285</td>\n",
+ " <th>C(Pclass):Age</th>\n",
" <td>2.0</td>\n",
+ " <td>5386.550285</td>\n",
+ " <td>2693.275142</td>\n",
" <td>1.352844</td>\n",
" <td>2.603528e-01</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Residual</th>\n",
- " <td>507660.125010</td>\n",
" <td>255.0</td>\n",
+ " <td>507660.125010</td>\n",
+ " <td>1990.824020</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" </tr>\n",
@@ -970,27 +1512,29 @@
"</div>"
],
"text/plain": [
- " sum_sq df F PR(>F)\n",
- "C(Pclass) 436677.926114 2.0 109.672659 4.283604e-35\n",
- "Age 4564.359068 1.0 2.292698 1.312222e-01\n",
- "Age:C(Pclass) 5386.550285 2.0 1.352844 2.603528e-01\n",
- "Residual 507660.125010 255.0 NaN NaN"
+ " df sum_sq mean_sq F PR(>F)\n",
+ "C(Pclass) 2.0 460882.452649 230441.226324 115.751681 1.699917e-36\n",
+ "Age 1.0 4564.359068 4564.359068 2.292698 1.312222e-01\n",
+ "C(Pclass):Age 2.0 5386.550285 2693.275142 1.352844 2.603528e-01\n",
+ "Residual 255.0 507660.125010 1990.824020 NaN NaN"
]
},
- "execution_count": 16,
+ "execution_count": 93,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "sm.stats.anova_lm(female_model, typ=2)"
+ "sm.stats.anova_lm(female_model, typ=sstype)"
]
},
{
"cell_type": "code",
- "execution_count": 17,
+ "execution_count": 94,
"id": "27783db9",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1013,8 +1557,9 @@
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
- " <th>sum_sq</th>\n",
" <th>df</th>\n",
+ " <th>sum_sq</th>\n",
+ " <th>mean_sq</th>\n",
" <th>F</th>\n",
" <th>PR(>F)</th>\n",
" </tr>\n",
@@ -1022,29 +1567,33 @@
" <tbody>\n",
" <tr>\n",
" <th>C(Pclass)</th>\n",
- " <td>269207.914985</td>\n",
" <td>2.0</td>\n",
- " <td>90.701809</td>\n",
- " <td>8.657391e-34</td>\n",
+ " <td>255902.065524</td>\n",
+ " <td>127951.032762</td>\n",
+ " <td>86.218789</td>\n",
+ " <td>2.149213e-32</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Age</th>\n",
- " <td>18986.128890</td>\n",
" <td>1.0</td>\n",
+ " <td>18986.128890</td>\n",
+ " <td>18986.128890</td>\n",
" <td>12.793652</td>\n",
" <td>3.857634e-04</td>\n",
" </tr>\n",
" <tr>\n",
- " <th>Age:C(Pclass)</th>\n",
- " <td>11620.048872</td>\n",
+ " <th>C(Pclass):Age</th>\n",
" <td>2.0</td>\n",
+ " <td>11620.048872</td>\n",
+ " <td>5810.024436</td>\n",
" <td>3.915039</td>\n",
" <td>2.062721e-02</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Residual</th>\n",
- " <td>663360.184028</td>\n",
" <td>447.0</td>\n",
+ " <td>663360.184028</td>\n",
+ " <td>1484.027257</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" </tr>\n",
@@ -1053,35 +1602,39 @@
"</div>"
],
"text/plain": [
- " sum_sq df F PR(>F)\n",
- "C(Pclass) 269207.914985 2.0 90.701809 8.657391e-34\n",
- "Age 18986.128890 1.0 12.793652 3.857634e-04\n",
- "Age:C(Pclass) 11620.048872 2.0 3.915039 2.062721e-02\n",
- "Residual 663360.184028 447.0 NaN NaN"
+ " df sum_sq mean_sq F PR(>F)\n",
+ "C(Pclass) 2.0 255902.065524 127951.032762 86.218789 2.149213e-32\n",
+ "Age 1.0 18986.128890 18986.128890 12.793652 3.857634e-04\n",
+ "C(Pclass):Age 2.0 11620.048872 5810.024436 3.915039 2.062721e-02\n",
+ "Residual 447.0 663360.184028 1484.027257 NaN NaN"
]
},
- "execution_count": 17,
+ "execution_count": 94,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "sm.stats.anova_lm(male_model, typ=2)"
+ "sm.stats.anova_lm(male_model, typ=sstype)"
]
},
{
"cell_type": "markdown",
"id": "ce321919",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"We won't include this last case, again, as we are not interested in the slope of `Age` for each level of `Pclass`."
]
},
{
"cell_type": "code",
- "execution_count": 18,
+ "execution_count": 95,
"id": "d27e6f05",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1127,17 +1680,6 @@
" <td>False</td>\n",
" </tr>\n",
" <tr>\n",
- " <th>male-female [2nd class]</th>\n",
- " <td>0.432770</td>\n",
- " <td>4.806509</td>\n",
- " <td>0.090038</td>\n",
- " <td>9.283634e-01</td>\n",
- " <td>-9.055761</td>\n",
- " <td>9.921302</td>\n",
- " <td>9.283634e-01</td>\n",
- " <td>False</td>\n",
- " </tr>\n",
- " <tr>\n",
" <th>2-1 [female]</th>\n",
" <td>-108.472618</td>\n",
" <td>18.421076</td>\n",
@@ -1177,27 +1719,24 @@
"text/plain": [
" coef std err t P>|t| \\\n",
"male-female [1st class] -19.279419 32.487685 -0.593438 5.536250e-01 \n",
- "male-female [2nd class] 0.432770 4.806509 0.090038 9.283634e-01 \n",
"2-1 [female] -108.472618 18.421076 -5.888506 1.224216e-08 \n",
"3-1 [female] -119.359262 15.928392 -7.493491 1.102750e-12 \n",
"3-2 [female] -10.886644 15.483529 -0.703111 4.826278e-01 \n",
"\n",
" Conf. Int. Low Conf. Int. Upp. pvalue-hs \\\n",
"male-female [1st class] -83.380352 44.821514 5.536250e-01 \n",
- "male-female [2nd class] -9.055761 9.921302 9.283634e-01 \n",
"2-1 [female] -144.749437 -72.195799 2.448432e-08 \n",
"3-1 [female] -150.727213 -87.991312 3.308354e-12 \n",
"3-2 [female] -41.378522 19.605233 4.826278e-01 \n",
"\n",
" reject-hs \n",
"male-female [1st class] False \n",
- "male-female [2nd class] False \n",
"2-1 [female] True \n",
"3-1 [female] True \n",
"3-2 [female] False "
]
},
- "execution_count": 18,
+ "execution_count": 95,
"metadata": {},
"output_type": "execute_result"
}
@@ -1209,7 +1748,6 @@
"\n",
"comparisons = pd.concat([\n",
" suffix_label(class1_model.t_test_pairwise('C(Sex)').result_frame, ' [1st class]'),\n",
- " suffix_label(class2_model.t_test_pairwise('C(Sex)').result_frame, ' [2nd class]'),\n",
" suffix_label(female_model.t_test_pairwise('C(Pclass)').result_frame, ' [female]'),\n",
"])\n",
"comparisons"
@@ -1217,19 +1755,20 @@
},
{
"cell_type": "code",
- "execution_count": 19,
+ "execution_count": 96,
"id": "80024bdf",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
"text/plain": [
- "(array([False, False, True, True, False]),\n",
- " array([8.61512929e-01, 9.28363378e-01, 4.89686469e-08, 5.51392265e-12,\n",
- " 8.61512929e-01]))"
+ "(array([False, True, True, False]),\n",
+ " array([7.32326021e-01, 3.67264854e-08, 4.41113812e-12, 7.32326021e-01]))"
]
},
- "execution_count": 19,
+ "execution_count": 96,
"metadata": {},
"output_type": "execute_result"
}
@@ -1240,6 +1779,17 @@
"corrected_rejections, corrected_pvalues"
]
},
+ {
+ "cell_type": "markdown",
+ "id": "5e3fa114",
+ "metadata": {
+ "hidden": true
+ },
+ "source": [
+ "Compared with type-3 sums of squares, we lost one effect whose test-wise $p$-value was very close to the significance threshold.\n",
+ "Of note, this difference comes from the type of sum of squares and not the correction for multiple comparisons."
+ ]
+ },
{
"cell_type": "markdown",
"id": "b7b20012",
@@ -1263,7 +1813,9 @@
{
"cell_type": "markdown",
"id": "35949307",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A"
]
@@ -1272,7 +1824,9 @@
"cell_type": "code",
"execution_count": 20,
"id": "47cf88a3",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"import numpy as np\n",
@@ -1286,7 +1840,9 @@
"cell_type": "code",
"execution_count": 21,
"id": "7a8dacd7",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"mi = pd.read_csv('../data/mi.csv', index_col=0)"
@@ -1296,7 +1852,9 @@
"cell_type": "code",
"execution_count": 22,
"id": "bea3243d",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1513,7 +2071,9 @@
"cell_type": "code",
"execution_count": 23,
"id": "d55877e4",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1555,7 +2115,9 @@
{
"cell_type": "markdown",
"id": "21e8879e",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A"
]
@@ -1564,7 +2126,9 @@
"cell_type": "code",
"execution_count": 45,
"id": "0b8cac58",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1586,7 +2150,9 @@
{
"cell_type": "markdown",
"id": "04a941a5",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"Note that one-liners such as the above expression will almost always be refactored.\n",
"We may need the transformed variable again, and therefore should reify it for future reference."
@@ -1596,7 +2162,9 @@
"cell_type": "code",
"execution_count": 24,
"id": "b0cc612e",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"logPA = np.log(1 + mi['PhysicalActivity'])"
@@ -1605,7 +2173,9 @@
{
"cell_type": "markdown",
"id": "adf75983",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"We may also append it to the dataframe, as a column, just like any other variable."
]
@@ -1614,7 +2184,9 @@
"cell_type": "code",
"execution_count": 25,
"id": "d357ea8a",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"extended_mi = mi.copy()\n",
@@ -1624,7 +2196,9 @@
{
"cell_type": "markdown",
"id": "fd32ce56",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"We cannot compare the skewness of both variables with a single test.\n",
"\n",
@@ -1635,7 +2209,9 @@
"cell_type": "code",
"execution_count": 26,
"id": "87455d20",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1655,7 +2231,9 @@
{
"cell_type": "markdown",
"id": "1067660f",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"...although they both happen to be skewed if we perform individual skewness tests."
]
@@ -1664,7 +2242,9 @@
"cell_type": "code",
"execution_count": 27,
"id": "645ef445",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1685,7 +2265,9 @@
"cell_type": "code",
"execution_count": 28,
"id": "f5984bf0",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1705,7 +2287,9 @@
{
"cell_type": "markdown",
"id": "478d1262",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"Note we do not need the explanatory variable to be symmetric or normally distributed for a model to be valid.\n",
"The point is mainly to make our sample exhibit a good coverage (in a linear sense) of the domain of possible values for our predictors."
@@ -1728,7 +2312,9 @@
{
"cell_type": "markdown",
"id": "c47489e6",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A"
]
@@ -1737,7 +2323,9 @@
"cell_type": "code",
"execution_count": 29,
"id": "bf32c7e6",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1762,7 +2350,9 @@
"cell_type": "code",
"execution_count": 30,
"id": "ac8ab2cc",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1785,7 +2375,9 @@
{
"cell_type": "markdown",
"id": "3cd48b4f",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"In the above example, we leveraged the expressiveness of `patsy` for Wilkinson formulae. The `I` «function» is a special symbol just like `C` for tagging a variable as categorical. `I` allows to evaluate a subexpression following the Python syntax instead of the Wilkinson formalism.\n",
"\n",
@@ -1800,7 +2392,9 @@
"cell_type": "code",
"execution_count": 31,
"id": "177c2ec7",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"# already done\n",
@@ -1815,7 +2409,9 @@
"cell_type": "code",
"execution_count": 32,
"id": "93e6b4ee",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"import statsmodels.api as sm"
@@ -1825,7 +2421,9 @@
"cell_type": "code",
"execution_count": 33,
"id": "d2ec65f0",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"name": "stderr",
@@ -1845,7 +2443,9 @@
{
"cell_type": "markdown",
"id": "c39bd1c5",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"Anyway, as can be seen in the plots above, we turned an influential observation (number 362 was above $0.5$) into a non-influential one, and similarly decreased the influence of several other points.\n",
"\n",
@@ -1869,7 +2469,9 @@
{
"cell_type": "markdown",
"id": "12d7e9a5",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A (with nested Q&A)"
]
@@ -1878,7 +2480,9 @@
"cell_type": "code",
"execution_count": 34,
"id": "82e9f8d6",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -1988,7 +2592,9 @@
"cell_type": "code",
"execution_count": 35,
"id": "faeab021",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -2103,7 +2709,9 @@
{
"cell_type": "markdown",
"id": "66ddaa3e",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"### Q\n",
"\n",
@@ -2115,7 +2723,10 @@
{
"cell_type": "markdown",
"id": "05844c56",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true,
+ "hidden": true
+ },
"source": [
"### A"
]
@@ -2124,7 +2735,9 @@
"cell_type": "code",
"execution_count": 36,
"id": "2588d8d6",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -2248,7 +2861,9 @@
{
"cell_type": "markdown",
"id": "c7f00836",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"In `seaborn`, a dot plot can be drawn with the `stripplot` function."
]
@@ -2257,7 +2872,9 @@
"cell_type": "code",
"execution_count": 37,
"id": "625111f9",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -2280,7 +2897,9 @@
{
"cell_type": "markdown",
"id": "1ba5537f",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"source": [
"Should we add a term, `BMI` brings more information alone than `Heart:logPA` for the same number of model parameters. This would be confirmed by AIC and BIC.\n",
"\n",
@@ -2310,7 +2929,9 @@
{
"cell_type": "markdown",
"id": "a6f37b2e",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A"
]
@@ -2319,7 +2940,9 @@
"cell_type": "code",
"execution_count": 38,
"id": "77e350be",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"# already done\n",
@@ -2334,7 +2957,9 @@
"cell_type": "code",
"execution_count": 39,
"id": "377783b7",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -2358,7 +2983,9 @@
"cell_type": "code",
"execution_count": 40,
"id": "e07bc434",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -2395,7 +3022,9 @@
{
"cell_type": "markdown",
"id": "0a822074",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A"
]
@@ -2404,7 +3033,9 @@
"cell_type": "code",
"execution_count": 41,
"id": "85d73982",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -2444,7 +3075,9 @@
{
"cell_type": "markdown",
"id": "7ce19023",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## A"
]
@@ -2453,7 +3086,9 @@
"cell_type": "code",
"execution_count": 42,
"id": "7dec1d91",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": [
"logPA = np.log(1 + mi['PhysicalActivity'])\n",
@@ -2473,7 +3108,9 @@
{
"cell_type": "markdown",
"id": "b588b5c1",
- "metadata": {},
+ "metadata": {
+ "heading_collapsed": true
+ },
"source": [
"## Q\n",
"\n",
@@ -2495,11 +3132,23 @@
"Compute the statistic $nR^2$ and the resulting $p$-value."
]
},
+ {
+ "cell_type": "markdown",
+ "id": "e4207c4c",
+ "metadata": {
+ "heading_collapsed": true
+ },
+ "source": [
+ "## A"
+ ]
+ },
{
"cell_type": "code",
"execution_count": 43,
"id": "db59adce",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [
{
"data": {
@@ -2525,7 +3174,9 @@
"cell_type": "code",
"execution_count": null,
"id": "1bac9cde",
- "metadata": {},
+ "metadata": {
+ "hidden": true
+ },
"outputs": [],
"source": []
}
diff --git a/notebooks/statsmodels_cours.ipynb b/notebooks/statsmodels_cours.ipynb
index 17502d5b7f7010d8514718e893777b16088858fb..1b95f5d5d641ed85e12ba7e8f4d4bda9ff5719e8 100644
--- a/notebooks/statsmodels_cours.ipynb
+++ b/notebooks/statsmodels_cours.ipynb
@@ -44,7 +44,7 @@
},
{
"cell_type": "code",
- "execution_count": 98,
+ "execution_count": 151,
"id": "d6590258-e1ac-4f62-8adf-3fa014376a22",
"metadata": {},
"outputs": [],
@@ -55,7 +55,7 @@
},
{
"cell_type": "code",
- "execution_count": 3,
+ "execution_count": 152,
"id": "7b4f63c1-9995-4516-aeb5-1fb4fe1db896",
"metadata": {},
"outputs": [],
@@ -84,7 +84,7 @@
},
{
"cell_type": "code",
- "execution_count": 99,
+ "execution_count": 153,
"id": "aec2f6b1-c4fc-465e-9434-b32333770e24",
"metadata": {},
"outputs": [],
@@ -96,7 +96,7 @@
},
{
"cell_type": "code",
- "execution_count": 100,
+ "execution_count": 154,
"id": "bd00bd2f-b7eb-47e4-8af4-187ce9473481",
"metadata": {},
"outputs": [
@@ -106,7 +106,7 @@
"F_onewayResult(statistic=2.3575322551335636, pvalue=0.11384795345837218)"
]
},
- "execution_count": 100,
+ "execution_count": 154,
"metadata": {},
"output_type": "execute_result"
}
@@ -126,7 +126,7 @@
},
{
"cell_type": "code",
- "execution_count": 101,
+ "execution_count": 155,
"id": "696ff83f-d827-4513-9801-51488e5a1df0",
"metadata": {},
"outputs": [
@@ -140,7 +140,7 @@
" 'C', 'C', 'C', 'C'], dtype='<U1'))"
]
},
- "execution_count": 101,
+ "execution_count": 155,
"metadata": {},
"output_type": "execute_result"
}
@@ -153,7 +153,7 @@
},
{
"cell_type": "code",
- "execution_count": 102,
+ "execution_count": 156,
"id": "c324bd0f-e769-4a0d-8e6b-d4d346e76f68",
"metadata": {},
"outputs": [
@@ -371,7 +371,7 @@
"29 81 C"
]
},
- "execution_count": 102,
+ "execution_count": 156,
"metadata": {},
"output_type": "execute_result"
}
@@ -402,7 +402,7 @@
},
{
"cell_type": "code",
- "execution_count": 103,
+ "execution_count": 157,
"id": "4a2a00ec-43c6-4299-82e3-15b807110828",
"metadata": {},
"outputs": [],
@@ -420,7 +420,7 @@
},
{
"cell_type": "code",
- "execution_count": 104,
+ "execution_count": 158,
"id": "21bf5c1d-caa3-4072-bc86-ce860b8f33c8",
"metadata": {},
"outputs": [
@@ -434,7 +434,7 @@
"Model: OLS Adj. R-squared: 0.086\n",
"Method: Least Squares F-statistic: 2.358\n",
"Date: Tue, 28 Sep 2021 Prob (F-statistic): 0.114\n",
- "Time: 11:16:16 Log-Likelihood: -96.604\n",
+ "Time: 18:26:09 Log-Likelihood: -96.604\n",
"No. Observations: 30 AIC: 199.2\n",
"Df Residuals: 27 BIC: 203.4\n",
"Df Model: 2 \n",
@@ -473,7 +473,7 @@
},
{
"cell_type": "code",
- "execution_count": 105,
+ "execution_count": 159,
"id": "638bdd6b-6964-4209-b762-2991fd3fb7fc",
"metadata": {
"tags": []
@@ -485,7 +485,7 @@
"F_onewayResult(statistic=2.3575322551335636, pvalue=0.11384795345837218)"
]
},
- "execution_count": 105,
+ "execution_count": 159,
"metadata": {},
"output_type": "execute_result"
}
@@ -507,7 +507,7 @@
},
{
"cell_type": "code",
- "execution_count": 106,
+ "execution_count": 166,
"id": "0fe841a0-21b1-4c92-a767-c4ab3abd37f8",
"metadata": {
"tags": []
@@ -517,9 +517,9 @@
"name": "stdout",
"output_type": "stream",
"text": [
- " df sum_sq mean_sq F PR(>F)\n",
- "Group 2.0 192.2 96.100000 2.357532 0.113848\n",
- "Residual 27.0 1100.6 40.762963 NaN NaN\n"
+ " sum_sq df F PR(>F)\n",
+ "Group 192.2 2.0 2.357532 0.113848\n",
+ "Residual 1100.6 27.0 NaN NaN\n"
]
}
],
@@ -1736,7 +1736,7 @@
},
{
"cell_type": "code",
- "execution_count": 148,
+ "execution_count": 167,
"id": "cc06b3ed-a066-4de3-884c-a7c82729c359",
"metadata": {},
"outputs": [
@@ -1744,15 +1744,16 @@
"name": "stdout",
"output_type": "stream",
"text": [
- " sum_sq df F PR(>F)\n",
- "water 15.552000 1.0 16.034261 0.000462\n",
- "sun 21.424667 2.0 11.044518 0.000337\n",
- "Residual 25.218000 26.0 NaN NaN\n"
+ " sum_sq df F PR(>F)\n",
+ "Intercept 394.218750 1.0 406.443314 2.139588e-17\n",
+ "water 15.552000 1.0 16.034261 4.623155e-04\n",
+ "sun 21.424667 2.0 11.044518 3.373296e-04\n",
+ "Residual 25.218000 26.0 NaN NaN\n"
]
}
],
"source": [
- "anova_table = sm.stats.anova_lm(plant_model, typ=2)\n",
+ "anova_table = sm.stats.anova_lm(plant_model, typ=3) # typ specifies the type of sum of squares\n",
"print(anova_table)"
]
},
@@ -1886,7 +1887,7 @@
},
{
"cell_type": "code",
- "execution_count": 149,
+ "execution_count": 168,
"id": "60de13e5-b798-4312-8d06-0ef79f480760",
"metadata": {},
"outputs": [
@@ -1894,18 +1895,19 @@
"name": "stdout",
"output_type": "stream",
"text": [
- " sum_sq df F PR(>F)\n",
- "water 15.552000 1.0 19.117394 0.000205\n",
- "sun 21.424667 2.0 13.168203 0.000138\n",
- "water:sun 5.694000 2.0 3.499693 0.046376\n",
- "Residual 19.524000 24.0 NaN NaN\n"
+ " sum_sq df F PR(>F)\n",
+ "Intercept 246.402000 1.0 302.891211 4.094979e-15\n",
+ "water 2.401000 1.0 2.951444 9.867817e-02\n",
+ "sun 8.041333 2.0 4.942430 1.593920e-02\n",
+ "water:sun 5.694000 2.0 3.499693 4.637649e-02\n",
+ "Residual 19.524000 24.0 NaN NaN\n"
]
}
],
"source": [
"model_with_interaction = ols('height ~ water * sun', data=plant_data).fit()\n",
"# remember `water * sun` is equivalent to `water + sun + water:sun`\n",
- "print(sm.stats.anova_lm(model_with_interaction, typ=2))"
+ "print(sm.stats.anova_lm(model_with_interaction, typ=3))"
]
},
{
@@ -1971,7 +1973,7 @@
"\n",
"y_low_daily = w['Intercept'] + w['sun[T.low]']\n",
"y_low_weekly = w['Intercept'] + w['sun[T.low]'] + w['water[T.weekly]'] + w['water[T.weekly]:sun[T.low]']\n",
- "ax.plot([x[0]-dx, x[0]+dx], [y_low_daily, y_low_weekly], 'k-d', markerfacecolor='w')\n",
+ "ax.plot([x[0]-dx, x[0]+dx], [y_low_daily, y_low_weekly], 'k-d', markerfacecolor='w')Q\n",
"\n",
"y_med_daily = w['Intercept'] + w['sun[T.med]']\n",
"y_med_weekly = w['Intercept'] + w['sun[T.med]'] + w['water[T.weekly]'] + w['water[T.weekly]:sun[T.med]']\n",
@@ -2069,7 +2071,7 @@
"daily_water_model = ols('height ~ sun', data=plant_data[plant_data['water']=='daily']).fit()\n",
"weekly_water_model = ols('height ~ sun', data=plant_data[plant_data['water']=='weekly']).fit()\n",
"low_sun_model = ols('height ~ water', data=plant_data[plant_data['sun']=='low']).fit()\n",
- "med_sun_model = ols('height ~ water', data=plant_data[plant_data['sun']=='med']).fit()\n",
+ "med_sun_model = ols('height ~ water', data=plant_data[plant_daQta['sun']=='med']).fit()\n",
"high_sun_model = ols('height ~ water', data=plant_data[plant_data['sun']=='high']).fit()"
]
},
@@ -2091,7 +2093,7 @@
},
{
"cell_type": "code",
- "execution_count": 136,
+ "execution_count": 169,
"id": "93ccdd87",
"metadata": {},
"outputs": [
@@ -2101,7 +2103,7 @@
"0.03098093333325329"
]
},
- "execution_count": 136,
+ "execution_count": 169,
"metadata": {},
"output_type": "execute_result"
}
@@ -2112,7 +2114,7 @@
},
{
"cell_type": "code",
- "execution_count": 137,
+ "execution_count": 170,
"id": "d1392464",
"metadata": {},
"outputs": [
@@ -2197,7 +2199,7 @@
"med-low 0.254253 False "
]
},
- "execution_count": 137,
+ "execution_count": 170,
"metadata": {},
"output_type": "execute_result"
}
@@ -2409,7 +2411,7 @@
},
{
"cell_type": "code",
- "execution_count": 126,
+ "execution_count": 171,
"id": "e6c92946",
"metadata": {},
"outputs": [
@@ -2434,43 +2436,45 @@
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
- " <th>df</th>\n",
" <th>sum_sq</th>\n",
- " <th>mean_sq</th>\n",
+ " <th>df</th>\n",
" <th>F</th>\n",
" <th>PR(>F)</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
+ " <th>Intercept</th>\n",
+ " <td>246.402000</td>\n",
+ " <td>1.0</td>\n",
+ " <td>302.891211</td>\n",
+ " <td>4.094979e-15</td>\n",
+ " </tr>\n",
+ " <tr>\n",
" <th>water</th>\n",
+ " <td>2.401000</td>\n",
" <td>1.0</td>\n",
- " <td>15.552000</td>\n",
- " <td>15.552000</td>\n",
- " <td>19.117394</td>\n",
- " <td>0.000205</td>\n",
+ " <td>2.951444</td>\n",
+ " <td>9.867817e-02</td>\n",
" </tr>\n",
" <tr>\n",
" <th>sun</th>\n",
+ " <td>8.041333</td>\n",
" <td>2.0</td>\n",
- " <td>21.424667</td>\n",
- " <td>10.712333</td>\n",
- " <td>13.168203</td>\n",
- " <td>0.000138</td>\n",
+ " <td>4.942430</td>\n",
+ " <td>1.593920e-02</td>\n",
" </tr>\n",
" <tr>\n",
" <th>water:sun</th>\n",
- " <td>2.0</td>\n",
" <td>5.694000</td>\n",
- " <td>2.847000</td>\n",
+ " <td>2.0</td>\n",
" <td>3.499693</td>\n",
- " <td>0.046376</td>\n",
+ " <td>4.637649e-02</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Residual</th>\n",
- " <td>24.0</td>\n",
" <td>19.524000</td>\n",
- " <td>0.813500</td>\n",
+ " <td>24.0</td>\n",
" <td>NaN</td>\n",
" <td>NaN</td>\n",
" </tr>\n",
@@ -2479,20 +2483,21 @@
"</div>"
],
"text/plain": [
- " df sum_sq mean_sq F PR(>F)\n",
- "water 1.0 15.552000 15.552000 19.117394 0.000205\n",
- "sun 2.0 21.424667 10.712333 13.168203 0.000138\n",
- "water:sun 2.0 5.694000 2.847000 3.499693 0.046376\n",
- "Residual 24.0 19.524000 0.813500 NaN NaN"
+ " sum_sq df F PR(>F)\n",
+ "Intercept 246.402000 1.0 302.891211 4.094979e-15\n",
+ "water 2.401000 1.0 2.951444 9.867817e-02\n",
+ "sun 8.041333 2.0 4.942430 1.593920e-02\n",
+ "water:sun 5.694000 2.0 3.499693 4.637649e-02\n",
+ "Residual 19.524000 24.0 NaN NaN"
]
},
- "execution_count": 126,
+ "execution_count": 171,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
- "sm.stats.anova_lm(model_with_interaction)"
+ "sm.stats.anova_lm(model_with_interaction, typ=3)"
]
},
{