Skip to content

Instantly share code, notes, and snippets.

@josef-pkt

josef-pkt/ex_multivariatels.ipynb

Created June 22, 2023 17:35
Show Gist options
  • Save josef-pkt/8a492802eee3896c47b646ee11beb1f6 to your computer and use it in GitHub Desktop.
Save josef-pkt/8a492802eee3896c47b646ee11beb1f6 to your computer and use it in GitHub Desktop.
Download ZIP
Multivariate Linear model in statsmodels https://github.com/statsmodels/statsmodels/pull/8919
Display the source blob
Display the rendered blob
Raw
ex_multivariatels.ipynb
{
"cells": [
{
"cell_type": "markdown",
"id": "c53d6f3b",
"metadata": {},
"source": [
"## Multivariate Linear Model - MultivariateLS\n",
"\n",
"This notebooks illustrates some features for the multivariate linear model estimated by least squares. \n",
"The example is based on the UCLA stats example in https://stats.oarc.ucla.edu/stata/dae/multivariate-regression-analysis/ .\n",
"\n",
"The model assumes that a multivariate dependent variable depends linearly on the same set of explanatory variables.\n",
"\n",
"Y = X * B + u\n",
"\n",
"where \n",
"- the dependent variable (endog) `Y` has shape (nobs, k_endog), \n",
"- the matrix of explanatory variables including constant (exog) `X` has shape (nobs, k_exog), and\n",
"- the parameter matrix `B` has shape (k_exog, k_endog), i.e. coefficients for explanatory variables in rows and dependent variables in columns.\n",
"- the disturbance term `u` has the same shape as `Y`, (nobs, k_endog), and is assumed to have mean zero and to be uncorrelated with the exog `X`.\n",
"\n",
"Estimation is by least squares. The parameter estimates with common explanatory variables for each dependent variables corresponds to separate OLS estimates for each `endog`. The main advantage of the multivariate model is that we can make inference "
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "9f00bd5f",
"metadata": {},
"outputs": [],
"source": [
"import os\n",
"import numpy as np\n",
"\n",
"import pandas as pd\n",
"\n",
"from statsmodels.base.model import LikelihoodModel\n",
"from statsmodels.regression.linear_model import OLS\n",
"from statsmodels.multivariate.manova import MANOVA\n",
"from statsmodels.multivariate.multivariate_ols import MultivariateLS\n",
"\n",
"import statsmodels.multivariate.tests.results as path\n",
"dir_path = os.path.dirname(os.path.abspath(path.__file__))\n",
"csv_path = os.path.join(dir_path, 'mvreg.csv')\n",
"data_mvreg = pd.read_csv(csv_path)"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "3cf0ef83",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>locus_of_control</th>\n",
" <th>self_concept</th>\n",
" <th>motivation</th>\n",
" <th>read</th>\n",
" <th>write</th>\n",
" <th>science</th>\n",
" <th>prog</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>-1.143955</td>\n",
" <td>0.722641</td>\n",
" <td>0.368973</td>\n",
" <td>37.405548</td>\n",
" <td>39.032845</td>\n",
" <td>33.532822</td>\n",
" <td>academic</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>0.504134</td>\n",
" <td>0.111364</td>\n",
" <td>0.520319</td>\n",
" <td>52.760784</td>\n",
" <td>51.995041</td>\n",
" <td>65.225044</td>\n",
" <td>academic</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>1.628546</td>\n",
" <td>0.629934</td>\n",
" <td>0.436838</td>\n",
" <td>59.771915</td>\n",
" <td>54.651653</td>\n",
" <td>64.604500</td>\n",
" <td>academic</td>\n",
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <td>0.368096</td>\n",
" <td>-0.138528</td>\n",
" <td>-0.004324</td>\n",
" <td>42.854324</td>\n",
" <td>41.121357</td>\n",
" <td>48.493809</td>\n",
" <td>vocational</td>\n",
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <td>-0.280190</td>\n",
" <td>-0.452226</td>\n",
" <td>1.256924</td>\n",
" <td>54.756279</td>\n",
" <td>49.947208</td>\n",
" <td>50.381657</td>\n",
" <td>academic</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" locus_of_control self_concept motivation read write \\\n",
"0 -1.143955 0.722641 0.368973 37.405548 39.032845 \n",
"1 0.504134 0.111364 0.520319 52.760784 51.995041 \n",
"2 1.628546 0.629934 0.436838 59.771915 54.651653 \n",
"3 0.368096 -0.138528 -0.004324 42.854324 41.121357 \n",
"4 -0.280190 -0.452226 1.256924 54.756279 49.947208 \n",
"\n",
" science prog \n",
"0 33.532822 academic \n",
"1 65.225044 academic \n",
"2 64.604500 academic \n",
"3 48.493809 vocational \n",
"4 50.381657 academic "
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"data_mvreg.head()"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "0b249078",
"metadata": {},
"outputs": [],
"source": [
"formula = \"locus_of_control + self_concept + motivation ~ read + write + science + prog\"\n",
"mod = MultivariateLS.from_formula(formula, data=data_mvreg)\n",
"res = mod.fit()"
]
},
{
"cell_type": "markdown",
"id": "4654d4bb",
"metadata": {},
"source": [
"## Multivariate hypothesis tests mv_test\n",
"\n",
"The `mv_test` method by default performs the hypothesis tests that each term in the formula has no effect on any of the dependent variables (`endog`). This is the same as the MANOVA test. \n",
"Note, MANOVA in statsmodels is implemented as test on coefficients in the multivariate model and is not restricted to categorical variables. In the current example, we have three continuous and one categorical explanatory variables, in addition to the constant.\n",
"\n",
"Consequently, using mv_test in MultivariateLS and in MANOVA produces the same results.\n",
"However. MANOVA only provides the hypothesis tests as feature, while MultivariateLS provide the usual model results.\n",
"\n",
"More general versions of the mv_test are for hypothesis in the form `L B M = C`.\n",
"Here `L` are restrictions corresponding to explanatory variables, `M` are restrictions corresponding to dependent (endog) variables and `C` is a matrix of constants for affine restrictions."
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "642ede51",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th></th>\n",
" <th>Value</th>\n",
" <th>Num DF</th>\n",
" <th>Den DF</th>\n",
" <th>F Value</th>\n",
" <th>Pr &gt; F</th>\n",
" </tr>\n",
" <tr>\n",
" <th>Effect</th>\n",
" <th>Statistic</th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">Intercept</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.848467</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>35.242876</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.151533</td>\n",
" <td>3.0</td>\n",
" <td>592.0</td>\n",
" <td>35.242876</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.178596</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>35.242876</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.178596</td>\n",
" <td>3</td>\n",
" <td>592</td>\n",
" <td>35.242876</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">prog</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.891438</td>\n",
" <td>6</td>\n",
" <td>1184.0</td>\n",
" <td>11.670765</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.108649</td>\n",
" <td>6.0</td>\n",
" <td>1186.0</td>\n",
" <td>11.354963</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.121685</td>\n",
" <td>6</td>\n",
" <td>787.558061</td>\n",
" <td>11.996155</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.120878</td>\n",
" <td>3</td>\n",
" <td>593</td>\n",
" <td>23.893456</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">read</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.976425</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>4.764416</td>\n",
" <td>0.002727</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.023575</td>\n",
" <td>3.0</td>\n",
" <td>592.0</td>\n",
" <td>4.764416</td>\n",
" <td>0.002727</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.024144</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>4.764416</td>\n",
" <td>0.002727</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.024144</td>\n",
" <td>3</td>\n",
" <td>592</td>\n",
" <td>4.764416</td>\n",
" <td>0.002727</td>\n",
" </tr>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">write</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.947394</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>10.957338</td>\n",
" <td>0.000001</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.052606</td>\n",
" <td>3.0</td>\n",
" <td>592.0</td>\n",
" <td>10.957338</td>\n",
" <td>0.000001</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.055527</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>10.957338</td>\n",
" <td>0.000001</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.055527</td>\n",
" <td>3</td>\n",
" <td>592</td>\n",
" <td>10.957338</td>\n",
" <td>0.000001</td>\n",
" </tr>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">science</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.983405</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>3.329911</td>\n",
" <td>0.019305</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.016595</td>\n",
" <td>3.0</td>\n",
" <td>592.0</td>\n",
" <td>3.329911</td>\n",
" <td>0.019305</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.016875</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>3.329911</td>\n",
" <td>0.019305</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.016875</td>\n",
" <td>3</td>\n",
" <td>592</td>\n",
" <td>3.329911</td>\n",
" <td>0.019305</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" Value Num DF Den DF F Value \\\n",
"Effect Statistic \n",
"Intercept Wilks' lambda 0.848467 3 592.0 35.242876 \n",
" Pillai's trace 0.151533 3.0 592.0 35.242876 \n",
" Hotelling-Lawley trace 0.178596 3 592.0 35.242876 \n",
" Roy's greatest root 0.178596 3 592 35.242876 \n",
"prog Wilks' lambda 0.891438 6 1184.0 11.670765 \n",
" Pillai's trace 0.108649 6.0 1186.0 11.354963 \n",
" Hotelling-Lawley trace 0.121685 6 787.558061 11.996155 \n",
" Roy's greatest root 0.120878 3 593 23.893456 \n",
"read Wilks' lambda 0.976425 3 592.0 4.764416 \n",
" Pillai's trace 0.023575 3.0 592.0 4.764416 \n",
" Hotelling-Lawley trace 0.024144 3 592.0 4.764416 \n",
" Roy's greatest root 0.024144 3 592 4.764416 \n",
"write Wilks' lambda 0.947394 3 592.0 10.957338 \n",
" Pillai's trace 0.052606 3.0 592.0 10.957338 \n",
" Hotelling-Lawley trace 0.055527 3 592.0 10.957338 \n",
" Roy's greatest root 0.055527 3 592 10.957338 \n",
"science Wilks' lambda 0.983405 3 592.0 3.329911 \n",
" Pillai's trace 0.016595 3.0 592.0 3.329911 \n",
" Hotelling-Lawley trace 0.016875 3 592.0 3.329911 \n",
" Roy's greatest root 0.016875 3 592 3.329911 \n",
"\n",
" Pr > F \n",
"Effect Statistic \n",
"Intercept Wilks' lambda 0.0 \n",
" Pillai's trace 0.0 \n",
" Hotelling-Lawley trace 0.0 \n",
" Roy's greatest root 0.0 \n",
"prog Wilks' lambda 0.0 \n",
" Pillai's trace 0.0 \n",
" Hotelling-Lawley trace 0.0 \n",
" Roy's greatest root 0.0 \n",
"read Wilks' lambda 0.002727 \n",
" Pillai's trace 0.002727 \n",
" Hotelling-Lawley trace 0.002727 \n",
" Roy's greatest root 0.002727 \n",
"write Wilks' lambda 0.000001 \n",
" Pillai's trace 0.000001 \n",
" Hotelling-Lawley trace 0.000001 \n",
" Roy's greatest root 0.000001 \n",
"science Wilks' lambda 0.019305 \n",
" Pillai's trace 0.019305 \n",
" Hotelling-Lawley trace 0.019305 \n",
" Roy's greatest root 0.019305 "
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"mvt = res.mv_test()\n",
"mvt.summary_frame"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "2f94f3ab",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th></th>\n",
" <th>Value</th>\n",
" <th>Num DF</th>\n",
" <th>Den DF</th>\n",
" <th>F Value</th>\n",
" <th>Pr &gt; F</th>\n",
" </tr>\n",
" <tr>\n",
" <th>Effect</th>\n",
" <th>Statistic</th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">Intercept</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.848467</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>35.242876</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.151533</td>\n",
" <td>3.0</td>\n",
" <td>592.0</td>\n",
" <td>35.242876</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.178596</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>35.242876</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.178596</td>\n",
" <td>3</td>\n",
" <td>592</td>\n",
" <td>35.242876</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">prog</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.891438</td>\n",
" <td>6</td>\n",
" <td>1184.0</td>\n",
" <td>11.670765</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.108649</td>\n",
" <td>6.0</td>\n",
" <td>1186.0</td>\n",
" <td>11.354963</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.121685</td>\n",
" <td>6</td>\n",
" <td>787.558061</td>\n",
" <td>11.996155</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.120878</td>\n",
" <td>3</td>\n",
" <td>593</td>\n",
" <td>23.893456</td>\n",
" <td>0.0</td>\n",
" </tr>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">read</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.976425</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>4.764416</td>\n",
" <td>0.002727</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.023575</td>\n",
" <td>3.0</td>\n",
" <td>592.0</td>\n",
" <td>4.764416</td>\n",
" <td>0.002727</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.024144</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>4.764416</td>\n",
" <td>0.002727</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.024144</td>\n",
" <td>3</td>\n",
" <td>592</td>\n",
" <td>4.764416</td>\n",
" <td>0.002727</td>\n",
" </tr>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">write</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.947394</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>10.957338</td>\n",
" <td>0.000001</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.052606</td>\n",
" <td>3.0</td>\n",
" <td>592.0</td>\n",
" <td>10.957338</td>\n",
" <td>0.000001</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.055527</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>10.957338</td>\n",
" <td>0.000001</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.055527</td>\n",
" <td>3</td>\n",
" <td>592</td>\n",
" <td>10.957338</td>\n",
" <td>0.000001</td>\n",
" </tr>\n",
" <tr>\n",
" <th rowspan=\"4\" valign=\"top\">science</th>\n",
" <th>Wilks' lambda</th>\n",
" <td>0.983405</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>3.329911</td>\n",
" <td>0.019305</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Pillai's trace</th>\n",
" <td>0.016595</td>\n",
" <td>3.0</td>\n",
" <td>592.0</td>\n",
" <td>3.329911</td>\n",
" <td>0.019305</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Hotelling-Lawley trace</th>\n",
" <td>0.016875</td>\n",
" <td>3</td>\n",
" <td>592.0</td>\n",
" <td>3.329911</td>\n",
" <td>0.019305</td>\n",
" </tr>\n",
" <tr>\n",
" <th>Roy's greatest root</th>\n",
" <td>0.016875</td>\n",
" <td>3</td>\n",
" <td>592</td>\n",
" <td>3.329911</td>\n",
" <td>0.019305</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" Value Num DF Den DF F Value \\\n",
"Effect Statistic \n",
"Intercept Wilks' lambda 0.848467 3 592.0 35.242876 \n",
" Pillai's trace 0.151533 3.0 592.0 35.242876 \n",
" Hotelling-Lawley trace 0.178596 3 592.0 35.242876 \n",
" Roy's greatest root 0.178596 3 592 35.242876 \n",
"prog Wilks' lambda 0.891438 6 1184.0 11.670765 \n",
" Pillai's trace 0.108649 6.0 1186.0 11.354963 \n",
" Hotelling-Lawley trace 0.121685 6 787.558061 11.996155 \n",
" Roy's greatest root 0.120878 3 593 23.893456 \n",
"read Wilks' lambda 0.976425 3 592.0 4.764416 \n",
" Pillai's trace 0.023575 3.0 592.0 4.764416 \n",
" Hotelling-Lawley trace 0.024144 3 592.0 4.764416 \n",
" Roy's greatest root 0.024144 3 592 4.764416 \n",
"write Wilks' lambda 0.947394 3 592.0 10.957338 \n",
" Pillai's trace 0.052606 3.0 592.0 10.957338 \n",
" Hotelling-Lawley trace 0.055527 3 592.0 10.957338 \n",
" Roy's greatest root 0.055527 3 592 10.957338 \n",
"science Wilks' lambda 0.983405 3 592.0 3.329911 \n",
" Pillai's trace 0.016595 3.0 592.0 3.329911 \n",
" Hotelling-Lawley trace 0.016875 3 592.0 3.329911 \n",
" Roy's greatest root 0.016875 3 592 3.329911 \n",
"\n",
" Pr > F \n",
"Effect Statistic \n",
"Intercept Wilks' lambda 0.0 \n",
" Pillai's trace 0.0 \n",
" Hotelling-Lawley trace 0.0 \n",
" Roy's greatest root 0.0 \n",
"prog Wilks' lambda 0.0 \n",
" Pillai's trace 0.0 \n",
" Hotelling-Lawley trace 0.0 \n",
" Roy's greatest root 0.0 \n",
"read Wilks' lambda 0.002727 \n",
" Pillai's trace 0.002727 \n",
" Hotelling-Lawley trace 0.002727 \n",
" Roy's greatest root 0.002727 \n",
"write Wilks' lambda 0.000001 \n",
" Pillai's trace 0.000001 \n",
" Hotelling-Lawley trace 0.000001 \n",
" Roy's greatest root 0.000001 \n",
"science Wilks' lambda 0.019305 \n",
" Pillai's trace 0.019305 \n",
" Hotelling-Lawley trace 0.019305 \n",
" Roy's greatest root 0.019305 "
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"manova = MANOVA.from_formula(formula, data=data_mvreg)\n",
"manova.mv_test().summary_frame"
]
},
{
"cell_type": "markdown",
"id": "ae0f0226",
"metadata": {},
"source": [
"The core multivariate regression results are displayed by the `summary` method."
]
},
{
"cell_type": "code",
"execution_count": 6,
"id": "8fa8e5ab",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
" MultivariateLS Regression Results \n",
"==============================================================================================================\n",
"Dep. Variable: ['locus_of_control', 'self_concept', 'motivation'] No. Observations: 600\n",
"Model: MultivariateLS Df Residuals: 594\n",
"Method: Least Squares Df Model: 15\n",
"Date: Thu, 22 Jun 2023 \n",
"Time: 13:26:17 \n",
"======================================================================================\n",
" locus_of_control coef std err t P>|t| [0.025 0.975]\n",
"--------------------------------------------------------------------------------------\n",
"Intercept -1.4970 0.157 -9.505 0.000 -1.806 -1.188\n",
"prog[T.general] -0.1278 0.064 -1.998 0.046 -0.253 -0.002\n",
"prog[T.vocational] 0.1239 0.058 2.150 0.032 0.011 0.237\n",
"read 0.0125 0.004 3.363 0.001 0.005 0.020\n",
"write 0.0121 0.003 3.581 0.000 0.005 0.019\n",
"science 0.0058 0.004 1.582 0.114 -0.001 0.013\n",
"--------------------------------------------------------------------------------------\n",
" self_concept coef std err t P>|t| [0.025 0.975]\n",
"--------------------------------------------------------------------------------------\n",
"Intercept -0.0959 0.179 -0.536 0.592 -0.447 0.255\n",
"prog[T.general] -0.2765 0.073 -3.808 0.000 -0.419 -0.134\n",
"prog[T.vocational] 0.1469 0.065 2.246 0.025 0.018 0.275\n",
"read 0.0013 0.004 0.310 0.757 -0.007 0.010\n",
"write -0.0043 0.004 -1.115 0.265 -0.012 0.003\n",
"science 0.0053 0.004 1.284 0.200 -0.003 0.013\n",
"--------------------------------------------------------------------------------------\n",
" motivation coef std err t P>|t| [0.025 0.975]\n",
"--------------------------------------------------------------------------------------\n",
"Intercept -0.9505 0.198 -4.811 0.000 -1.339 -0.563\n",
"prog[T.general] -0.3603 0.080 -4.492 0.000 -0.518 -0.203\n",
"prog[T.vocational] 0.2594 0.072 3.589 0.000 0.117 0.401\n",
"read 0.0097 0.005 2.074 0.038 0.001 0.019\n",
"write 0.0175 0.004 4.122 0.000 0.009 0.026\n",
"science -0.0090 0.005 -1.971 0.049 -0.018 -3.13e-05\n",
"======================================================================================\n"
]
}
],
"source": [
"print(res.summary())"
]
},
{
"cell_type": "markdown",
"id": "e3d3815b",
"metadata": {},
"source": [
"The the standard results attributes for the parameter estimates like `params`, `bse`, `tvalues` and `pvalues`, are two dimensional arrays or dataframes with explanatory variables (`exog`) in rows and dependend (`endog`) variables in columns."
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "4aa98cb0",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>0</th>\n",
" <th>1</th>\n",
" <th>2</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>Intercept</th>\n",
" <td>-1.496970</td>\n",
" <td>-0.095858</td>\n",
" <td>-0.950513</td>\n",
" </tr>\n",
" <tr>\n",
" <th>prog[T.general]</th>\n",
" <td>-0.127795</td>\n",
" <td>-0.276483</td>\n",
" <td>-0.360329</td>\n",
" </tr>\n",
" <tr>\n",
" <th>prog[T.vocational]</th>\n",
" <td>0.123875</td>\n",
" <td>0.146876</td>\n",
" <td>0.259367</td>\n",
" </tr>\n",
" <tr>\n",
" <th>read</th>\n",
" <td>0.012505</td>\n",
" <td>0.001308</td>\n",
" <td>0.009674</td>\n",
" </tr>\n",
" <tr>\n",
" <th>write</th>\n",
" <td>0.012145</td>\n",
" <td>-0.004293</td>\n",
" <td>0.017535</td>\n",
" </tr>\n",
" <tr>\n",
" <th>science</th>\n",
" <td>0.005761</td>\n",
" <td>0.005306</td>\n",
" <td>-0.009001</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" 0 1 2\n",
"Intercept -1.496970 -0.095858 -0.950513\n",
"prog[T.general] -0.127795 -0.276483 -0.360329\n",
"prog[T.vocational] 0.123875 0.146876 0.259367\n",
"read 0.012505 0.001308 0.009674\n",
"write 0.012145 -0.004293 0.017535\n",
"science 0.005761 0.005306 -0.009001"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.params"
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "0439bc39",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>0</th>\n",
" <th>1</th>\n",
" <th>2</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>Intercept</th>\n",
" <td>0.157499</td>\n",
" <td>0.178794</td>\n",
" <td>0.197563</td>\n",
" </tr>\n",
" <tr>\n",
" <th>prog[T.general]</th>\n",
" <td>0.063955</td>\n",
" <td>0.072602</td>\n",
" <td>0.080224</td>\n",
" </tr>\n",
" <tr>\n",
" <th>prog[T.vocational]</th>\n",
" <td>0.057607</td>\n",
" <td>0.065396</td>\n",
" <td>0.072261</td>\n",
" </tr>\n",
" <tr>\n",
" <th>read</th>\n",
" <td>0.003718</td>\n",
" <td>0.004220</td>\n",
" <td>0.004664</td>\n",
" </tr>\n",
" <tr>\n",
" <th>write</th>\n",
" <td>0.003391</td>\n",
" <td>0.003850</td>\n",
" <td>0.004254</td>\n",
" </tr>\n",
" <tr>\n",
" <th>science</th>\n",
" <td>0.003641</td>\n",
" <td>0.004133</td>\n",
" <td>0.004567</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" 0 1 2\n",
"Intercept 0.157499 0.178794 0.197563\n",
"prog[T.general] 0.063955 0.072602 0.080224\n",
"prog[T.vocational] 0.057607 0.065396 0.072261\n",
"read 0.003718 0.004220 0.004664\n",
"write 0.003391 0.003850 0.004254\n",
"science 0.003641 0.004133 0.004567"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.bse"
]
},
{
"cell_type": "code",
"execution_count": 9,
"id": "1d2be72e",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>0</th>\n",
" <th>1</th>\n",
" <th>2</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>Intercept</th>\n",
" <td>4.887129e-20</td>\n",
" <td>0.592066</td>\n",
" <td>0.000002</td>\n",
" </tr>\n",
" <tr>\n",
" <th>prog[T.general]</th>\n",
" <td>4.615006e-02</td>\n",
" <td>0.000155</td>\n",
" <td>0.000008</td>\n",
" </tr>\n",
" <tr>\n",
" <th>prog[T.vocational]</th>\n",
" <td>3.193055e-02</td>\n",
" <td>0.025075</td>\n",
" <td>0.000359</td>\n",
" </tr>\n",
" <tr>\n",
" <th>read</th>\n",
" <td>8.192738e-04</td>\n",
" <td>0.756801</td>\n",
" <td>0.038481</td>\n",
" </tr>\n",
" <tr>\n",
" <th>write</th>\n",
" <td>3.700449e-04</td>\n",
" <td>0.265214</td>\n",
" <td>0.000043</td>\n",
" </tr>\n",
" <tr>\n",
" <th>science</th>\n",
" <td>1.141093e-01</td>\n",
" <td>0.199765</td>\n",
" <td>0.049209</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" 0 1 2\n",
"Intercept 4.887129e-20 0.592066 0.000002\n",
"prog[T.general] 4.615006e-02 0.000155 0.000008\n",
"prog[T.vocational] 3.193055e-02 0.025075 0.000359\n",
"read 8.192738e-04 0.756801 0.038481\n",
"write 3.700449e-04 0.265214 0.000043\n",
"science 1.141093e-01 0.199765 0.049209"
]
},
"execution_count": 9,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.pvalues"
]
},
{
"cell_type": "markdown",
"id": "fb79f12d",
"metadata": {},
"source": [
"### Prediction\n",
"\n",
"\n",
"The regression model and its results instance have methods for prediction and residuals.\n",
"\n",
"Note, because the parameter estimates are the same as in the OLS estimate for individual endog, the predictions will also be the same between the MultivariateLS model and the set of individual OLS models."
]
},
{
"cell_type": "code",
"execution_count": 10,
"id": "5c884936",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>locus_of_control</th>\n",
" <th>self_concept</th>\n",
" <th>motivation</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>-0.781981</td>\n",
" <td>0.759249</td>\n",
" <td>0.575027</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>0.334075</td>\n",
" <td>0.015388</td>\n",
" <td>0.635811</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>1.342126</td>\n",
" <td>0.539488</td>\n",
" <td>0.432337</td>\n",
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <td>0.426497</td>\n",
" <td>-0.326337</td>\n",
" <td>-0.012298</td>\n",
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <td>-0.364810</td>\n",
" <td>-0.480846</td>\n",
" <td>1.255411</td>\n",
" </tr>\n",
" <tr>\n",
" <th>...</th>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" </tr>\n",
" <tr>\n",
" <th>595</th>\n",
" <td>-1.849566</td>\n",
" <td>0.920851</td>\n",
" <td>-0.318799</td>\n",
" </tr>\n",
" <tr>\n",
" <th>596</th>\n",
" <td>-1.278212</td>\n",
" <td>-1.080592</td>\n",
" <td>-0.031575</td>\n",
" </tr>\n",
" <tr>\n",
" <th>597</th>\n",
" <td>-1.060668</td>\n",
" <td>-1.065596</td>\n",
" <td>-1.577958</td>\n",
" </tr>\n",
" <tr>\n",
" <th>598</th>\n",
" <td>-0.661946</td>\n",
" <td>0.368192</td>\n",
" <td>0.132774</td>\n",
" </tr>\n",
" <tr>\n",
" <th>599</th>\n",
" <td>-0.129760</td>\n",
" <td>0.702698</td>\n",
" <td>1.020835</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"<p>600 rows × 3 columns</p>\n",
"</div>"
],
"text/plain": [
" locus_of_control self_concept motivation\n",
"0 -0.781981 0.759249 0.575027\n",
"1 0.334075 0.015388 0.635811\n",
"2 1.342126 0.539488 0.432337\n",
"3 0.426497 -0.326337 -0.012298\n",
"4 -0.364810 -0.480846 1.255411\n",
".. ... ... ...\n",
"595 -1.849566 0.920851 -0.318799\n",
"596 -1.278212 -1.080592 -0.031575\n",
"597 -1.060668 -1.065596 -1.577958\n",
"598 -0.661946 0.368192 0.132774\n",
"599 -0.129760 0.702698 1.020835\n",
"\n",
"[600 rows x 3 columns]"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.resid"
]
},
{
"cell_type": "code",
"execution_count": 11,
"id": "c3212cfc",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[-0.36197321, -0.03660735, -0.2060539 ],\n",
" [ 0.17005867, 0.09597616, -0.11549243],\n",
" [ 0.28641963, 0.09044546, 0.00450077],\n",
" ...,\n",
" [ 0.6252098 , -0.23716973, 0.11864199],\n",
" [-0.3024846 , -0.29586741, -0.47584179],\n",
" [ 0.77574136, 0.2878978 , 0.42480766]])"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.predict()"
]
},
{
"cell_type": "code",
"execution_count": 12,
"id": "69744e7d",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>0</th>\n",
" <th>1</th>\n",
" <th>2</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>-0.361973</td>\n",
" <td>-0.036607</td>\n",
" <td>-0.206054</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>0.170059</td>\n",
" <td>0.095976</td>\n",
" <td>-0.115492</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>0.286420</td>\n",
" <td>0.090445</td>\n",
" <td>0.004501</td>\n",
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <td>-0.058400</td>\n",
" <td>0.187809</td>\n",
" <td>0.007974</td>\n",
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <td>0.084621</td>\n",
" <td>0.028620</td>\n",
" <td>0.001513</td>\n",
" </tr>\n",
" <tr>\n",
" <th>...</th>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" </tr>\n",
" <tr>\n",
" <th>595</th>\n",
" <td>0.185458</td>\n",
" <td>0.036897</td>\n",
" <td>0.034498</td>\n",
" </tr>\n",
" <tr>\n",
" <th>596</th>\n",
" <td>0.330408</td>\n",
" <td>0.097329</td>\n",
" <td>0.489407</td>\n",
" </tr>\n",
" <tr>\n",
" <th>597</th>\n",
" <td>0.625210</td>\n",
" <td>-0.237170</td>\n",
" <td>0.118642</td>\n",
" </tr>\n",
" <tr>\n",
" <th>598</th>\n",
" <td>-0.302485</td>\n",
" <td>-0.295867</td>\n",
" <td>-0.475842</td>\n",
" </tr>\n",
" <tr>\n",
" <th>599</th>\n",
" <td>0.775741</td>\n",
" <td>0.287898</td>\n",
" <td>0.424808</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"<p>600 rows × 3 columns</p>\n",
"</div>"
],
"text/plain": [
" 0 1 2\n",
"0 -0.361973 -0.036607 -0.206054\n",
"1 0.170059 0.095976 -0.115492\n",
"2 0.286420 0.090445 0.004501\n",
"3 -0.058400 0.187809 0.007974\n",
"4 0.084621 0.028620 0.001513\n",
".. ... ... ...\n",
"595 0.185458 0.036897 0.034498\n",
"596 0.330408 0.097329 0.489407\n",
"597 0.625210 -0.237170 0.118642\n",
"598 -0.302485 -0.295867 -0.475842\n",
"599 0.775741 0.287898 0.424808\n",
"\n",
"[600 rows x 3 columns]"
]
},
"execution_count": 12,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.predict(data_mvreg)"
]
},
{
"cell_type": "code",
"execution_count": 13,
"id": "28f414f4",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>locus_of_control</th>\n",
" <th>self_concept</th>\n",
" <th>motivation</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>-0.361973</td>\n",
" <td>-0.036607</td>\n",
" <td>-0.206054</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>0.170059</td>\n",
" <td>0.095976</td>\n",
" <td>-0.115492</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>0.286420</td>\n",
" <td>0.090445</td>\n",
" <td>0.004501</td>\n",
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <td>-0.058400</td>\n",
" <td>0.187809</td>\n",
" <td>0.007974</td>\n",
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <td>0.084621</td>\n",
" <td>0.028620</td>\n",
" <td>0.001513</td>\n",
" </tr>\n",
" <tr>\n",
" <th>...</th>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" </tr>\n",
" <tr>\n",
" <th>595</th>\n",
" <td>0.185458</td>\n",
" <td>0.036897</td>\n",
" <td>0.034498</td>\n",
" </tr>\n",
" <tr>\n",
" <th>596</th>\n",
" <td>0.330408</td>\n",
" <td>0.097329</td>\n",
" <td>0.489407</td>\n",
" </tr>\n",
" <tr>\n",
" <th>597</th>\n",
" <td>0.625210</td>\n",
" <td>-0.237170</td>\n",
" <td>0.118642</td>\n",
" </tr>\n",
" <tr>\n",
" <th>598</th>\n",
" <td>-0.302485</td>\n",
" <td>-0.295867</td>\n",
" <td>-0.475842</td>\n",
" </tr>\n",
" <tr>\n",
" <th>599</th>\n",
" <td>0.775741</td>\n",
" <td>0.287898</td>\n",
" <td>0.424808</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"<p>600 rows × 3 columns</p>\n",
"</div>"
],
"text/plain": [
" locus_of_control self_concept motivation\n",
"0 -0.361973 -0.036607 -0.206054\n",
"1 0.170059 0.095976 -0.115492\n",
"2 0.286420 0.090445 0.004501\n",
"3 -0.058400 0.187809 0.007974\n",
"4 0.084621 0.028620 0.001513\n",
".. ... ... ...\n",
"595 0.185458 0.036897 0.034498\n",
"596 0.330408 0.097329 0.489407\n",
"597 0.625210 -0.237170 0.118642\n",
"598 -0.302485 -0.295867 -0.475842\n",
"599 0.775741 0.287898 0.424808\n",
"\n",
"[600 rows x 3 columns]"
]
},
"execution_count": 13,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.fittedvalues"
]
},
{
"cell_type": "markdown",
"id": "a37e6b56",
"metadata": {},
"source": [
"The predict methods can take user provided data for the explanatory variables, but currently are not able to automatically create sets of explanatory variables for interesting effects.\n",
"\n",
"In the following, we construct at dataframe that we can use to predict the conditional expectation of the dependent variables for each level of the categorical variable \"prog\" at the sample means of the continuous variables. "
]
},
{
"cell_type": "code",
"execution_count": 14,
"id": "59592e70",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>prog</th>\n",
" <th>read</th>\n",
" <th>write</th>\n",
" <th>science</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>academic</td>\n",
" <td>51.901833</td>\n",
" <td>52.384833</td>\n",
" <td>51.763333</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>vocational</td>\n",
" <td>51.901833</td>\n",
" <td>52.384833</td>\n",
" <td>51.763333</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>general</td>\n",
" <td>51.901833</td>\n",
" <td>52.384833</td>\n",
" <td>51.763333</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" prog read write science\n",
"0 academic 51.901833 52.384833 51.763333\n",
"1 vocational 51.901833 52.384833 51.763333\n",
"2 general 51.901833 52.384833 51.763333"
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"data_exog = data_mvreg[['read', 'write', 'science', 'prog']]\n",
"\n",
"ex = pd.DataFrame(data_exog[\"prog\"].unique(), columns=[\"prog\"])\n",
"mean_ex = data_mvreg[['read', 'write', 'science']].mean()\n",
"\n",
"ex[['read', 'write', 'science']] = mean_ex\n",
"ex"
]
},
{
"cell_type": "code",
"execution_count": 15,
"id": "d73c0ce4",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"predicted mean by 'prog':\n"
]
},
{
"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>locus_of_control</th>\n",
" <th>self_concept</th>\n",
" <th>motivation</th>\n",
" </tr>\n",
" <tr>\n",
" <th>prog</th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>academic</th>\n",
" <td>0.086493</td>\n",
" <td>0.021752</td>\n",
" <td>0.004209</td>\n",
" </tr>\n",
" <tr>\n",
" <th>vocational</th>\n",
" <td>0.210368</td>\n",
" <td>0.168628</td>\n",
" <td>0.263575</td>\n",
" </tr>\n",
" <tr>\n",
" <th>general</th>\n",
" <td>-0.041303</td>\n",
" <td>-0.254731</td>\n",
" <td>-0.356121</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" locus_of_control self_concept motivation\n",
"prog \n",
"academic 0.086493 0.021752 0.004209\n",
"vocational 0.210368 0.168628 0.263575\n",
"general -0.041303 -0.254731 -0.356121"
]
},
"execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"pred = res.predict(ex)\n",
"\n",
"pred.index = ex[\"prog\"]\n",
"pred.columns = res.fittedvalues.columns\n",
"print(\"predicted mean by 'prog':\")\n",
"pred"
]
},
{
"cell_type": "markdown",
"id": "68752e2f",
"metadata": {},
"source": [
"## Outlier-Influence\n",
"\n",
"This is currently in draft version. \n",
"`resid_distance` is a one dimensional residual measure based on Mahalanobis distance for each sample observation. \n",
"The hat matrix in the MultivariateLS model is the same as in OLS, the diagonal of the hat matrix is temporarily attached as `results._hat_matrix_diag`.\n",
"\n",
"Note, individual components of the multivariate dependent variable can be analyzed with OLS and are available in the corresponding post-estimation methods like `OLSInfluence`."
]
},
{
"cell_type": "code",
"execution_count": 16,
"id": "829d60f4",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([3.74332128, 0.95395412, 5.15221877, 0.82580531, 4.5260778 ])"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.resid_distance[:5]"
]
},
{
"cell_type": "code",
"execution_count": 17,
"id": "6d3c7fef",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[0.36844484, 0.05748939, 0.06050103],\n",
" [0.05748939, 0.4748153 , 0.13103368],\n",
" [0.06050103, 0.13103368, 0.57973305]])"
]
},
"execution_count": 17,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"res.cov_resid"
]
},
{
"cell_type": "code",
"execution_count": 18,
"id": "2ec68e92",
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt"
]
},
{
"cell_type": "code",
"execution_count": 19,
"id": "679a98f5",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"[<matplotlib.lines.Line2D at 0x1d9d826e2b0>]"
]
},
"execution_count": 19,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABguUlEQVR4nO19eZwcRdn/t3pmj+wmIdcCgRCSIPcNkUOUM0gA5VBRFH29XlFfb/yhIL6evIoXoqIgH0Q8ABUFFVAg3IdcCbkTQjhC7mRzH5vdnZmu3x/d1V1dXVVdPdM9O7Nb389nPzvT011XVz311Pd56ilCKYWFhYWFRfPBGegCWFhYWFhUByvALSwsLJoUVoBbWFhYNCmsALewsLBoUlgBbmFhYdGkKNYzs3HjxtFJkybVM0sLCwuLpsesWbM2UEq7xOt1FeCTJk3CzJkz65mlhYWFRdODEPKG7LqlUCwsLCyaFFaAW1hYWDQprAC3sLCwaFJYAW5hYWHRpLAC3MLCwqJJYQW4hYWFRZPCCnALCwuLJkWiACeE3EIIWU8IWSD57f8RQighZFw+xbOwyBYLVm3F3BVbBroYFhaZwEQDvxXAdPEiIWQfAGcCWJ5xmSwscsM7fvEUzv/l0wNdDAuLTJAowCmlTwDYJPnppwC+AsCeCGFhYWExAKiKAyeEnAdgFaV0bsblsbCwsLAwROpYKISQDgBXAXi74f2XArgUACZOnJg2OwsLCwsLBarRwPcDMBnAXELIMgATALxICNlTdjOl9CZK6VRK6dSurlgwLQsLCwuLKpFaA6eUzgewO/vuC/GplNINGZbLwsLCwiIBJm6EdwB4BsCBhJCVhJCP518sCwsLC4skJGrglNL3J/w+KbPSWFhYWFgYw+7EtLCwsGhSWAFuYWFh0aSwAtzCwsKiSWEFuIWFhUWTwgpwCwsLiyaFFeAWFhYWTQorwC0sLCyaFFaAW1hYWDQprAC3sLCwaFJYAW6hxbINOzFj0bqBLoaFhYUEqYNZWQwtnPrjxwAAy645d2ALYmFhEYPVwC0sLCyaFFaAW1hYWDQprAC3sLCwyAELVm3FpCvuw4JVW3PLwwpwCwsLixzwoG/8f2hxfk4AVoBbWFhY5AhK80vbCnALCwuLHEDqkIcV4BYWFhY5IkcF3ApwCwsLizxA6qCCmxxqfAshZD0hZAF37UeEkJcIIfMIIXcTQkblWkoLiwFEb6mCtVt7B7oYFs2KHElwEw38VgDThWszABxGKT0CwMsArsy4XBYNBpqnJabB8Ynfz8QJ3394oIth0QCYdMV9uOwvc4zuJXVgwRMFOKX0CQCbhGsPUkrL/tdnAUzIoWwWDYQhLL/x5NINA10EiwbCXS+uSnV/o3PgHwPwb9WPhJBLCSEzCSEzu7u7M8jOYiAwhOW3hUVVaAgOXAdCyFUAygBuU91DKb2JUjqVUjq1q6urluwsBhDuUFbBLSxqQJ5Dp+pohISQDwN4B4Az6FAmSIcIrAD37ACkHmqVxaBAPXpKVQKcEDIdwFcBnEIp7cm2SBaNCCu/vTaw8tsiLWiOBKSJG+EdAJ4BcCAhZCUh5OMArgcwAsAMQsgcQsiNuZXQoiFgBbhdhVikQz0m+0QNnFL6fsnl3+RQFosGRp5aRLPAtoBFNbCxUCwGHK6VXlYDt0iFethLrAC3MIIVXvWhkZau244HF67NPyOLuiHPbmPPxLQwgpXf9WmDM3/6BAB7BqmFGawGbmGEoeIpOn/lVpQqrvQ3uwqxqAaWA7cYcAwF2fVa9w688/qn8H/3LZb+bgW4RRo0/E5Mi6GDoSC8NuzoBwDlGYaDvwUs8sCA+oFbWABDwwuF0UQqzYnKmRULCykaIhqhhQUwtPzAVQNvKKxCLHKA5cAtBhpDQXYlVXEINIFFhrAcuEXDYEgIcL+OqoFnNXCLatDo8cAthgCGgvBK5MAHfxNYaJDWldaeSm/RMBgSAtz/L3LgTKDX0xfeHQpW4yZDta8/z35jBbiFEYaA/FZSKOxrPWXqUJgwBzssB27RMBgK8oR52igplDqaMStDocGbDNW+EbsT02LAMRTcCAMNPEaheN/rqoFbn/OGQ3oO3PqBNz2mXj0D//27mQNdjJoxFCjZgANXUSh1bASrgTceqtbAMy1FFFaA54wNO/rx0OJ1A10MY/zXLc9j0hX3xa4PBU62kQJ2DYX2bjakfSWWA7eoO554uVt6fSjIk1ADj448J6BQrBfKUEa1NKLlwC0GHI2kneaGgAMX4F+op0ytWAHecGjEIWByqPEthJD1hJAF3LUxhJAZhJCl/v/R+RbTYqDRgH03cyR6odRxBNeTA9/WW8Kh37gf/3llQ93yHEoY6GiEtwKYLly7AsDDlNL9ATzsf7cYxBgKnCxVaOAD4gdeRy+URau3YWd/Bdc9vFT6++sbduL255bXr0CDBA1xJial9AkAm4TL5wP4nf/5dwAuyLZYFo2GoeDWFm7kaYCdmA00YV74q6fxtbvnDw0aTYPqd2JmWw4e1XLge1BK1wCA/3931Y2EkEsJITMJITO7u+UGssGKwWSIGhJ+4P5/R7mRp36oJweuEjAVl2JXfwVbekra+4YK0o6BesRCyf1QY0rpTQBuAoCpU6cOqS5QGgRqKyHewB0KgzfUegUNHAPghTIADS4KnC/+eQ7umbs6+D4EuoAWjTgGqtXA1xFCxgOA/399dkUaPOgvN78AZ2jEzps1lLFQmBdKHV9nI3ih8MIbGCKeSBqkrX0j+4H/E8CH/c8fBvCPbIozuFCqmL/ycsXF9/+9GJt39udYovQIDXhDd/CyNqgnjdQA8juGRixTPVHtBDag0QgJIXcAeAbAgYSQlYSQjwO4BsCZhJClAM70v1sIKFXMVbaHFq/Drx9/Dd++Z2GOJUoPMgCbWAYOvhuhcJW1QT2boBHbeyjYQXRIrYHnUoooEjlwSun7FT+dkXFZBh3SUChMWy/VoOYs39iD3TpasNuwlqrTEBFqn4MfyeFk6+gH3oDqbgPOKXVF1V4o2RYjArsTM0f0p9DA5eazdDj5R4/inJ89WUMKcQyEC91AQXWgQ/B7HZqAtXcjCfCwD8R/27ijDzc/+dqQ6B9pUQ8/8Ny9UIYy0lAoWQ2AVVt2ZZIOgyfM6JDQvlgdHVGtCbbS598IXms3prYro1C+fOdcPLakG8dNHoMjJoyqf6HqiZTvhI3pRvQDtzBAqZz+zdVj1k6FAYgDkheSJkkmoGPxwIPf8yhVFCxwViOGk5UVaesuz0c8jcG+WZHWBlCPFrECPEekolBUgZQaBI1oVEuLpCoEP8fcCOtnCQgEeAPOmIOhD9SCtNVn9w90LBSLKpHGiJkUSGmgEIiuQTB2k6oQnEovXCd1XIWQOtI1IhLbR3ZtEPQLU6StqtXAmxyMAy+q9mZzaFQNfDAZMY0pFDEWCvu9DhI83DTUeO2ta75GUzzyQNoxYDnwJkcgwAspBHiDjQTGBzeeOEmPpDok7bSsi0Y1ABy4cZfTFGkQzO+JaMQqWgGeI5gAb4m5NcSRhRthHhjIJX3WSKoCE5rigqmem5kKLK8GisLAqj0Y+kAtqJ4Dzw9WgOeIft8yXzDSwOXL90bBYBi7ScYkFQce/p5xgSQI/MAbsMF1JRrIbnv37JWYdMV9mbvQ1op67Fy1AjxHJAmE6L3e/0aT3wNxHmReSKqCq6Cx6mnIDbT9BggnG7+vMfvA3bO9oFtL123PNZ/UboQ0+j8PWAGeI9K8uMALJaeyqHDZn+fgvb9+Jl4eYfJp0LGbKVzFhFtPGskZQMoqqe/J5pRG6BZ166NpKZR8ShHBkNyJub23BEIIhrflW/00M/ZAaeB3zV4lve5SgGd+BkMgI1MNXHWqcX0EeOP6gTdqHwg8pXIuX2o3wuAB6weeKQ7/1oM47JsP5J5POg3cw19mrsSOvnIu5UmDYLk8ALGw84I5B644Ui2XUkXR0NEfJUVKq29QSnH/grUop9jklgQ26eXdZKmNmJYDb26kUaL4zvHNfwx8SFk3Kr8bVPdKh0QN3JVvpgqX6Nm1wrbeEt7762ewfGNP5HrBH5EZyrfMkAWF8vDi9fjUH2fh+kdfyaRMQP1CHVgOfIghzYDnO8fmnoE/1CHcGdrAGmFKJNWgIkxaseczbIIHFqzF869vwnUPvxy5PhDHt5lCJ8BMNfGNO/sAAKsz9Bip12azBnwlzS3AD//WA/jdf5YNdDGUSPO+G61ziJx8o3ogpEFSHUJXzuj1PLbSh3R7NLOBNGImoQGLBIA7cCPnfNJz4HYnphKUUmzvLeOb/xx4ukGJFI78/D2NICzFIjRAkWpG4k5MZTTC/LTi+GTRuEbMLOufZVPmQXHJkH4rfU4F4dDEAnygS5CMdF4o4b2NUDXRpa4B5UlqmPuBR6/rDjSovjDyy07Agasz29FXxgvLNmVYGDNkUX/VYRk1pekn+exrmzDpivuwdmtv5nlUAxr8t14oMTTiElNEKi+UBquOGNipUV3IUiFRgNfPiKny+w9iz2iy+uztL+KiG5/BlgawlTRCx2VeKL97ZhkA5Da5VbuVPk/UJMAJIV8ihCwkhCwghNxBCGnPqmBJaMStxiLSceCNVR/RJXpQaOCJboTsk+hGyCiU7MskThaMA9f174WrtwFIF644C+iUpoEMAZHLCikDsP7WkBw4IWRvAJ8HMJVSehiAAoCLsypYEhrtZcmQZpUQ5cCzL0taUF82DC4jpv53Rluoov9muQoJ3ewFI6bTuBy4rv1S88M1loVHHrSMDINOA4e3k3MYIaQIoAPA6tqLZAZLoeQLKujgjVa+amBsxFTIg1y8UGIaeCO7EWaAPGStwmaRNao9Ui3PN1m1AKeUrgLwYwDLAawBsJVS+qB4HyHkUkLITELIzO7u7upLKqABFZQYUlEoVT6XF0SDXiMKlLRIPtBBfj2PVYgqqXoeHpEWjdoH6kXepK6+/0Ce7VYLhTIawPkAJgPYC0AnIeSD4n2U0psopVMppVO7urqqL6mARu1MEaTwA414oTRA3cQyNECRakaSqyYNBpzi+Uw1cLm2Hx7okF1eWUFLodSvGDHUi39PLb+rfTAFaqFQpgF4nVLaTSktAbgLwFuyKVYyaANuNRZRz06dtdAXt9I3xYSZgGQ3QvmEm+8qRB53pa7hZI17qvq+geweBicWZoJq/cAbUgOHR52cQAjpIN4UeAaAxdkUKxm8lf7aGS9r7hw41JMDz7qPiLsS6zFAZy/fnGmQIxG8oJLVx1WoTOFGngzLovIDT3GkWtavJEmRzdSImWHhVa6YA43ACyXHPGrhwJ8D8FcALwKY76d1U0blSgQ/q/384aX1yjYV0nTqWmfprDtJaMKsjx/4otXbcOGv/oMfPbgkv0yo9GOAJA08Wz/waNoMJht5ModhVroimZY2Fxtmo1IogQaeeVEC1BQQm1L6TQDfzKgsqdAMS/q6jsHMKZSoBp53XdimlLkrtuSWR6IXiqsQ4Oz5HNpAFD1haNTkzLISW+YEisRuwP4P4HBU7ZxNi007+9G9vQ8H7jlC+nv6cLLsucakUAYUTSC/q/dCqaJumWvgVP8dAHr6y5lRHq1FryvmuTmFr4NsULFJShRUuURkVAbOYn7gBklkW5Sa7htIw3tWlMk5P3sSZ133hOaO6miihtzIM9BoBg08VTjZyPI+fd2ybg7Wvjq/5EO+8QA+/6fZmeTXViwAAPryFOA8By75XUWh6J6pviwe4oGzPJhw4FmNAdP+pstvYL1Qskln7TZ9DJX0GjiN/M8DTSzABy7vWW9swrfvyTYKYq0vOetOEtfA5en/a/7aTPJrKXqjsJSnEZPKP4vXxJ/y8MRJmiRMvFCyKk4WGngtbbOtt1STBi96oeTFiKcuIePAc/SYa14BPoAS/N03PIPfPr0s8b5wCeV9WL6xRxnIvvG8UITv2SYfA9NEc6VQEn5XauB5xANX7fpMsfkjMwplABJhCseKTT044lsP4tYa4vo37FZ6///9C9fmFmCreQV4A1AoiQcECL365B89irdc84hBujUVKxPEjJg5T5isrZZt7MkxmhxPocg4cBq7D+A0uhxejEr0GFEoGb0TU+03Dwpl+SbvSLkHF66rMoX6HQSe/ki18P6nlm7IujgAmlqAD3QJksdzOj/wBqgQh4CjrdOBvvwy8/pHsjsvkUcShcLYm5gCnkM0QlVSaSiUrGDshaI1YpqlIbr8ZSF7s3YjzGos8sm0teQjaptYgA+8wEsqQSovlAThkub5LCCeTlPPA2OHt9Xk3arOI3HC1WvgeXDgovBh1428UOrMoehuS902MYqu+spk5UbIoPLBr6WK7b6RPms0rQBvBI01qdOm6dR8n6nKCyUnI+ZAHBjb0ZpPZ09qo4BCUT1fFw08mQPPekKppxdKTLZmsNM3awZFRV+lFuBWA1dj0FEotXqhZG7E9BIMtnbnzYFzyXfWQQPXbaVX7cTMVgOXGzHD3XvJwjIzBdxUA5fcJxrq0yLc6ZuMj/z2edwzNx6x2smYQlFq4KnDyYb3t1kNPIrGoFCyX1JX/Xw2xYilx8ZG3icg1YVCSfi9nho4g+hBEVIoyZnlwdUm3FnFL/r70sjex5Z043N3zI5dz9qIWc6KQuHub7caeBSNcGKJMadqklaKdHV5ZYVd/RVs3VVCgRnwhPbOOr8IhdKWE4WS4IWi0ib5eDA7+sqYv3JrLuXz8vBgRqFkm2cStPkZJqKsV6YUSm0SvZJDLN+WQj6iNh9Vpw5oAAU8UwoFEeFSRVmqeEaHi258Bv0VFwfsMRxA3KiWl9EUADpa8uLAuc+S8jOlQKUluhT41B9m4alXNuCl705Hew3lFG0MInSbP2rrKbKymKWTxUYepYtmDcjaCyU7Dpx7UznJq6bVwMUOU6q4WLm5J9M8lm/s0Wr6iYfkGuZDKcV989ekKJksjZoej6Hfl9ihF0o0g6wpLD61vKLLRThwye8qP3D+91lvbI7cWyvi+3i8dI3CyWasgSdtiNEJetOysOF01+xVmLNiS/Cukw+crp/Glh0Hrv6WFZpYgEe/f/OfC/HWHzwaRLWrFSs29eDkHz2Kn2jCmyYtYU373BNLN+DV7p3cg2bPRTOr4hkDqAx4mXPuXPr52TeShET0P0MgZBKMoOlKok9A5weeOYWSUvhGng3ifZimEd5569OvG8eb19W1Xhx4WvB1yovxbWIBHm2Rx5d4521u7y1nkv46P7DNs69tVN6TfMai2VvbuqtkXjBVWXKS4CqjWvZeL/mlLc8jnokyHjj3DGvnWo26KgolaG8jLxT1Pbc89Tr+MWeVaWkM79Jp4GZp8N3IISQM4Zv4nPoO0QulZj9wgQNfsakHu/orVQezAvJTSpqYA893ScVS17koJZWAxj7IUczgTKjchJ5CYGXud8595jv7jY+/ij888waevuL0TPPQUiiKukUmgBpDtgS0hbiRh7V3jcGsvnPvIgDA+UftrU1j444+fOqPLybm5RcuhrQbvfhxSwgxfs7EqBvmYZamCmJff9sPHwUAfOMdh6RKpx5KSdMK8LydUJhmoJvNEwexwgulv+wG8a+BuACvbiNPPgj8knPWwPn0+bSv+fdLmeVhukxX+4GHAiszDVxxXZd8lht55qXwqNFRKMZafEQDN6+DzqjriOOnxnapKDJjk6IJ/jZrJe6dF9q18tLAm5ZCyX1jif9fZ9gxNWKKHeqE7z8c+V4syH2B0yDvFUneXihR7TgnOgjySSK8Jp9wZXaAWvtfUBaRQkmR/uadJdw5c0VN5UgD3XvZ3FPCLx99JTGGC9+GDiGJhmPZcyJik6A2pWSUI8pEdal9+c65EWrUauACxBeatSGDhhI8+Z6E38XbNu2MGloLTu3zaG4auP8/bsTMNkc+/bzm5mUb9IbiwI1Qkb+K5qkFKgXBRMP/7B0vYktPCVMnjcHkcZ2ZlEcHXZGuvGs+AOCEKWNw7L5jlPdFOHBH7bopQtseQhPW+m7KFV6A15RUmE4jeqEQQkYRQv5KCHmJELKYEHJiVgVLQt5eRazBdfR0UkcxiWsBNDgHztzaBKmauZDl0stjufmv+WsSuV7VqfT8Bh/Td5oE5STh/2ASjXBLj6fh1RJDPY1gMblzWIteJxQ5cBPKCNDTlaKdqppXw/dv/nNWfTGvQx1qVf1+BuB+SulBAI4EsLj2Iukxe/lm/OLhpcqGzSzAT7DC1RsxV2zqwc1PvqZNg9J04UHT1ODu2SvxyvodOdIOHmJGzBz9wPOYjBaujnK9+njg0euBHSBDCoUh5oUilKWaNPLC3BVbEt10RTpQRJRCMW/HdIHh0r8b/iQovq/Xe7drWlQtwAkhIwGcDOA3AEAp7aeUbsmoXEpc+Kv/4CczXo41bBCzI2MfTq0RkwIf/M1zuPq+xTFaBODdvfS+pbWU+Ut/notp1z6eXw/xETNiZp0+N2hkk0OeW/fFayqfd0o5r4tavVASeJocdnPXjGtnvIz3/fpZ7T3JETrDzxEOPCFvrVul8NPds1el7i8RAZ6HBt6ARswpALoB/JYQMpsQcjMhJEbEEUIuJYTMJITM7O7uriG7KGIceMZxq5mGphfgNPA71/sVU+0LzEK7zW28M4FSRz9w2TusdWKOBY2S3KMSJgGtwf1QqxcKS0tlgKvXgQ5pq7Fk3Xbt80nviR8HhP9egxFTXE09tqQbMxalO+GnxM2YuXDgDSjAiwCOAXADpfRoADsBXCHeRCm9iVI6lVI6taurq4bsoli1WX62ZNYUSpIfOOE+S2/w09Jp4OLLraYG+fmBexCLz5e5p7+s3fCUJh8v7fjveUdDBDRnYvrIkgNngi6+kYdGfm82qFYmr6zfASDatrwfeFJttSseycPbUm7o4+0IfNvXSk2a7jStFrUI8JUAVlJKn/O//xWeQK8Lvv73BdLrWXX88ExIjQCn+hfEUyi6cmVxEHvt8cTlz4fap1oDv/yv83DxTc8qD2w2QdQLRUahVJ20FNIVk+JINZk3Ua0aMk2YLOoVLjnrbGTlfnDhWky79nH8a/6aSLs5hGTCgct+SWsW6CtXgs9ZcuBBNM9GE+CU0rUAVhBCDvQvnQHA3NM9Y2QddD80YqqR6IXCUSh6AV57mWutdtLzMQqF+/zSmm0AgJ19NYQxoNKPyvxrhZZCURnIeU+FjCiU2GQhlMUEdZL1RpC1y+I12/3/2yL1LTjJqx4GrQDPoAGiGnj4uVZ5wjYZNSIHDgCfA3AbIWQegKMAfK/mEtWIrNx1jDhwyf2R3/1LLtULIJ12a4pau4fqeXZdx9OHEeWqR5IRsx4aqSqPkDbhdmLWOKEkebxkPWFVXIq/z16FchbLPQ0opaCU4nXO5z4YSxA4cM6ImQRteF1JEmk9c/o4AR7hwGtsLuYhnFfvrUmAU0rn+Pz2EZTSCyilm7MqWLXIiitlHUbLgfN5aSgUwIvpocyryjLTBKFXbVrR697/vlIFG3f0Bdd5+ULEe8uV1IIiasRU0xtpsWDV1oB/VeUX5iv/LaRQuKV1jQM7zEtc2fgcuOZ1il3SpP/8+olX8cU/z4mFLeafzMIdseICtz+/HKf9+DG8sGxT9Ech9gkhIX1oel6pDLJf1mztxXo/IJ0JIhRKhl4ojEJpRCNmQ4H1vczcCIV0pffQ8A5Zrvw7+81TryvT0dETOkSMLbVSKMrr3i8PLV6PY69+KHYd4OwA/rUDv34/3n3Df6rOXy5cq6vgO37xFKZd+7jE31q2YqLS32jwO6+NZzNhxigUtmrTGr2j303K8pgfrXOYcAhF1oKl4lLMXr4FAPA6HyJZkl90K70+XT2FEr/2oweW4LjvPRz/QQFeA49y4BlRKDmR4INGgDOk7ZCPv9yN+xesjV1XeQlE8wo/S7VG0+WhSuVLfM7oturKkFSUiAYeb6S5KY8dS9pKH48Q9wh++7R6UqwGbiA8o9dlhtzaOXA/TYV7ZqqNK4ry8ljs2yn4IGpA9kt7j0Lxv0gM/OJGHlPBFjoVSPLMoBZ9Ci+UWseYkwG9qE0/p3QHDGk18A/f8jw+9cdZseuu2AsliMb7jf726Evrcet/lhmVIRYoyugpkTc2fEgBteuc4nqKe9PmL/cQiV5bsWkXvn1PDTZzjZav08DDjTy1CnBlMQBE+3FPfzmSXxKFIisaE1B52xIqNHS1jMXpFsrm8G6EiRp4mIaILKrUr+LAaxS9hcCIWVMySgwaAc4MaWk0I522zvy2dWFKKOW8X4Q3dEsK7bBaYRCl4PPpIeKEuGj1NuzoKwvGqHh50oNKPoXgi7FU2ExSW24hlMGsJFpxrVQdn1d/2cXn7piNNzaGlAPLq1Rxccg3HoiEMhXLJ/Z5rcE8YxtmnM6Je3DRyO98vyHG4zVcEccHZNKqe8WmHizfqD9uUaWB1zo5OJYDNwN7rb2livY+HroxaBIPPCsNuFqtqB4auFi2c37+JD53+4sKv/fqC8G/C9mExg/0M3/6RNX5MEjLn6AVuzTqkVJb/iEt8+xrG3HP3NW46u4FsY08rD/rwsaKwkEU4PNXbg00zCw8nnRwXRqUJxg73Hc+P4eobQGxdFkakt+Snn3bDx/FyT96NPj+wrJNsTbqK6n8wGtrICcT5UaTfj7JZo/5K7fij8++kXjfx26daXwupm5WDGb8hGBWDLW86LiLntlzEQFede7xtJKuL1wd9eeVnRmZFlGeVFKOlBLz0ZfWY+6KLcF3EwcLVyHBZa1c68CWVYeQKF0DhP2QP7RAVCpECk7sT++8/iku33wpFJeG03ho3Pa/gyg38iRpqLqf01TpsSXrcdGNz8TozX6uEfld02kn6pWbe4LjGAGeQsmn3ZsmHjjfCZOwcWc/RnW0Jt6nezlsEOhCdXucaPT+alAthcI/ljwAKBav2Y5D9hop/90gD4Yxna1RTjbIQ1sELaKTUTwh08D/DB+99QXt7/zRZT95cAk+8bYpag48MDiG17LzA5dP3qw/MWFSdAgmXXEfvjht/ziFIpRFv+dAvJKtYKm4oRFTpvzwuXkn8pin6z0jS9O8DgtXe8bc9dujLoZ9JZ4D5zbypHzPb/3Bo5HvrLyWA08C9151vts8dC9ex7kFz0cojHRvKLKrL+ZGaJZWJP+Ee297bjnO+fmTePqVDYlpRcopuT6mMzo5BpEga5Dg/JOyZAIPkYwGAsvj4cXr8KvHXsW371kYUigKDdylNLMj1fj6ROrOTSxAGCWP9enrHloqqUu0LDqhk3eMFZkGzvDTh17GrDc2R+9Pu5Ve5oWSokrd2729DF3D23Dbc2/gG//wQnLwGnhUMTJPWwamADbkgQ6NCtPzEXQvJ6RQNM8jel5iGvDLNHHTxobt/dLNJyLSdDQWbGoDtxmHh1IDl1RsdGdrJL8sQvnyQmhHXzm2CYOlnfVSNOSaXU4DF8uG4DobiFltnIpr+95/1u7MI4IXhjEKRSiLLnBa1u0npuZS3uUvTq3N4Wgt/t5kLxT1eDStEaU06P9jOltx1d0L8PtnPFqW18BZ22/dVcL8VencYUUUMqAXdRikAtxQA9c0qnEwK1THcUUt3dFn127r9WJ8JyCapz5/Fq9c1J7DMqjyiF8b3dES3chTZRuo8v/rrJWxTRhZR+mTaYmVQJjIV0RRL5Tw9519ZTy8OF340sAPnMpXUqyaoWao7odik+jeQ95GzAq3pEhSpFxKwzZP3Inp/Zf6gRvWoa/sYr2vgYsyQhbM6kO/eQ6fuV1/ihPDsd+dgf/+3czYdbuRxxDV7ALWxug20MAjm09SumeVXbnbUhqk8UJhAryoIPXTUCgFYUs0QzX1eNevnsbB/3t/4gAOhWvqLDzEjt2KC2n2SsRqyKgVvq5fu3s+Pv67mXhl/Xa8vG47zrv+KexICOwVbCGPrGTCMrL6lgxCEojCIY0bYa1iRRwfLqWctsw2saj6Vlj2RA1c27fMarGzr4yt/jF0Yhv1l93gaENW/nkpNqNt3NmPhySTuN3IUwVMNUHdXaHRxPv8qT/MivF3Uc+J6jVwEz61VHFx2DcfwN9nr5Lmn5QCE+Bpd1zK7i9zrmKA2hfeBC8u34JdpUriAM6LQglCIXCasCoHvs58OVjgph19Ffzw/iWYt3Kr0tYQT8uMQtFBbBNt5MuMNfA4hcIZMRNc6LyDTszyqQiTQjQdszT4/Qtim/WV3SDMQJbachhO1mrgWrzazW+CiP++vbeEa2e8HFkq8Y06b+WWyKDjjZgbdvTh/oVr8ck/RHds8ppF2vcTcVUy6DDd2/uwo6+M7/87PHa0Gg1cNbhVj0td+mhUp0ryxKGUJm6+SWoC0TujVty/YC3unLki4uoWlEEh5Fwqj0YYeuFQ401N4VZ63s2OlSSsZ0ihhAmu2xa1Y6TxQokZPDO3KYTjImlVTKn5qfTipCD7LQm8ABefKbsu2lo8cZgl21GtjcwUg0aA85B14FueWoafP7wUf3x2eXCNf4nnXf80Lrn5ueB7OOOHHVGMd01p9XHI08Zb2NbrLf06W0PPT/65NVv1hymUJRqsiReNdFMN5yoWuVdRjztnrsSZP30CTy1Va6VJ2+eDga5gFPrLLuat3KJMX8TV9y3G5X+dF80vQQP3Ji6JBsctk0P5Ek3lE7+fiU9zIRtY1Sho5NZgomJeKGUTDlwQ4FoOXP9srZBq4Np7zVZWeiOmWR129JaVikDFpQG9mKWnDiHEawergZtDJgxGdbQACA8fUN3HEAgPApT8z7uEXZ61UChliXDSYdMOT4PubOMEOPfcR36r93sOnlFo7aoimFAobKSqVhLsRPhX1qu1cFnzRV279AP96vsW4bzrn8ayDfEIeEA01gUPh/MSCISqkIWMWlGVQ6WBz1i0Dv/mgqbJNEF+Iw/rE6yf6HcEC99TuBFmLVdcN2re1uXh0lCQJq1C9Vvp1c/x/XTF5l0Rl9Bo+t6mG0KyndSYAmg18BSQaSAj2j3B9xofaF7TqGEsFKKMba0LZpUEvsO63NI7kj5XwE09TIAXuN/T5QlEvSf4x3WGptg1V0GhiIKd/S5sZnjutY34j8ARywZNSSLAVdol23W5ZVdJ+rsqHjtPfySHk+UnXWly2p27kTT9pOIGU7+e/g/9BkZMUSinOX9VFtOmFnheNWbp8Ry4rszevX6ast80z5VdiuG+0jN3xRZO4xfLTVFwiG+gz1ID92SI9QNPAd7SPnPZJryyfkeggW3mttlr3a24N1xSGJI8TtQDpRSPLVmv1ABFiBp4QdLb+XwZhz2c18Cr6GjRsK382t08jbJAofAbeWRFEpfS77vpWXyAo6uAMJhQgfM96y+7MQ5RaYRleZlWQlI2ladLaFiURyOU7URNak5+16eMM45TKMlpAV7slLs5Q7cIUdhnrRlW+IkwIW2X0qAdk47jC8qdkgP/v/sWBx5BL63dFvYjyaRXdIinsLkUl/15jr7whmAC3GrgKcB36Pfc+AymXft4oM0VOQGha9NgaUdpxOWPR1Qjo/jIb1/AqT9+zKiM/Ll7FUojsS4YeA10o0+hDGutUYBH/M/D62lS4rlgMW1ZmUwisvX0e/QUP5H1V9zYQR2qJLJQmlQubcHmHcg5cFbkcsVNbcSM5+WnJVAoJmkBwLUzXsZNT7ymuVf9bFps2NGHXf2iXYhydZOvaMJ7w/w395Rw9HcexH9eldtJ2H19JRfzBfc+nXbLxzzpryhsGPDeveMQOI431u7STIJpQECAjGkZHoNegDP0V0JKhEHXprzxX+XKRRGnB0xRFgSpTAPn82UaCi+Ab39ueeyZJPAUhOhFw3xkk1Cu0Mgqh/dCkVEcTqBFqxuJ2Rd4N3VPA49qvEn2giVrqws1yx9EoA7sFX6W1bPsyqkwGUI/8OiqJWKXcKmhH3j4efNOfSC3+Eae9IKFrQanXv0QlglhWiscvRZ69ajKEqWiNveU8IP7lyjvBbxJ/Z3XPxUYrH/79OuJLpsMlIb9Vua5UyAehWLS5qbwNHCk05BSoGYBTggpEEJmE0LuzaJAWUA2xsWYEoC+8zIN2aVUyUPWEguFF84Vl0p3rfH5MoHRV65gR18ZFZfiZs0xbSqoYh27lOKtP3zEKA1RA+/1tyGL3il3vbgSTy7tNgrowzQ5fiIrcRp4MoXiXf/K3+YZhSFg4N3YknzkReHKwMrID/wkzpO9W/4uQkikH5VdqjS+8uAnE97ILYPJ4Q9J0B2Xx/cBE88S8Z62olwkiffdO8872/Pb9yyKuVXqysYgFo2tgh1CjNrcFATRo+OyRhYa+BcALE68q46QaWklCceqc7diY9HVaeAJGplpGSuunELhqRuWfF/Z29Bz5V3zYveb5KXS9iiA7b16HjIsV1RQL/I9e8QBedlf5uJDv3k+kHC825gIRqHw7dDHceCsfa+d8bL0eT7Z5ZvM7BCsTOx5pRdKkEf0nTGwVUK5QhND67KDfnvYiopGDXT8Yx59R4PflHXgyjKstaC5M258rUawvK6x81AKTgMPJ0ceHX4ZqUSAt7fIyy/elybuf5hG+A5lFErB8fpfX4YCHIQ0rhcKIWQCgHMB3JxNcbKBTEjItGjdLjfebY33Qnm1O9Tu+KdNdsxF8ha8UAoyDrwcFxi7fEH3l5krjfPqV502EqFQzMuv8gOvuPLJk3fV6y3JBwerV8yIGWyccXHSNY/grhfl3CRfno/dGo9JoQLfLdThZKO/A8C37lkU09TKLr9ikLfnRTc+AyCcsGLtzn0tm1Io3DOdCQJcJriyhMyIKdZxr1HDMLqjBS71ysMb5lUauIzyqCYCKEsm5gdOfQrFIVVNDiowDTwnBbzmeODXAfgKgBGqGwghlwK4FAAmTpxYY3ZmkPVJJsD5jqAbHLzRrMQ983buNBi+8/M7PE3A5630QnHjFIroi26CfuWJ2+E9aTqYx3PGH+AHCI/wVBKK7b1ynl1qxCy7geq5q9/Fqi36zUrVIAympDYsMqkqXt5VqkQOCf7Zw6/g9Q1m9E2PTxmJafLfKy5VekDx4N8pb+SWIe5GmJh8Kjz9ygYsXrPdT1ueuMO51lVcT+tmniIqDVxMyqVmbSM+E0aTjP5WcT0lyiH6FYY0XU0jEoLMfct5VK2BE0LeAWA9pTR+IjAHSulNlNKplNKpXV1d1WbH5Zt8j5xCobHfdB2A37lYrqg02LA8aZdd/CzvUvkGhcjhqjUIcH5yqYW3Z1Bq4AqKhOfAtyvcxUIjJs+BUym/LEO1wyP0PAknn5hQDTTw6HUx4NniNduUKwwROzkNnCXLukBLIdyuz+qt6/d8m4sabPy4teizacIxmODJpRuCkK0queYQAuK71jHqgqGt6Mh35cbqQYNVmylc7h3zQveWp14PBHhPfwVrtvYqUpBDR58SxG0bWaIWCuUkAOcRQpYB+BOA0wkhf8ykVBrINFURssZiA4FvbJ1QcDkBrhL0fDZ9hgOXoVeIP1yQvAlRSweAXoNOW6q4kbjL0QNbw/tEvtUUvKeBeF02ebI39sdn38BWxUabXVI3wkoguGSGJb4rVDtAohx4KMx5sG/iddUGL5aeDowD5/MFvIHOaKSy6xr5gYtKBY8kt8G8uFmA35QUvU4ICc7DdGl09fnXWSsx+cp/JUZYrLjplRne3sDLge/cu8h3JKjOA0XnGRXUNXWqZqhagFNKr6SUTqCUTgJwMYBHKKUfzKxkHNipGYBnpDn7sD2198tmxFJaCiXgwL1z7mTgxVhaCqVXOERVvpEnLnhNOu2PH1yCC375NBb5x0fJvFkAROKKpKJQJMYnIKrh8GCri/Xb+/AHP4C+CGawVXHgvZL2zWDjIOeFEq4q4hq4fNndr9GO9SFdKXpKTAMXfPMRhvx13ZC+0wla/jdR8MUFtjA55UXOQm47AEIKxXW9viRbfc4R4trIKJT0ApyPvRL/reB7oaRFkvLT6F4ouUNs0vcfp+fSpRy4r8nwnh3ak0v837q39+Hq++RONpRmRKEovFBKVVIoTHB3+0tZ2WkjgGDENC04dEZM/UYeQH0iUBC6gOuRvBeKbIXDD/xqxwcfFTBMTH8vQ3BajmQq0S2re8uViHATbRHFAqeB+5OETttfv60Xv3h4KVzJCii+81L/PUsE4Q+EMgxrKQS8MKWQGvDnLN+iLWe1FArrZ+LE1dNfQcEhVWnKSaELvLpWkbABMhHglNLHKKXvyCItGUThljRLspfNH8sVGDE5oShbnrIXy17Kxp1qH9MIhZJSgO+KcOByL5RfPBKef8gEggnHGnLO3jOyoFCA6FKo7mFFoWzegIzfrxbg4We+nfg82WQV9QMPOXCZF1E2Grj/P2LclVMofKQ9r3yusiA6DbyHEzxUyI+ChgcLuOGEq/Ny+vUTr+EnM17GzDc2B2mxVWrcWBf9LivmXS+uxKf+oDVtGWnuKo+PSeM6fSOmeg9E2XVx3UMvBwHQZF4ju0pmbq9BeTgvKfH9bOkpeeOmCkGrM2JSMA48fbomaAoNnBcgnltO+NuksR2x+1mD8sdySTlwScOLnU43cKKbWdJSKDytIZ+Unlwaj09uAt7rAwD6eLrGd796dMl6pUeKCNG3eOHqbfjgzc/H7lNRKFHDpHwFxDZOOSKF4rdLn6R9Ixx4lSyj7Bi1OIUS/o8KcHWeOq1sB+dvz/tCE+Jr4P4ypOy6qYJabd1VCtI6YA/PMSyJMpFNuJf9ZS7uX7g2dp2HSX8MVzfReyeP6ww0cJfKV58bd/bjuoeWBv0sRnm4FLv60ylNfFgMMb3NPf2+Bp6+H+njryPg+/NAUwhw8QUTQUsTIWtPdl+EA5dozcHJL/5/3dFYLg2Xz2kF+C6RQlGok7wx1RSBBu5XT9TA/zFnNT762xdwGxcb/WO3qsPRDpO4dsmonIqr10YAQYBz764s0cArnF+1bIVjGvlPB9kWfbGpeQ48rRFc/AyEkSVZmnysa4qQUqi4YWgCk5goO/vKgYbNPFmSOfDEZKUw2bim2jQzmWngft0dQiK+4EBImQUCV0INVcOBh+7BCgqlivbQtQWF10+HNAfOD2riO9szyGY/qRuhxIgpC1Il8nY6Ac6fwJKWQhG1YhmFAoS7I9OcuUlECqXMC3BgtX/4w/JNoXGW/yyiI2FzSJh2MgceEdoSe0Ts3eraN+KFYlTEGGRLalU4WbF+rPyyN1eO9LNoerwnDkX83QZuhJRKJxgV+BNnCswQmkShcDeksd+ZtLdLPe8iMc89RrbBYRq4643vR758SuQe5hTA+oNsIkorwCsaCgVA1aFkte+GUl8DT52sEZpCgIscLP9VppnIlit9gREz/K1for2XhResa3j+J97IdvmdcyP0hwyiF4qK19+yqz+4xxRh8Cjvf0mIu5L2FHnV5goRqo08PPqVGng8Vk0llN/YKAnSxLeYLtvjJ49R/hZM2ArbAJ8448BZev1ldR+JRJsUGoUPnRqb9GgotMoVeXAwFbb09AcCOdDANUbM17p3YDbnbpoGJhPK4y934+Bv3B+ED2AY09kW+IFXfCVo95HtkZUeM1AWlRMRNXKpjT7jhn7gkuKnNWL+4uGleOSldfq24Hze80CtOzHrAp5C8XY2RZfZInRuhHyHlln2gxk6xRIRiLoR3jkreZu7SKEwzwMRW3pK2HdsMjXBQwzfGjm+zaUoFNMJ8KT4GgwVKvcD5/PhaYd3Xv9UeF2igXsrHO/7PXNXx9I19QNXta1XHkaZ8WlF72Ff2eTe5gsapoGrIhOGn6P9LMqBh/l5NjTKCS2qPDhCho07+7HbMO/kqZZCmAYP/vvpP3ncPHEBJuODhYYVV3djOlp9vp9G/N55xYyND/bbdUIMnAo1swvw4BUZ6UoxJYXyE79MT1x+mvKeAkFQ1zzQFBq46CPNv+iDx4+M3S/nwOMauIzDTLNkpZyGmN6NMLy/5NJgySuCnTCTZmknHoJQFrbkszKbGkbbi+YaeNLZljyds3JzuDWelcURKBTdsj4SWVJTLp3XEpvETdqXrRja/d2OOhc/2XmeDDwtRxHX/kM3Qppq4t62qwzXbzPWjqKgzYqL/atBLJ4WRZ8eOawYcOBsAw0QfffMU6dY8DY3iTt4XZdqXStF7LVbu9Ibi0Ezz8cQeb9aDy4n8LjJA80hwAuiAA+/3/ihYzGlqzPyu0yIsMEX1QbVmpNJR+dvqcWNsFxx0SLQRPd/8W0AQr40TdgH1j5McJcjmkfczTAJLYoAQyIqCi8U/prKd5cJQ/5V85ONDKbjTWVf4PONxmcXOHBhJcM08B19FSxYtVVK4z380vpgN6z4OxPgntGM80IB8TfyhPRHGgqFUQSOv/vPS0O8R/18Gtn+nXsXxa61CtuJVf2G7U5kXktOoIGH7ymkUIhU06641MiwCwDvOnpvTD9svPSIPh4ybxgVdnIHWchYgDDNBo2FUk/ENfDw+8j2Fpw4ZSwA4IQpHjcpi1QmO+FEpoEzv2+T6IL8S9HtxPzkyVNi13gOvCwYMWd+fRrGdLQCQBD8KVXUOD8pmeeNy2m1poZRcXJRQRWNkBdCOxUCnD3Ht4NXVnXeO/sruPBXT+Nvs1ZqVXDRhsIj2OkYMWJGQYN7vQZj8Ub+351z8Y5fPIV12+KxM2Yv34ILfvk0AIkG7lMoLQXmicHlRSmKBeZGKDcKq1D2BX6BhDsKxT6uW8rXKmPaWgQBrlFp2e5El4YeWAWJBl5wHGkYBRVdJ8N+uw+PGRJlj+r6iQjeEK2bFNnuziFtxORfLEHcWs4aPuT94lqP7LgsmZCeft2TfhoGGjjA+Smr36JMCPEae7niBmUHgHHD24JBzFwd0wzkhxat856VaJf81mXTNHUcMg+VF4oqjjYP/hDp8N5kLXv28i147OVu7RJVq4FLjNsqVztW9nZfUDHBsm2XfkOJygulteBEtncTQZBVUlIoz7y6ETc89ir6K24gFMW8V27ehTc2yqPt1XrwrhhIqyihUL57wWHBZ9enUJiCFtHAS5wGLhHgbgoNvOCQWB+Qtauun4jg37luIik4jtXAxYYVOU3GHxe4pafo1yx72TojiCkHzqCjUGT9QjzQQRSSLRwPaloesSyh62R0Kz3LybRTDW9rMbpPtcXejF+Oe6Fs3VWSep/Ey1fQapYmFAprX1FgXHX3/ECYsAm/TbAJJNWvIigK67d7q7zWohPwwAw8heL5gWuTjoBvK9aO4vL+qVc24JQfPSbdt1CrjBEVFZkG/qET9g3KF/iBS4yYLNxusaCgUCg1jsFfdEisbEnurklgk3B7i6N9/0XH8/myGrgPQgjEiZ0JP956L7rxyQSgrgPwPx2wx3D85ZMnxu5ReaGIkHUM3gBTqtDY8o1p5CUJd2+K/kCQR7VLGnw2S6etxcEvP3BM4n2/eer1SAxzBtU8ufeoYVy5/PMDuXd7y9OvG5UvKS60bmAyCoUJiY7WQsQ+cRt37igT9u0CVZD0bkQvFEa5tBQcUIR2g2Anpv/u02rgPEI7iPz5bZKokLXKGLEZdJyy4/i7UF1IKZQIB67UwM04QBbnm4esT/L5J9Ep23xqs7O1GGnjokPwpWkH4EB/J6zj01m1rm5UaAoBLg7AuAbOKBRGDcTT0G3ukYHXXPYZ3YHxu7XH7qEw80JhfeHoiaPw24++GWcctHvMzUxcbrIOFLg/VvH+Re2dEMZTxzc16UApxagOMy18rSSWsk47ftfRewftQ5BOC2IoV1zt8BCNazxKwmpleFsRvSVvC/s1/34pcm9gxIxp4Pryie3Ma+AuDdvnX/O97esRDbxaAe5XWaWkzBNOdgeQQczTaAK6sjuE4JXuHXjhjU2cG2GcQikoBHiF1kah6OLWA+qTgRjYBDistRCpZ2dbEV+Ytj8+etIkAN67ZJEX80BTCPD4Rh5BW/V/D7jDBL9cBpm2yMDL9q4RbVIemM9GF7eZLd9aHAenHbg7igUi7EiMUygFh4CQUFOvZiCL/C7jXNNq9a5rvktP1s6qsveVPe5fNoDToKSgbhhOPWh39bN++7K26vC3dM9buQU3Pv5q5N7AD9zQKwcA5q/cijO5U5yA8FT3loID0Ph7YP2dGSWrQUihyJ/nA6Ux1KolikXVrXAJIXhjY48fK4S5EYa/sz5adBypclRx9REaeRQMKRReyLcmvONev0xFJ7p7k+1aZnVi43hIc+DxaITCDX5jMU3L5EQPIOFMTK7jd41okxpkXM4RXLarMyxvtIMWC050G7mEQiGEoMVxuHjQVQhwwU+5tehg1ZZdQXxzYw0c6p2iIqTLXUU2fWXvSDI2WMuueT48PA1cnskP3n04DthjuPJZljfT+Ea2ewKc91Hn8wHSCfBv/nOB8jeRQmEocnFMqqVQgoOWFUrKa91xQ2atMoaN07ftPw6AfoXLv+VA2EnefcGRHzKcxohZlFEokkfTCPCK4pQutpuUN0iTIe8HzsdCQVxTe3W9dw7hQeM93kkmmDZJjGEmBzoAnleIbBm+esuuYCDo0hI5vqJDohRKxQ14Tx7FAgk0w2o08H5OMALeBDdj0Trc8fyKxDLzEL1BdHy4dHOUQjL0laLeN4Dc4JuEJGNWR4t6w3G/QFF1jWgDgJj2DYTtaBpaAABeFOJa82AUivhuw2iE1WvgbMyo7AOyo+1qFTIsz7cfuqe3etRSKPHPssmb58D5MZjGjdAhJDY51KqBRyNphp/bAgEeauCL12zDjEXr8FRCeI1q0BQCXKQXxBfNXsbUSWP872bp6oxf/EsZ0V6UUijf5zhSJrjest/Y2H1MM2ExSIqOE6dQJJKrpeBwG4u0VZGC97AoOCS2kjHdfEQ5T4GiQ3DuEeOV90p3JSqEUH/FRUsxvvIQkSTUdRQKAdGGAhDLO264J8AX+odiRO71X0JLwUkV+EmFVn+XobhiDDlwt6r3DoRtlmbir3W7d6Co+IY7vVITNmBAoUledIHbyMOPQdc1P9S4WIj3fVldeSEv2jlEBAdDQNTAoyEMZPXMEk0hwPlG8GKhRH//7gWH4ZcfOAb77+4tlXU7o3jIlvsA8LYfPhJ5KR2tcgHOo1RxUXAI9hkdj0/Oyht6y5BYJD5Z+i2FcBBUs5QO+F3XRdGJayGm53jymy2SIOcr1WVvM9DAdZt5AJ9CUWVB9NEURSHANHAZ+M1G4jbxT54yBUdM2E1bThEtBU8DF5uHDfRNO0tYvCY+kZggpFDyWrzHwYRkwfHeo6kAZ2WVvXuKcJzySo6ngZty4E5M6ZOt2nghrzN8A4IbMNf5mLLAxitfZtOonmnQFAI85gcufB83vA3nHjE+eEkqwSxC5fq3YtMuQYAXlHEdGEoVT8sVNUog7KzsZcaMmJW4F4p3vxM5iCKpU8XL5D/rc+ziADE9x9Nz8dML0dBrRm7E7Gwt4Pyj9or9JlIosq07SXOH1lgGPeUhChmmgetASNzHuaOliBMlqy8digUndiIPuw4A3713Ebq3q0+E0iHUwM3dH2qmUAIB7oAQkmDE5J6TbORhcF0a9FPWVxzib6U31MALJN73ZRMbr+AkUii8ay6XFosbFB4RaAW4ZCOP/j5TakB3PBk/qDrbCnAkrkgiig6RCCQDDlxFoRTDQeBSqt2aLAMLd1qukUJhQfcBtTcKW0HItK7+iovh7UWMbI+7IorxMuQauL58Jdc13sjzpWkHRJ8VyqvTwPn0xEmhpUhSbcVuLTheXSWGyjTpqMDeh06ZEftTrUbMQsD7ehOnTvvn3ynTXWTjq+xS/H32ar+8TvA/7U7MmAYu3YkZfk7mwMM4Q3xa7UwDp2y1Fj5jGtUzDaoW4ISQfQghjxJCFhNCFhJCvpBlwXjoYqFEy+T9N9XAH1q8TvlbWaBQgOSBVXCIVEsONHD/t2LBwdZdJUy64j48+9pGzwtF8lzEC8WlxkGlGNhJ7t5OTyc2QExPEeKjLqrABpcszXLFmwBkk1vciCnRwCXXDt0rjEJpqokBwJH7RGmONBQKQ4FIBLjjKCNKylAseBYRGYViGrpAm75fFt0kLZ60VKsG7nAaeJI3Ef+7o9HASxU3GKdsddtacFIZMYsS5UVmq+HvSfI04vdYSDVwzg2SobM1++jdtWjgZQBfppQeDOAEAJ8hhBySTbGiiAoeotTImKBXddqkWZWHK1AoQDIvptLAiUihcPV5dMl6lF1Xql1HvFCqoFB6fE+DsoIDNzkgGRA0cIUoZ2X71WNx741SxfUEuIReEt+J7N3Kcoz44Gs28iRRP7wGXnQIRg1L3rDkOEQauClpgv/UKfvhvVMnBHkRf4eeyg+8FrB+KFNmWJN0CseY1aqCR42YSfdKBLjkIT70LqMxW4oOKq65F1WhEC1PwZHTO8UUApxNHmXXjXDg7RojZkNp4JTSNZTSF/3P2wEsBrB3VgXjkRQLRbyuFOApBCAfVyLQwBM0o4LjmFEoXDo3PfEaXCpfPrYUHNy/cC0+c/uLcF0ZX6zHzr4whkfRiU98pkdSuVQflxvQT46ligvHkbd/a4KHkeoaD+1SXfwupMULgc62YvCudXAIicVIl61wRFxy/MTAu6HFp1DEaIRAsgeECZhCIBsL7PxJkZOtWQPnKJSkd0YEgeo9H7+PF+Bs3LANaaYaOB+dEfDiucs28fH3JCl7gYNAhUr9wEODd/hMmv0DpsgkRULIJABHA3hO8tulhJCZhJCZ3d3dVaUfjYVSPYWSlkNm6GwLB50OxQQjJtMg+GUVm7xlaTNa5b55a+BSmmoFAYSdv+JSFArJHL4KFTd5QOomt/6KF3HOhEKRauCSa/zQ1XmhiM+KSfEUyvC2ovQAZxEFJx4PxUQDbys6oSdSgQAgPoUSLbyYdjVg7SoT4MwWIdJFNXPgESOm+XO6jTz88XPFQAMnqLg0oBeTIFIobS0FqQbOj4+ksc6er7g0OLcWCLXs8IzXMJ2k1WA1qLmnEEKGA/gbgC9SSmM+T5TSmyilUymlU7u6uqrKI+JyBPmL9ssCh6ijDKbVYBnai2YCXOZe5pXL/51zI5Q9K4KPw11x0xsxd/aXQSnFy+u3o8WAl1SBGmjgurYplV1jDlzuBx6/xoyWnqslTdwG/vkz9sdXpx8UqwevgbcVHbS3JvcRR8aBG2jgbcVCsAop+kKOSoyYBSc5rSQUAw08vsoa3ekJ8AmjO3Ak5/pY61b6UIAnH47AK1nBRh4ZhcIJR9b/mRFT5mFzqST2Ph9eF/A1cImM4Ns8qf2ZEbNUcSOxwVm/cCUaeB6oKXlCSAs84X0bpfSubIoUh6jdEU2pHULQr3CPq1aAOxLqQ54+kQpZptmEboQyl0E5hcLg0uoolNueW44Fq7bhtQ07qxfg4ISoIgltwKiKC8eRt424qpB6oWjK1l4soOwma+CXnXkAPn3qfrHfeQ21WJAboUXIvFBMKJS2llADb/F52Yob30pPiLw907ihseeZoLzw6JDdPH7y2KAMd//PScH1ajXw9x+3D+7/4ts4CiVZWZCtDGTNxw4AueUjU4PphRkxSxUaoyX2GR1GuPzhu4/Ax06ajOMmj4lSKC0FKf3C38NPHDLw4R9kArxC4xp4HqjFC4UA+A2AxZTSa7MrUhxJ0Qgjvzny2AlAMoWS/LuBBi6hOdhsXwwoFL2wZojsPKO0CgFexnwu6ly1Wp3JRh6tBu56z8vuEQWhLBudLGhv9ZbD/HAcN7yVS0+Y/IXvvNeMQ+JBj2QoEBITHK0GFEprwQn6QMEh2H1EO9Zu641RKATx020ASCNiqlAUKJRpB++BX37gGPzwPUeg058IRNfSagX4mYfsgYP2HKk1Yh4yfiT+fOkJwXdeA2fZ6vpne0shaKfWohNEahTfA//+xnS24hvvPAQdrcXoLsuWgtRuwisTG3ZE/e+fvuL0yHc+yNxmzl7G6K+AA8+BNuFRy/RwEoAPATidEDLH/zsno3JFYOoHzn5Tc+D66l7zriMi3889fDx+/7Hjgu9JA7ToOFItKTgsIMJ/RiHbbs6Xt7/spubAd5Uq+PPMFcF3frAePXGUcTre4cLVc+A6CkUUSlIjpqTdWXMNaynElsNpVhqmhlweREKhFA3cCB2HRN7hpHGd2N5bjgkLQuQGr9EdrfiPIEhaFLQcu87GgkOAc48Yj/dO3Uc5Diiq2/ErHqjifY2W5/gpY3D8lHCjE09zsnepe2/DWgoRe5Hrh5MV30PEu4WrJp90e4uD1zfEg3kN57xyxINEuoa3ReQOz7/zcZZkRsy9Rw0zsq1Ug1q8UJ6ilBJK6RGU0qP8v39lWTgGvnPqjJiA9xuvgU87eI/gc6IRUhBC+4zpwMkHhLy9iQYuO70mNGioOXDZpMOvCPrKlaqNsAx8vt9856HGz3khP/X3aDVw5kYoKf+eolaZkkIZ5mtT/JJYNHrzUB2XBpgLfo9CEYyYRcfI/Y+/Z/I4L+zCIy+tj90nm6zXbO3FmM7W6H1+u4sCPzRiehMUPwHLKDzA4+OrCZ7F6hR1DZTfw/CKH4AOCLl3rQBvLQRUU0uB+Dsx3dhKpRAR2nJOW+UNwitfG3dEBbhDEPFQ4v3IH3s5dM44bvIYAMAJ/mR13OSxeOIrp2H+t96urFstyJlizwYxTwIivw54SxYmDH/w7sMxcUwYmyRpI4y4nV3U7JIEaLFAMKI97obGH9fl3Rcvh8wnm7+vVElPoYhg5f/saW/CWEEQ6ODSZA2cHxRi2kvX7/DcCCXtL+7ONHUjZIOeUSg93Cnh/GAd1REti871TMxGveM37urXIvENl4HfTXjIeM+A2CMc9Fx2aST9yeM6AQAn7jc2thplbdomMaoCYbwb/jHWD2Syupqol6Hx0vvvSryWtKsTGk1HBk8DZ4ZrL4bMkrXbY+9B5l8ufla5aXa2FfH9dx2O4yePiY1jQghOPTBU5njvpYpLMe3gPbDsmnMxwY+FdPIBXVj8nek4bvIYFByinDRrRXMIcE4HIwh9OlV8KRPgE8d0RjQuMXCSCLEDibu1kl6Cp4HHBXhZ4MBlHVXmLSCeBi8zbF125gHYY2Ty7kEAaC2G3GeqTU0cB64aYnzZrpeEm12wapvRBCQrl27uGNbioL/iRiZA1r6XHD8RJ/uxqRm0MaqFfFQCReqFUnRw5iF74EfvOSIQuDLwnhR77taOcw7fM3ZPqRzldoe1FPDE5afh6gsOiwlGJoxiBwoLfuCy48JEUU1pdcGviqIApxRrt0VPZhJf/XXvOyo45YnlqFvADOM4cFbXskuxXshHtQqQ+Z2L6Gwt4v3HTcSfP3ki/nTpCfj2eeEq1SHAjy86Ehf48XzEfjRacmJVHht3RDSHAI9pRt6FYyaOjt1bcAiWrNsOwBssvACX+WjzEDVs0c80UQN35Bp4WeDAZenIDK+iwJOtIC7yd/bx+Pzpb8LPLj4qdp1lq9oxqoJLk+kFflcfb0TkIeY5QjLZiRtkPOg5cBHMcHTu4eNjKwetBi7ko6qzd8hAdMItOgQdrUVcNHWfyFmfIoqBBu6lvduweFuVKlF7h0spJo7tQHtLISbk2H3ihNIqUCiRsKZKDpzGDmA2AROIrK1lNIyogV9w9N748tsP9PKV7FoU0R6hUMK0tgneInw2qlCuqomC7fcAgH3HdkaCrzG7x/7+WZeiH/nuhkpU1sh+c34OkGlG//jMSZjcFdd0IppGIXpidJLQSoqXYMaBqymUkAOPpyPVwGOeDvHnCoTEBM/uI9uxm2RLeLDUVbg7qmDiBz6cm7hE2oIhFsGvLS58ZTSEbiOPTMtxBIHCwzTAEp9HLH1CsEo4sYfvGzofaJ5CAeR8bMl1I9f5PizWaXRHC5ZvkmjgTtSIyT/GPsoolK2Sw46TEKws/YRlhlCZfUBcCSQZMV2OQlFBDD0tv67QwIWxq3uPz7y2MfK9yyCKZR5oDg1cooEduc8oeXQ77uUWHRI52zJR6xT6nTjLygQvj4JDIoIsSEfkwCUd45Mnx32URW1eehCrZIt8e0tBGkKVj8mStQbOa9OyVQgQp0dk29ZlGrjYXD9935Hh/ZJ6hka1eBlkGji7X7ydCaJPn7of3nNsuNJxCMGV5xyMi9+8T3AtIsA1TcUmsZC7lgjwcpQDV805C799ViClRAHOzmLsC7xQZIItbtBdvdWbmGQ+8zJMHteJA/ccEeQJyNtYRlvEgtQpuiTrr7pdy0EaSg5cnS+DGGxKHlhNnm/XCHMXzyzRHAJcwWXJwHdk76WHnSnJS0DciSbu6EzSWvvLrtRAEnLgciPm505/E/YZEz8IYoygye7sj2vpss44TCLAf3zRkcG9RUXMFhVMYqHwKw9+pXDHJ0LfX3HCldEfUg1cEK0nTBkbCFdZGrrgSDIOPEhDqCSbeKcfuifGcrRQwSGYPK4Tl50ZhqblPZh0vr9MCQi9RyT9RaOB8+hsKwZqtJgOO1O1XyLAGcRkSxUXF9/0LABIDyaR4TcfnhpMRm/Zz7M3TBwbf1aWPxPqrBwqBxj2fthY5hWBz5+xvzIf1WfmQbTnyKjQFT2L5JvKwovjhrcG53/mEevbBM0hwFPcy7/cou9uFKSTkJC4M1fUwJOEnuhNEKTDKJRgG3W0IKp0RwsCnPe0YHBIfH0yrNWJdcb3HDuBO1Q5XVwUXgNXtSG//OQF5/7cgcKjBe+UThmFIqEUxKI6hAQ8q5QDT6mBMxqG3T5136htpa3FiQxc9roiZyhy71C2RBfpLy2FInDgMsEmHmsnm/haCoTTwMPrqoiSa7aGBsEJo9U8Pg9eMH70pEl4+orTcdCeI2P3ybbpFwQKJZ629z+Mse19Z22296hh+PhJk4U048+L5ezw+2o8EmayDYS/1FdyA3dB0/bKGs0hwCX8nQrx8zL5Z/VPi5qOeKJ30lZ6tinkq9MPilxnQqMl0MCj6ag8QkYJlm0WXZCH48Q7XnuxINXswq3O0fsfuuxk/NeJ+0rLAKTnwHnwJ9yIoVqHSSgUWbllA4u1aSwkKvgJxIwDZ9oTe+y2TxyPeZzfblsxajwkknbk36k4cRw8fiRmf+NMAOEKgNk3ZAK8P0ahRMv80nen4+cXHw0gFH6y1WVv2cX8VVsjZfY+x24FEDWkq96nCCK0i86Aq3qWadbiRLW3LxTZJM3agW3+2tZbitEufH9Q+YGznahJnlhSAc597i1X8OlT9sMjXz4lMG7WG00hwE1CfMrQ4jiCAUh/Px9rAYgH+086Vq3X18BF/vCwvT1/3/38MztFLl2lgYsGUZkGLtOk21sL0oh2qo1Ek8cNlxpfGfh44CqI/OGdnzoRD112cuSauMzslCw7ZeUWm90hIT8t0+ILgmDgIdPAGd3EBnxbsRChe9qK0eh6wRFgXDvy71B8JyPai4FgYbRcQKFIVhAlYYOKWI32loLkcA79WaSyfqKiLL75zkOMY5LXEsGQtTdLQtTS99otKsBZXlO6vHG0vbccGUv/d+FheOubQrdRlRshUxzYOzt4/Eg8c+XpsfIlNUGp4h32zcozEGgKAb7XqGH43oWHAzAPyXj+UXth4tgOrQVfBBvwjO/ca1SUI0tyQ1Rty77k+Il46LJT8OZJ3i4tcfu4GBObQfSwkGrgJGrEPPfw8ThqwiipMVB19qDK+MrguqHmoVrFiB4lb540Bm/aPaqViO0vTsx3/c9bpBq4SG0VnJBCkU3ubFDLtG2dBq56u21FJyIoAm8erj78BLaXoIXyNA/TwFuLcuMju6eDe8Zkc41scuch9UJRkBcHjx9pvCu1BvkdtrciEabNMwqFjU+2gxWITkyXHL8vCCE4YYo3zvixz+rT0VoIFAemBIxoL2L8bvGVA0v7gog7oVHV6oamcCMEPM7va3fPx7SDdze6n3FTfOdPmlEnj/Nm0v85dT90jWjDKQdE80ryQlGdcEMIwZt2D2dpkS9TaeATBEPSeUfthRuEE29ETemKsw+C48Q3mgChxijTmkZIPHoYTDTwamI9MO356StOx86+Mg7YYwQe2hE/5k4WsIjZImUrB/aaZO5ssk1PbBJQVbGtpYCDx4eTUXh0WNwwBgCXn3Ugbn9ueTCh87+xd81cLeUCnEYMgbrzPhlY7Pep+47G5wTDHhA1IIfUhTytjtaC8bFuKgOrCYJy+BJ8lzB+GIXCNuCx19nRWsSItiIuOWFf6Urhpv+aitufW44jJ4wKrrF3NbK9hYvZ7btYKstHMOcbZ0b62K5+80Oi64GmEeC7DWvB8187IxYLQoXW4KVzGrji3lMO6MK6bb04cM8ReP6qM9A1vE1uiErgzHgN/G+ffgvunbc6sFLzEL1QwkHdgi09oR9u14g23PKRqfjYrTMBAJe//cC4AC/Ig+fLBAPTGGUbLWSbahhMDjXWnTZy3+ffGtOsT5gyJqCaeN5U5Qe+7JpzMemK+7x6OARXnXsQrrxrfiyWykF7jgjd2ST1PPfw8Rj+0SJWbt6Fr/99AQDOiKmoXFvRwVR/9QSAO7yX55V5YV7AeUfuFQQS4ye3dx29N7q39+FjvvFNtuIoVVxMGhvucdAp4KyK7FCBT5w8BaccEI+7r4vQKKKjtWg0aaSBLDlWDvYb21V53fuOwrbeUuBEwGgn/piy+d8+S5nXyPYWfOqUKI3JXs9uw1qC92ay61Tc07C5p19x58CgaQQ44G1QMUWbENYRkG+FX/yd6RGqYneNP6doVNTh2H1H41jBm4HHBUfthb/P8U7bZhPDk185LbYjk5VnRFtRubFgbGcbVmzy/HeDyG7cvSN9eoRdky3Jp05Sl9Wl8RjsDokKFt0xYIfutVvwedbXp3llVmx8ED1vWP5i3hcePQEXHj0Bi1ZvC66/cNU0jGgv4hO/9yY8WT29mBbeyuo79y5Cf9kNBKxKrLUUnIgxljWtaahQfjVULDj4zGlvCr7LJqyy6+IgTuM30XIvmjoBv3z0VRy1zyjp73y7Mj/9Mb5Qby04EZfZjtaC0YHX7z9uIqZowgYkQVwJ/OS9R2LR6m24wI9dfvfslQDCkL9u0LeryMx/drdhLYHWziaINLTIlgYT4E3BgacB04SYBv5JbibeT2JsSBOvQCXc01jeGa67+OhAO2fDc0R7S0RQAKF3g25X2E0fOjb4LOM1WSxjNn/JBMKE0R24+b+mStOnEgpF1FZNAjkBnuBWCW/AM/j+5KIjI9fE8sp8egFvxdLeUggGaBJ3zARwIMANBzJ7F7p3wr8HGZ3FINtdWypT7D6iHXd+6kQACRq4n89Zh+6JZdeciz0USg5f1umH7YnvXXg4vjTN82P/j2DA62wtGrmZfmna/spVy68uOSbxeZGLP3Sv3XDR1HBzFNN+mS87WxUkrSBk2O5TTCOHFTkbCaNQzNPbUsVO1Twx6AQ4A7PuHzNxNK4423Pr02nEJhDPEAQ8b4u/f+Ykyd3JOM5flm/eqZ7VGb2iG1C7j2zHvj5nKhvsjN8OKBSFRFBpNrLbxeKkOWwgCe8+NhrfRSyveLqKCN2OwEi6bBOMPwmIA/mu/3kLrjrn4Hj6BpKen3N0Alw2mTJtmAXFqpbOePIrp0mvE0LwgeMnBuUaN7wtwvMOa417uYi48Oi9peOB4ZzDxyt/C8vh/VdVj7mdslXp+Ud5mrnOY0qFbb7gHclRKNVEXmQU55ETdsMNBpNU3mgqCiUNeM3m0rdNwYlTxuJIxfLSFLJ4B8yz5HsXHo4juPMFTXDpKVNAiDwgFQNzXUw85dv/rxvsTFtS9VuVNnX5WQfG8vHupbjns2/F4SnrnRai0TdqPIwLx/e9eSIeWrweh+6tLxfT7A7zKR5xAjtm4mhpwLSk8x5F6A4pFnnY/bo68TV/0mATRbWnmbPdvTJ3TRH86e+txeTj4b59/qGpDumV98uoG6EIRvswCuWqcw/Gl87cv6oof2xz0VmH7hn0Yc943Z+KQpkyrhNzVmzBDR88NuZtNBAYvAKc6/SOQwLhPbqjBZt7qlsG6bTMDxw/MXV6bcUCPnt63GOAh+4gZB6iTy0ATDt49yBOBRBOAqpTV8QcPvKWSfjwWyZh8rjOyAD30vL+TxpntuW6Wtz4wWOCMxzFvAG598uZh3ixmZNw7L6jMeuNzYEh1HQpnfZoOp0GLvrPP/zlU4PPoztbcflZB+Lsw+IhZxm+Ov0gfPFPcyJeTjzmf+vtVZ2FytxQJ4wehpWbd2FsZyv6Ky6295bx7y+8TRqHKC0SNfCOqAZecIjWW0qHwyfshkXfOQsdrUW4LsVXph+IyWM78enbXkyVztUXHoZLTti3IYQ3MJgpFIXWwrxYbrjkmMCgZorRna0Rnpj5m+aJ8HRrUw08vHbzh9+My88Kd4WyRYnq1JUjfLerGz94LJ772hn41nmHBsv4dr89vzBtfz8/5pUicuP6+qTF9MPGB1vwWdwLPs9qtVMAuP0Tx2PRd84KJj1jDtzgPr6FdXEy3jxptDT0L8NnTnuTdqPI2/bvwqz/PVO52W1Ee4t0t6oK7/VXg6M7W3HDJcfg1o96RwqWXRpUylR43XDJMThkfHxbvSlGtrdgv65OXPPuw6tOgwdrI8ch+J9T3xRMBmn6bEdrsWYqNksMWg1chbGdbXi1eyd262jRGtNUeCvnFnjLR96cZdGkYIP/+KTJImlXBDyvgT+/sALv8ONo3PjBYyKbg8Z0tio112LBifzmBNpTmN/fP3MSdtfworXisjMPiASQAtLTGTyY50xgHDMcySYabRg5j+CiY/dR3kcIwflH7Y0v/GmOUd5544fvCQ3IZx8+Htt7vdVqxaX46vQD8b//WGgcuOnsw8djweqtWLRmm/T3sBVVNhkSWZFkDbZKHyWJyd4sqEmAE0KmA/gZgAKAmyml12RSqhzBNPDNO6ujURi3/uUzD6h6i38ajB3ehge+eHJAVXzylCn49eOv4cqzD4qc1ynTwEVM6RqOed8K/WenH5ZsaFLBkfDpKhe2avDo/zsVb2yMHzybB1ibmU4F/GpoSlcn3nnEXrF7mGfF9y48vC4ns9SKX3/oWMxZsSV2nffY+NCJk/ChEydllmdA+2Xrcm6MN08ajf99xyGRUMHNhqolECGkAOCXAM4EsBLAC4SQf1JKF2VVuGqw96hhWLxmm5J3ZNvkqwlcD3hagQm/miV4HvvKsw/GlWfHPSP26xqOV7t3avnWLPGFafvj6vsW53ba9uRxndqjyXj8/P1HY38FB2wC5qopi6vC45iJo/Di8i0RL5RHFBoiM3ibcrYTx3Rg+aYeo3vzwFmH7omzDo1z7SyE8tH7VEcbnHfk3vjlo69KlYXQjbA2fP3cg6UuwkkghODjb51cY+4DDEppVX8ATgTwAPf9SgBX6p459thjad7Y0tNP/z57pfL3zTv76OV3zqE7+0q5l6We2N5boo8sXjfQxWhKVCou/emMJXTTjj7tfdt29dN5K7YYpbmrv0zveO4N6rqu0f07eku0e3uv0b31xvyVW+i2Xf2Zp/voS+vovl+9l37y9zMzT3uwAcBMKpGphFa5fiGEvAfAdErpf/vfPwTgeErpZ4X7LgVwKQBMnDjx2DfeeKPKqcbCwmIwoeJSXDtjCT520uSq7FFDCYSQWZTS2E67WrxQZJRhbDaglN5EKZ1KKZ3a1RWP0WBhYTE0UXAILj/rICu8a0AtAnwlAN68PgHA6tqKY2FhYWFhiloE+AsA9ieETCaEtAK4GMA/symWhYWFhUUSqvZCoZSWCSGfBfAAPDfCWyilCzMrmYWFhYWFFjU5MlNK/wXgXxmVxcLCwsIiBQbtVnoLCwuLwQ4rwC0sLCyaFFaAW1hYWDQprAC3sLCwaFJUvROzqswI6QZQ7VbMcQA2ZFicgYStS2NisNRlsNQDsHVh2JdSGtsJWVcBXgsIITNlW0mbEbYujYnBUpfBUg/A1iUJlkKxsLCwaFJYAW5hYWHRpGgmAX7TQBcgQ9i6NCYGS10GSz0AWxctmoYDt7CwsLCIopk0cAsLCwsLDlaAW1hYWDQpmkKAE0KmE0KWEEJeIYRcMdDlSQIh5BZCyHpCyALu2hhCyAxCyFL//2jutyv9ui0hhJwlT7X+IITsQwh5lBCymBCykBDyBf96M9alnRDyPCFkrl+Xb/vXm64ugHcmLSFkNiHkXv97s9ZjGSFkPiFkDiFkpn+tWesyihDyV0LIS/6YOTH3usjOWWukP3ihal8FMAVAK4C5AA4Z6HIllPlkAMcAWMBd+yGAK/zPVwD4gf/5EL9ObQAm+3UtDHQd/LKNB3CM/3kEgJf98jZjXQiA4f7nFgDPATihGevil+8yALcDuLdZ+5dfvmUAxgnXmrUuvwPw3/7nVgCj8q5LM2jgxwF4hVL6GqW0H8CfAJw/wGXSglL6BIBNwuXz4b1g+P8v4K7/iVLaRyl9HcAr8Oo84KCUrqGUvuh/3g5gMYC90Zx1oZTSHf7XFv+PognrQgiZAOBcADdzl5uuHho0XV0IISPhKW6/AQBKaT+ldAtyrkszCPC9Aazgvq/0rzUb9qCUrgE8wQhgd/96U9SPEDIJwNHwNNemrItPO8wBsB7ADEpps9blOgBfAeBy15qxHoA3iT5ICJnlH4AONGddpgDoBvBbn9q6mRDSiZzr0gwC3Ojw5CZGw9ePEDIcwN8AfJFSuk13q+Raw9SFUlqhlB4F7/zW4wghh2lub8i6EELeAWA9pXSW6SOSawNeDw4nUUqPAXA2gM8QQk7W3NvIdSnCo01voJQeDWAnPMpEhUzq0gwCfLAcnryOEDIeAPz/6/3rDV0/QkgLPOF9G6X0Lv9yU9aFwV/aPgZgOpqvLicBOI8QsgwenXg6IeSPaL56AAAopav9/+sB3A2PRmjGuqwEsNJf1QHAX+EJ9Fzr0gwCfLAcnvxPAB/2P38YwD+46xcTQtoIIZMB7A/g+QEoXwyEEAKP01tMKb2W+6kZ69JFCBnlfx4GYBqAl9BkdaGUXkkpnUApnQRvLDxCKf0gmqweAEAI6SSEjGCfAbwdwAI0YV0opWsBrCCEHOhfOgPAIuRdl4G23Bpad8+B5wHxKoCrBro8BuW9A8AaACV4M+3HAYwF8DCApf7/Mdz9V/l1WwLg7IEuP1eut8Jb1s0DMMf/O6dJ63IEgNl+XRYA+IZ/venqwpXvVIReKE1XD3i88Vz/byEb281YF79sRwGY6fexvwMYnXdd7FZ6CwsLiyZFM1AoFhYWFhYSWAFuYWFh0aSwAtzCwsKiSWEFuIWFhUWTwgpwCwsLiyaFFeAWFhYWTQorwC0sLCyaFP8foUuCgQHXMdYAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plt.plot(res.resid_distance)"
]
},
{
"cell_type": "code",
"execution_count": 20,
"id": "aeeb9710",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"[<matplotlib.lines.Line2D at 0x1d9d8372100>]"
]
},
"execution_count": 20,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD4CAYAAADlwTGnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABhUUlEQVR4nO19d5wcxZn2807YXeUsEJJAQoggshBCZBMtwBgcD3wGjAPGhrNxOJ8w39kYHDjOgA3GcIDBYBtjbGyQkUwSOUuAEBKKCKGcs1ar3Z2p74/u6q6uruqunumZ2dHU8/tJO91dXamr6q03FjHGYGFhYWHReMjUugIWFhYWFrWBJQAWFhYWDQpLACwsLCwaFJYAWFhYWDQoLAGwsLCwaFDkal2BJBg4cCAbMWJErathYWFhUVd466231jPGBsn364oAjBgxAjNmzKh1NSwsLCzqCkT0keq+FQFZWFhYNCgsAbCwsLBoUFgCYGFhYdGgsATAwsLCokFhCYCFhYVFg8ISAAsLC4sGhREBIKKJRDSfiBYR0STFcyKiW93ns4horHu/hYjeJKJ3iWgOEf1EeKc/ET1NRAvdv/3Sa5aFhYWFRRxiCQARZQHcDuAsAGMAXEhEY6RkZwEY7f67DMAd7v1dAE5ljB0O4AgAE4logvtsEoBpjLHRAKa51xYWFnWIxeu249UP1te6GhYJYcIBjAewiDG2mDHWDuAhAOdJac4D8ABz8DqAvkQ0xL3e7qbJu/+Y8M797u/7AZxfRjssLCxqiFNvegFfuPuNWlfDIiFMCMBQAMuE6+XuPaM0RJQlopkA1gJ4mjHGR8kejLFVAOD+HZy49hYWFhYWJcOEAJDinnyMmDYNY6zAGDsCwDAA44nokCQVJKLLiGgGEc1Yt25dklctLCwsLCJgQgCWAxguXA8DsDJpGsbYZgDPA5jo3lpDREMAwP27VlU4Y+wuxtg4xti4QYNCsYwsLCwsLEqECQGYDmA0EY0koiYAFwCYLKWZDOBi1xpoAoAtjLFVRDSIiPoCABF1A3A6gHnCO5e4vy8B8Fh5TbGwsLCwSILYaKCMsU4iuhLAkwCyAO5ljM0hosvd53cCmArgbACLALQCuNR9fQiA+11LogyAhxljj7vPbgDwMBF9BcBSAJ9Lr1kWFhYWFnEwCgfNGJsKZ5EX790p/GYArlC8NwvAkZo8NwA4LUllLSwsLCzSg/UEtrCwsGhQWAJgYWFh0aCwBMDCwsKiQWEJgIWFhUWDwhIACwsLiwaFJQAWFhYWDQpLACwsLCwaFJYAWFhYWDQoLAGwsLCwaFBYAmBhYWHRoLAEwMLCwqJBYQmAhYWFRYPCEgALCwuLBoUlABYWFhYNCksALCwsLBoUlgBYWFhYNCgsAbCwsLBoUFgCYGFhYdGgsATAwsLCokFhCYCFhYVFg8ISAAsLC4sGhSUAFhYWFg0KSwAsLCwsGhSWAFhYWFg0KCwBsLCwsGhQWAJgYWFh0aCwBMDCwiI1MMZqXQWLBLAEwCIxZi3fjBGTpmD6ko21roqFhUUZsATAIjFeWrgeAPDsvLU1rolFV4NlAOoLRgSAiCYS0XwiWkREkxTPiYhudZ/PIqKx7v3hRPQcEc0lojlE9G3hnWuJaAURzXT/nZ1esyyqATvZLWTYIVFfyMUlIKIsgNsBnAFgOYDpRDSZMfa+kOwsAKPdf8cAuMP92wnge4yxt4moF4C3iOhp4d1bGGO/TK85FtUAUa1rYNFV4egA7ACpF5hwAOMBLGKMLWaMtQN4CMB5UprzADzAHLwOoC8RDWGMrWKMvQ0AjLFtAOYCGJpi/S0sLLoQLAdQXzAhAEMBLBOulyO8iMemIaIRAI4E8IZw+0pXZHQvEfVTFU5ElxHRDCKasW7dOoPqWlhY1ApWLFhfMCEAKn5O/syRaYioJ4BHAFzFGNvq3r4DwCgARwBYBeAmVeGMsbsYY+MYY+MGDRpkUF0LC4tagVkeoK5gQgCWAxguXA8DsNI0DRHl4Sz+f2KM/Z0nYIytYYwVGGNFAHfDETVZWFjUMSwHUF8wIQDTAYwmopFE1ATgAgCTpTSTAVzsWgNNALCFMbaKiAjA7wDMZYzdLL5AREOEy08BmF1yKyxqArvbs7Cob8RaATHGOonoSgBPAsgCuJcxNoeILnef3wlgKoCzASwC0ArgUvf14wFcBOA9Iprp3vshY2wqgBuJ6Ag4oqIlAL6eUpssKgyyVh4WFrsFYgkAALgL9lTp3p3CbwbgCsV7L0NjE8YYuyhRTS0sLLo8rAiovmA9gS0sLFKDFQvWFywBsLCwSA2WA6gvWAJQI4yYNAV3vvBBratRHuxkt5Bgh0R9wRKAGoCHzL3hX/NqXBOLWmDe6q3YsrOj1tWoCGw46PqCJQBlYFtbB7a1JZ/Iu80cscZAJWHir17CBXe9XutqVARpDO1xP30GF9z1Wgo5WcTByArIQo1Dr30KALDkhnMSvVfcbSiARamYu2prfKI6RBpDe/32XVi/fVf5GVnEwnIANUDRrv8WNcCl972Jp+asrmwhdmzXFSwBqAGsqZxFLfDc/HW47A9v1boaFl0IlgDUALuNBGh3aYdFarCbm/qCJQA1QL0TAHsgjIUO9T62Gw2WANQAXAlcq4X0dy9/iEVrt9em8AZHrcwkq1WuXf/rC5YA1AC1nCSMMVz/+Ps4//ZXaliLdMEYw2MzV6Cto1DrqsSiVjvkapVr/QDqC5YA1AAeB1CTsp2/O9o7a1B6ZfDaBxvw7Ydm4udT59a6KrGo1fJYrXLt8l9fsASgBmBF5y/VQAaUpg9CV5nsW11nvDVb22pck3jUygekaiIgRTEL1mzDjl27z4Zjd4IlADVALS0lasl9WNRQBFS1coIlFYoMZ97yIr5uzU+7JCwBqAFq6QjGF6ByuI+uRjzqSexcK+JftT6SyuEbjtcWb6hSBSySwBKAGqCWirLdmQOoh5PKascBWCsgizAsAagBOAdQSyXw7mTLX0+Lzu5uBdRVyrUwgyUANUDX0AFUjgK8umg97n91ScXy16EeiNru7ikrL/i8vdY8tGvCRgOtAVgNd+HcAimN9V83qb9wzxsAgEuOG1F+IUb1qEoxqaBW+p+q+QFIBK6evk0jwnIANUAtw0GnoQPoqjvtrlovEbv7Tng3b95uB0sAagCPA6iBFqAWYSgYY/jVMwvwwbrKhJ+oJ7FK7RzBaqMETmuzs7m1HRf97g2s3db1fT3qCZYA1AC15QCcv9UkPptaO/CrZxbi3+9+o6LlWCug2pcrczhplfuX6cvw0sL1uOelD9PJ0AKAJQA1QS3ZZFYjDgAA2gvFmJSl5l+RbI2xaUc7Vm7eaZS2ZsHgqlVOSAmcdv71w+3VA6wSuAbwxnANNqy1UEJWukgv/xoxAMfeMA1tHUWjo0FrxwHUNgRFuaXXg36nHmE5gBqga5iBlo+ka0ql53Ct1oi2DnPOZrcPBlchDqAeNv4rNu/E6yl7PHcWirjnpcVo76wM92wJQA1QW0cwLgIqJxRE19qO1ZNYoHbB4GpSbF0s3GnhlF8+jwvuej3VPB98cyl+OmUu7n5pcar5clgCUAPUcsFiNSQ+lUYtoqsmRc0+fcJyf/iP93DrtIUlFKNmAcptdx182ors0re1OVFUKxVN1YgAENFEIppPRIuIaJLiORHRre7zWUQ01r0/nIieI6K5RDSHiL4tvNOfiJ4mooXu337pNatro5bhGIq7MwWoA9QsGFzCch98YylufnpB8nI0nsAWpaHSRhuxBICIsgBuB3AWgDEALiSiMVKyswCMdv9dBuAO934ngO8xxg4CMAHAFcK7kwBMY4yNBjDNvW4IeB+1Jn4AcMuuHqq1660Hmra7i2LkYhpJBFQJVLr/TDiA8QAWMcYWM8baATwE4DwpzXkAHmAOXgfQl4iGMMZWMcbeBgDG2DYAcwEMFd653/19P4Dzy2tK/aCWcyINHQCHndvJsdufByD7AaSev7qc3R2V2iyaEIChAJYJ18vhL+LGaYhoBIAjAXBvoD0YY6sAwP07WFU4EV1GRDOIaMa6desMqtv1UUtHsDRYyqTvVloMUE9rQe3OA6iNJ3Ba5YoL4FUPvYMjr386lXy7Oir91UwIgGq6y/WKTENEPQE8AuAqxthW8+oBjLG7GGPjGGPjBg0alOTVLotaLli1FAFVWudRD4rC3Z8DqHy5j85cic2tHVi9ZfcPC8H7M1MrHQCc3fxw4XoYgJWmaYgoD2fx/xNj7O9CmjVENMRNMwTA2mRVr1/UIh5PuOzqFV7pRa+eFI215P6qg8pHAx3YswkAMHPZpvQz72LwxnaF5qsJAZgOYDQRjSSiJgAXAJgspZkM4GLXGmgCgC2MsVXkrDK/AzCXMXaz4p1L3N+XAHis5FbUGWrKAfAD6VPIa1Nru1mZVWpwHTAAu78SuMJWQAxAcy4LAOgo7O7E1O/PJ2evxrptu1LPP5YAMMY6AVwJ4Ek4StyHGWNziOhyIrrcTTYVwGIAiwDcDeCb7v3jAVwE4FQimun+O9t9dgOAM4hoIYAz3OuGQC0tMdNcjP/+9go8MXt1bLqKh4LY/deBslEzLqkCxWbcVatQy8O1q4z5a7bh/VWJpOdGMIoFxBibCmeRF+/dKfxmAK5QvPcyNOscY2wDgNOSVLYroVjG4KuFGIYjbXn8mx9uxMRD9owps0ocQB0oAWomAqqRGWgl1uiM+507G4AAiC2sxOi2weBKRDkTuZbDtlCDSHR+V1WmzHrhABat3V6TozKBWiqBU7ICEoYO/1koViY+TpeC0KGZCmxwbCiIElHO5qPau8BdnQWMmDQFf3treUgBvWNXJzorFKaZo/JKYAddff//qd++gj+8/lFNyt6djoTkC2Haw3bphtZU5OxpcrxiTpWwBLIEoESUxQGU+Oriddsx7qfPYNUWs9jzHBt3OMra/31ynuCF7ODgHz+JKx58W/vuzvYC7nvlw1REXhVHjSlAXB9tr1A8FxMk2Yk/8/6a0supghko37ykzQGc9L/P4eifPVN2PhUb7pYAdB2URwCCi7Ap/vj6UqzfvgtTZq0qqVwCeZyLyE4+OUc/4f/niXn4yT/fx1Pvxyt7dai8ErhryIAKXaQeKiSp2lcfmFF2OVtaOzBv9daKfJuurgNIc8MjZmVFQF0IScceYwy3TluID9fvKHlBLFWeKg4ivks1HUtbdnYAAFrbC8rnJvlUy++h1mGqu7KNf9V0AG5Jn/rtK5j4q5dS3w0zJnIAXbO/K1WtSoxuSwBKRNLJvm77Ltz89AJc9Ls3fFFBldcrotI9gcuZyH78lvCzu19cjBGTppQV7rarLANdWSdZtVAQbjGL1++oWBmc0HdVDiBNk1sxr0wFlACWAJQIlnSyu99xV2ex6scyikpSPxZQ/GAqFJly4Uhubqlv8O9dqxiupzDFJ257CeMleW2trUB3FxFQVy+3q3MAabY5KAJKL18vz/SzbAyUOtkZ86m66nsu29iKHz82O3JwJ12AxUVfzDZOaTnqh1Px6MyV7rvh/EwRdf5BqcHpZq/YirXcYqOLrANxC1IXpg8VQyUc0DwdgKEn8MsL1+PyP7xVNS4oVR1A4MpyAF0G4kdOOrCikl/553dw/2sfYfaKLaVWLbI8UR6fhIilIQJSPnP/pqHgqrUZaFdRRqugq9rMZZtTDaoWsgKqYJeYjt9L7nsTT8xZXTWRUaXabDmALgSRACQZV0TRA6TDPVYuW6Hwf8F6V3ZC8AUxqpxaBsZLG11VJBGF829/BR/75XOp5RfyA0gtZx98zJiagVabMFsroAaA+GGSi0T0cni+iOSy4WeljitxYRJDQVRaaRml/A3Vp4z9e1eJBtqldQARfdTWkd5ACHMA6fYJA/PKSLqjr5aVVpr7ALH/KrFJsgRAg85CEfNXb9M+L5UDAKJ3RZ3uqpyN+NpJx4G4y/Z+g0pesEx1EDx3v0xV3Xja8mdNrbmILrz+1+5IyJTyFcecxwEY6gB4qqp9n4qJgCwHUDXc8swCfPxXL2LRWjURCOyqDb64mCJK7MHzTdPkS0WsiCovsjDb/XExUTnlhO+t3dZWlvdyKejKIqCqrX3Sx6jEN+Dj2ZQDMOFE00TllMDpwxIADWYu2wwAWL1FHRskKAKKzy+QJiI9H9Rp0npx0Rd345WWjZrsvPizNBYKLkZas7UN4382Db96ZkHZeSZBlyYA1fIDkK4rIRbjWSbt7+qJgCpTjuUAqggv4JTmYwatgOLz4+kZixaJ8EEdlWXScRDUAfj6h8pzAM5fEyVwWVZG7t+2zgI+d+ereH6+c7jcs/Ore8ic9QQOf0dTU80k+fM5mZS4VOv7pFmKSLgzFVitbThoDbgVjm5nGrCnNxhYQRGQlFeRgchZlDkHkOaOreBxFb4fAKEcXwZT2Ws8MZP1BOXgnaWbsXRjK+asdA7OqMSOKQpdmAGoon5CEgFVoOCkOgCOanVBpayAKhHqxHIAGnAlrG6XHNQBxEOMwSMvoPv+cCq+/9dZgXyjFpMduzpx9d/fM44uGbD88c1uEi0KpQxpmQNQrcdFr72lTxqxfYCgR6kyAejKIiDVF6yGjqQSfcKt15JaASX23i8RqXoCC7+tH0AVwZWwul2yuIgbcQCCIsqXyftf9JG3lwOAF5tflScv895XluDPby7FPS8tji0X8NtACMYCSmNymoxJE0ewspTAmhj0lfKl0CEpEVu8brsXbK/SUFWtkvJ5jqR9ouMuSZEmaTjoqomAKlRMJU68swRAA752mIiAkugA3De06TwOIGJsJz3ARfQ7SEMHYDrAfaIXTwHSEHl55YFzAGVnmQhJ+/PUm17Ap377SoVqE4SqZpXYnYeUwAl33bphIN7mhKvr+gFUSARkOYDqIedqXFSDrLW9E9c//r53bbJ4qXQAhDCB8XQAEUSCPzGVCYplBK2AjF5XlJ+eDoDXJ00zUNWZB9VAKRN/8brKRc0UoeQAKkEAyuQATMCrnbT+1RLQpaoEFqOBWg6geuAiINUAvuvFxXh50XrvelNrB8697WUsiQiBK8rBxSw7pK2+ZwVkMIo+2rADIyZNwbuuyaq+bOevQ3CEskpWAidL5+k/FAQrDSVwyPbcva66CKgrh4NWLEuViI1Trh+ALrVaBNRFOYAU+1WsstUBVBE8EoNq0Mimbf+avQrvrdiC3z6/SJsfH7QBM1AK58UnpVIHIP14fsE6AMA/3lkR2ZaCigMowxHM9C0W+qtQRDJ9e0utD8+q6krgrmwGqqhaJZTAIRFQSjoAESVzAFX6PBXTAVgroOrBUwIb7Op8wxr9BxIHhbjY6eyko8a2LwIyA/MrGKhrxR3BDGz8ReV4qdAtZJU4QCOyHl2YAKiQFgcgjiO5CyohpknqCezlXS0CUCFhk9UBVBHcDFReXGav2ILfPBfc6YtpGGO484UPsKU1aN2hi4nfqZEbmEXQNBsR4i6s9HDQ6kn+0cZWjJg0BS8tXBd+R3q3EiKgYpFpOYxqK4GrHXoiCZQcQEoropiN/A2SWwHFp+FEJS0R0LKNrdjcmuxAouhyUsvKBoOrFbIaM9BHFeIW8YO/smgDbvjXPPy/x2ZLacIiIIC0uxjV7ly+ZTogAjoAIQ/PQYyAye+ujAx+F6iHMMnf/HAjAODRd1aG0xV5+oi83EbNWLIJq7bsNCpfRGeRaSdcVEC9SqAr+wFUUgfAtBclWAFpRkvwQCLnr27zpIOOAJx443M485YXE+VVSjmlQMzJKoGrCIpxBBNREHbV7QXn8PTtbUEOILBLEn53uDMkJ21XTcQm67ap4xTJUFoBEQniIOBbf34HH/+VfhKI3IZKnKXabfPJHDUh+KPrHn8fJ9/4fEQr1CgyphU1pSkC2tVZiE1TbzqApJ60+rz1+aQlp1eNuTR1AGsN5xLgOGLKHL5pOUkRVAJbAlA1ZN2ekRevqGMNo6CLyMl1ALLFSrQOILq8bW0dOPK6p/DqB+sDZQf8ACByAKqzB1jktVAZAOrBaSLfF/ulPel2EdG7WN36f89LzkH0pv4Uz89fiwP+3xN4Z+mmyHRdeP1XoiKOYNJ1Wrth1UFGaegA2jriCbuME298Dodf91REOVYHUPeICwUhggkLuoytbR0YMWkK/vDaR14acYB48f9DBKD0QTR31TZsau3AzU850TD9WEBSOOio3bth8R4HoBhJsnxfSTwj8v7Da0vw5JzVkeUXikzbVzoz0BufmA8A6DDcAb+00CGkM5ZEE4AuLQJScQAp2a2KWUcpgZ+bvzb2+Mkgp6wSgzLP3JZztne9+AG+9sCM2HqqxsmmEmT/G3dEv5OqDkDo3UoQABsMTgPfCkjmAMJfIWqxXrXZGfB/mbEMgDPARdFLh5YD0OcZN8B4Ff0F2L8fOBwmYJ9vrrwLhsFw/qo5ABaohzov/bP/fmwOAGDJDedo0xSLTJuHTkneHhFuQ5mP+zeO8+rSIiBF3UtguNR5i4t2xDi69L7pGNKnBa9dfZo+L+F9xoSxrBABcQ7g51PnGdVT9b35Yt6UzZTEgapQKSugmomAiGgiEc0nokVENEnxnIjoVvf5LCIaKzy7l4jWEtFs6Z1riWgFEc10/51dfnPKw8rNO7HDDbDmWQHJIiDFe/6CHH6qUlSpzEBlAhC5lsQRAC+PsAzeVwj7kUFVjYoUQSl1AHoLn0q6YIpK4JAOIGbCmBIAvhmIS15rK6DOQlErflDdNlGixokzWts7Q4t2sE7BG6sScACBkOtgofvy5mzEpCn4/Ssf6vNW3Nu0w5Hl9+6Wj6xXEqTpECj2RyVMGmIJABFlAdwO4CwAYwBcSERjpGRnARjt/rsMwB3Cs98DmKjJ/hbG2BHuv6kJ6546jrvhWXz2ztcAANlsEj8ABavq/lXZ+Yt3OjRHQKrmnR9aIXpSyhyAGA7ajwUUFA3JkBdHpvvtcQCK+jJ1XmmiyJi2PxTHKgffNZyoPJu49b2W6//arW3Y75p/4Y9vLFU+V1XNpP1Rbfpoww6M+dGT+Mv0ZdpyEjuCaX4Hw687f1U6gIdnLNfnreIAXBFQn27pCUPSHO9iE2vFAYwHsIgxtpgx1g7gIQDnSWnOA/AAc/A6gL5ENAQAGGMvAtiYZqUribmrnFjyOg5ABdGcUrZ1Vw1SlTeubLFiYjnjX6vFVPICHAgHjejde6j8gGxW1RYVB8CU9U0TnYIISK5CnBWQ6eLk9WecCKiGFGDZJseE9m9vqRdA1eJnwgFEjcOVrnjzEaFMY+MBDXT+Jv5YNqubCqrP09buKIGbc9lEeVUPfqVrRQCGAlgmXC937yVNo8KVrsjoXiLqp0pARJcR0QwimrFuXdjZqFLI6JTACcUlInHgMDEDjSQA+uICVeTp1CKgZHqGOHZcOThlAhRT71Lg6ADU7YibMKYLtkoOraxLDXUATa7Zms6ySckBJLRek9Hb3TVvFBSpIQ4gsR+AumyVSXHS7mYMeHj6soC/i2oTUK4VT6ocgNh/NbICUhUrt9AkjYw7AIwCcASAVQBuUiVijN3FGBvHGBs3aNCgmCzTg04JrAIfnKpO4BNSXIxEz9hSdAC6Hb8iIYDgIPIPhPHDQasGrHxP1w2+ElhRvPTXBLdOWwjG9Iu6jEKEDiDOEUxXxuJ12/HFe95Aa7ujD5J1KjpUggBs2L4LV/zpbWxriz43IOfKu3ShRZQ6ACFte6eGcEQ0iY/pzTuEuknpk8cCir6v4mBNUWQMP3hkVsDfRZVHuZxcmoyg2H+1Cga3HMBw4XoYANnt0yRNAIyxNYyxAmOsCOBuOKKmmiG0qLp/kxz2ogIXAYmLEQs8NzcDNR3vvsjCgTiIRAskzztZkYd8elLczkslbuHpkuyobn56AWYt3xIw0Xx10XqtDb4oAjKpkwjd4vTzqXPx8qL1eGXRBicfSaS2bGMrPlREfq2ECOi2ZxdhynurImXbgM9BdiTYcovt/9z/vaZMYxKSZJtwMl3ICsiwT56dtwYjJk3B2q2+kjgoAorYzBhCudEphudAudZcafoBiHWu1YEw0wGMJqKRRNQE4AIAk6U0kwFc7FoDTQCwhTG2KipTriNw8SkAs3VpqwHZJpxfyTJ8VTwblZhHfibayYvHQ+rNQPV1jZtTGUlkIR4CI5btTSZFfiElcIAdD0MtAWKBegCO4jDO+abAWMAk7wv3vIFP/fZVZdrTb34BLyxQH/4et2PSLdjt7jfJu7tq3jae/MQbn8Mpv3w+9F4lOABf/BTHfTh/9aaM0TtdXUhxsYu27+rEt/78jmc6adLcsBm1Ot2DrvJ6plAP1aaDm1MDya1tVPVVDYFSrHjEEDGm+4B/vLMcIyZNifToF+tcEw6AMdYJ4EoATwKYC+BhxtgcIrqciC53k00FsBjAIji7+W/y94nozwBeA3AAES0noq+4j24koveIaBaAUwB8J61GlQJZIaaLOa4811YQ6YTzdQmA8KKYpScCClkBlb6Y8HrwBVhcE25+ZoGbxt/pJBYBKdJHegK71wXGcPL/Po8r/vR2ZP0zRNiVwEPz7aWbjdOq6iejwxWH5F25epQfgPid0rKrF+F9y5jhwMdpEhGQCccijoOH3lyKye+uxO1uMEQTLlVOo9PL8L7eJYiixDeVC3UK4iX//I3SRUsL1mzDVX+ZKZZk9N6f33SI2QfrtmvTiHWphBLYyPbJNdGcKt27U/jNAFyhefdCzf2LzKtZeYQ4AMbvG9hKRzxT6gDE51oRULJylPWS5KaOnbhzT+QGlHLQhBxAlA7At9t2rl9URA4VQSgtLIS8oCn7MLBgq3uUf3OPAJB+ES7HKsUEHjcX8/V52brxqnpb1f6QFY9C5BJlFit3gVyGbhfblFMQAI3hgVwfU5hsdIDkIqCd7cHNiikHEHfsbJK8SoUNBeEivHioOQAVVB+QjyGfAxCf+RNJJwJK84xc1aErjPnKU5OJfN+rSyIHqtKU1Ntdebka1ZsI2NVRflwgVR+u2eqz27qJ7llmuSIgXwcQvYBUwhHMP5kuOl0sATDkAEILuEp/5InEFO/LZRgaLHArJlEZLVbPVHwTBVMCUCyy0KKeqBzDiukiDuvyssHgKgjZfI53e1gHEEZUrBvVDl/83gXhuSrEgglCuzbPYczNK7QQO88iB5707KMNrfjnrJWhfDh0hHL2ii14aPpS7XsqEKgkDiAstgoXuFhgt+WJurO9gI5C0dMBMOm7qqovllGJUBC+BVJ0Ot7/uvhGKuKlJADSddAU00EUQYw7EjKeAxAWXqb8qc07DqrkqmH2z3dX4qAfPYGFa8zCo8vz3rRWvB+j9X0iATDMOAFsLCAXHaHdo/PXJGKkaFsvUwg+IYM6AL6wUOB5WuIE3yzS5WIkTsB5lpz13NbmWHuo48qod5OfuO1lIU/O+USP5DgOQKewNBEBfSBY78jPD/rRExg/sr+3i+afPmoRFu9VJBREJUVARrofFv4tKcWjypGrE6cDaA/oAIIcq6quphFdeXrdPfHR/724GACwYvNOjN6jV2y+8ng2DjGiOXRKl1clrIAsAXAR4gA8ubmBCCjig/MJKZokiqk7BUcwcQErZynRWfCIdxmixVtKsZaXXzi9ejcZvyCrkCFCa0dn6P7T769BLku49L7pyvfkOiuDf233nZZUdX7zw40YMaB74LkvhokWmaS1/m9p7cDvXl6MQ4b2wXT3wJ14DsD5qz9gSPVOvAgowDV6679+4YoSIQF6AqDSAUSJgHIZZ8P0b3e9rsxPBdVc9s1A/WfLXa9q088Z4gAMX8wa+BoFQ0EYVigBLAFwoTMDNTkPIGqkcIuSrIIDAAQ/gQyFZPQA8PislRjUsxnH7Dsgrgmhd3l2KocvRwegr7jykYKQcKgWnjiLEG3ZYErHpLiQv/Jio1qgdI5E4g6Tj4WoeEiqPNLyA7h+yvuhkA7lOqEZc20R0Tx5HTIRHMDKzcFT3eRvoNvEqnQAUUrgXNYxYnjro03qDBVQcwD69KXqAUwJgC8CipqHleUArA7ARUFYiAF/4JqEjPdPBKPQKuGLePx74vcWzUTlxWlbWweufPAdb5djOrB4Ot993s2zGExj4uQTyFdRf5P0cr3iXNoZC+4ETSEToTjTQXEB3CmYnbZ7IiD+Xf16ReWXlhXQToUJbFzWseInxWMzoh3+HaUDuO7x9wPXSTkAnRJYrmouk0nc36r28vqpsmotkQAUGcPsFVs8T3IdogipV78KmwFZAuCCi2r4Tl3LAShWryjZPV9MfOrthzmQ4/OLH7tYBA69Vn/qUBR8HUCwTvKOKkp8qpQPRwxGI5NCwwlbZKwkAmAiAtIt2Nvb/Mnq6QAknYVaDh1dXlqIyzlOAa3kXkwIgCiWdH/KjnFRMOUAfD8Af9E9+mfPeEr7kAgoS0blR9XFyZdv9MLPVIR4zdbocNaAEyL7E7e9jCsffCcyXTZCtMhhzUCrBG8nLh0FGY67E35XFWyNwyMsGVKmESdVOW7uqvrwHHwRULBcU9bTu+f9NRUnyPXSFhdKV5IfQIjgqPP2fwsEQAhnwMV2otc0oF5EK+0I5pcT/TyWATDlAKDvQy/mVQmRcjniTBnl7/7yovXKeuUkkakJxLpsdgPYeXND0Rc7pR3820s34ZifT4sth3P9ry/eEJlODjh57m0v4zsBh7LKbioASwA8dMocgNvvRixYxC6QEwD/PgUmlerAdrF8E+hELfKBMCEroIQ7D++e4pmJQtE0KmihyEqyqDHRORQ1CzYnALmMb5n14oJ1mLV8s5B/ZTgAE84o1goopr+SWG4F8lVsbpLEyQqLgNTpeF6y9ZeXXCpK1pnp8P7Krd5vkeAdcd3T2Nza7rVJxQHIIqAFq83MQnm9dAH2OETjgoVrtuG9FVvwDyGkBJBsHSgFlgC48GT1mSC7b7IOibJi2eKG51tQLPSOCMi5R6Cg0416v2nUFrkORRUHgGjiprZ4CeYrwkwEpK+z/F4ass94qx2BA3BFQJ1FPw7R/a99hE/+5hUlByXWlaPUOhuJUmLSlFK2ymQ00pdCGLfCZSTCIiA1BeDpZNGfjttwdADRZbd3FnH2rS9513IfbW7t8DdHirVaVgLrmBe5bl5YjpgK8gOLCkWGM255UZnGcgBVgqwE9ja7Bh9AtG+Xk/NdgDgYRCsdcdEPcgPJ6h+sDy/IyfPWZ53YLeIh4KInsDIPxYSQOQsRJjso07NSiywd5ZeSiymq+1uMaClDtqrSlVHqZDU5mCVutY3VASgeq5zG5DtRSuA0OQCeLuAIBj2xyWfjOQB5DMl1yWb8cCiq8SaPX53/Sqgcw7HL+zHK1NwqgasEfjRjTlLMhILBKd6N0gHskuTJgD+wGBN31SxIDFLSAazf4Yc+kAeaqQMKR9QJX6qBukOSoUZFTZXTpeFV+8biDfjt84sC98RqimWISmAZnlI9Jhhcscjwt7eWY2tM7H4ZYkCyTTvaSzpcI26jonqq4gDkbFRB0kysVzhMHcH4WNT6MUjXjud8dNnhnXm4vXwMmGxg5A+jkxIYEwC3I6P0XVYEVCV0Sh67vOOTiICc9MEX+CHz4gATFwiRBY0yexPrFAeRexH9D4JcSPSuUfUoqnzVoJdloHIS3QHejLFY9tkEW9s6ceMT86U6+PmKdebOPyr4C3T4mVjN2Su34vt/fRc/+OusRPXkbb33lSU48vqnsXRjayhN3LePU0CrCISaAEibBOGS/0yiBA6HgiBsae3Qjo3whktdFrcakiHSF3l8yxugjMABKJXAkjOiTLx0HLHp2M16HID+41kRUJXQKVnr8G4XrYFmr9ii3L52Cjtb+YNx5WLBG3wMW3d2enmKu8ugyV06HICYS9DTONoRTB0mgL+rSK8Y9DpTTj6pr/3n+8rnRRav1CwVujNnF6zVK/jE0OAX3BU8OEXsQ37OwYrNemKiAh8bz81zzjVQEQD5W63d1oZlQrrwYlcM2KGLT7/xx7cwZdYq5c5T/uyFIguZYpqeUaCqFxFw+HVPhZz6eDp5kX76/dWYNneN0gxUhYCOJ0Y0I84/pRloSAQUhE5KYBwKwl193/hQf2S6JQBVghy3X7agefDNpfjEbS/jJUUoY0+0gfAE4qIFcYBt29Xh5i0SGGmBLuO7izJrcQB1BHQA6YqAVOlLseXneVVK9il+B7GMqMBfPFlnkeH1xRulZ34efPOQ5FQuuU46yCnG/2waTrzxOf+5lMeVD76DMT96UpnBv2avxhUPvo2OTtU3DuKxmStw6k0v4Ln5a33/FW9XHltt7Rh7YUFwHvF+3LIzKD57bv46fOX+sAd4NhO/dMXpAEQuWDXeZB2JvPfjr8ivmoSPAfy15l+zV2vTGGZVMmwoCBedUlROvkvnaybflc1ZsTX8rvCV5IWQKxd9tpB8DgBMYCNRsg5ATsrrLOYvpssQtwIK51UsMmSkyKQcfKeqkoWr2N44MzgdChUkADodwMYderm9d3BOjNkkF40kFV+ZcDs8TVtHwUgH88Sc1W79GIhI+c1MrIDecQPvLVyzTRABAZc9MANPvb8mtt7yoqsTtfH26bgneTxyXV0+S9oIqHLZqp06i/i2snJeFgHpfIVM564cAl6FNI+XVMFyAC7EuPx3vbgYU95zTrTkg6h/9yYAaoWNuLOWx9F2vtsXTHO4DkDchcsioHLWP68kph6M3I1e9cw/JSycL4+YaqoEli06OERzWRUYS0cJrMtb9TvSJyJCQRkksM5FJTgAPj6P/tkzOOhHT4Tz0IgPuSLzzucXK/KMFwFx0ahj4ca5ZBgt/nK9ohCXTH7e7IaOOGBPfbROWecr78wLRX+jofoGcacB6vSE4ji54/kPsMAwrLQKVgRUJfz6mYUAnME9VWDJ+KDv5xIAFUTnlZAOoE3mAICtLpvLEGQjA7vTUg4mlerAmHpiZTJO4arBNW2uI4dW7YiiFjYTJXBceo5isXLmb2K3FgK/9eXx8MCqbxII7OcuMKYiALnsKOsoXvY2yVpJFlU6v/3n23d1orNQxJtLwnJmpQ5A4hT8UCZ+vnHhvEWYLmBx31uuF48dFNXXMjcR9nEQ5p9iqOo4Czm/KE7jf56Yh2v+8V7k+5FlVNC7HLAEwAMfLLqDWfiAU0Hc6fryd4e93K6wL/fi6jME2Ehx4ChttA3XFdk0UUaWSOsIdvkf33LeU5QVZa2QRAkMRO96p81bg/eWb9E+LwcBK6CYfpKhWg+CkV2Lgb8iomLImExynVjJMzMObB6YZ6p5zT9mR5x8prgv3RL1BKIDYxSC4y86reod9fPgNbcC0vVLW0cBn7nj1cC9cLBAn+vWcQB/mb4UD093zu7ViYB0jmAce/fvoayjCaNYaQ7A6gAkNOWCUQY9FjFigWgv+JNQ/F67OouhSdZZFE0cmTCI1LtJjufnr9WWHwpbzPhfpiQanMhFm4GqOACmfaaSh+o4AEL0wsAPy04LXA4OSKa2woWJGEbNAfi/ef+oFtbP3fla6F6ystVptrV1oiWfDTwvFBla8lm0thfwzNw1+rOPpe+zdmtbiPB3eByAb3cfF89HzIPrIuIQ1wdyvXwOoIiMwJ3wcaiK5KnSAegseQCn7f/1iLN7//zRw7VK4Dgz0AE91dKDJE6mlYLlACR0z+eUct0oAsBFQM6A8u+rdv8FYXEIcABSGfIu8kv3Tdd60spVE801Vbu2bMbhAFRja1CvZqeeioe8Tqpa8EEfJID6cLqVkvGroPqegCw2ia+PStwQPEfA5QCEv+fe9jKen79WadrJYSLu0okj+BiTOZtu+axfbx0BkLag438+DR8IR2YCgggIflvjOIBSFq343XAwz2aXA+gosIBPgOrsC69eUj/s2NWJafP0G6tYsZS3NsjvBW/oRKdG541XeJo0NAEYMWmK5ynqezcGRTFRdsIcnA2XD1lRiUA6i8FFx2cGgu+qJrxuMOisEEQCkxMsDrKZTMjslOOI4X21ZXUW/HxDdVNMvCgdwCE/flL7LG3I1lWqk5hMRBWq/hK/p6cDcO9tau3Aeyu24Pt/fTcyX7OomuoKembGCg6AY74miJlKBzBnZVD0xr8htxwD4kWRpRCApCIgzgEUikECsHDtdjz6zgolsV6wNkjcrv77eyGdigiZcIZFQPyvLAKS8jG0UlLBcgAVAl+wbnxivuSQFRxsJiIgvtN1REDRC2Aw4JuoxJNs9g0sNLy26AgA/HZlAwQA7nkA+l3SDgX30hFhBvru8i3oLBQDT7SOYBU42SgKsoJUDvcBGIphlGIxv43cGqyjUMSDbyz1ApEViixGwRtf9qwVWzB3VdgEWfQpEfNrzvtTWyd+Uo2xDcKRmUDwPAvf6iXZYq2C3B9xfSCX6YmAikXkJaewq/4yUzn3/vnuysD14nU7QmlEyFy4zhFMJl4ysdbFejLpp0ozyg1LAArSouD/ZoEFTqfoEeFxAAiKgJR21gFzvaDNvphcpXDV7ZJ07DMTCFJOcJzJupNZzq8pl/Em/Df/9HYov7igZdPmrQ0MWB0HUGnb5nB5/u8iYwIB8O+bLMLyTo6x4NGV/Ht3FBgenbkC67Y5cZg6i8GQHDJMyl68bgfO+vVLofucAyjKHEAuG0orQ8VlcosnL02naAVk9t1M2iP3R9wrUUrgnCIsxLrtu0L3kqIgf2/puai/EyFzDmI/dxSK+HC9Q3h0/fSLf83FXS9+ECijUmhcAiBNGA7ZHJP3f9Sg9ndG8SIQPjhyriJWjAUUsAJSlKerwc6OTnz/r+9i7bY2tx7MS6/kALIUcjwDgJZcJtLUUycC4jbZTnx1/6GOA6ju8i9/X3/BiDoMRAX5ezIWFKOIFjNymOhMhNNPOZOcjyfxW7YXiujRbEIA4uVefPEi+N8trr4m7ZGdoJJG9hTNQPOKvl2+Sa9zMYU8B3XGFiHlsnQtbub++9HZOOWXz2PTjnYt1/l/LyzGz6fOU+adNhrWCkirAJRk8VGu4jLkEMtRUf5yWQos0Eyqh4oD0E2Sqe9xr0/gps8f7nMVAkEK6ACIUFQcCdmSz0baPvMFQ07x1RNH4vbnPgjFTtFxAJUe1DIC37ro7x5NuDsR8hGBRYkDEDkkUbbcKZhlqlBOf6iUnsff8KzRuyae2t4YJgqIK6Ng0pykBEDeTDR5SuAi8goT7ajgfibIkGJhl6qoGz8hDkC4ftZVOrd1Fow4YSsCqhDEjyRzALLIQPwbhWLRzAwScEUyghhGViCX4gfAZaHCPtTLM6gDIM9PQYRDAPR19giAVI8+3fIAgFYpTIHOCqj6BED8zbx+ivICVaEtRAAgiYD8fLYKMW2KRRZpOmlKgHq3hPdrvA2tu5IfYJ7EY1k03Y2rrsnCJouA4saEPJZEJXBOyQHEE4CoEpsU3LAuUqrc3FAMImFc7PKU6mRoBWRFQKnjidmrsFZwzJGtRMRO57JynRlfkyB/lCNsRk0wzgEUPAIgiYAScAAc8s6WCcRMnCT8NCV5ADbnMpFES/RfENGzOY8MAW8t2eSd4QroRUCVHtSh8iQCn5V0AKbOSm2dYQ7AE5FQ8JttbPWVqZ0xBED+1Kru2W9wTwx0TXRFdBYZlm9qxW+eWxR+KQZxnq4iiHzlv+77/cep+7nP4/PLSorbuCGh4wA6JSsgDq5/iUKUKUJzLqvgANSinrgDYUTOkBOyzmKYA1ehSxAAIppIRPOJaBERTVI8JyK61X0+i4jGCs/uJaK1RDRbeqc/ET1NRAvdv/3Kb0482joKuPyPb+PCu9/w7oUUs5L89k9vLMUfX1+qzE8MSytzD9EcAAXMQGXCk8QKiIPvivhOZduuTnzitpcBBCdcJkOhA2gAnwPQOg5pOADndCZHCSyG+d3VWVSKPqrMAEgE3V8worxAVdjZHvwmRcbQXnAmdEsuG1ASt0ln2yYRAcky5FGDeuCwoX2U46lY1G9O4pA0ZhHvJt0On5ueGukABIL4j3eW462PNkWml7kv0cpJRQDk9EnRrOAA5I2CzipK/p4dCg6gUAhz4CrU3A+AiLIAbgdwFoAxAC4kojFSsrMAjHb/XQbgDuHZ7wFMVGQ9CcA0xthoANPc64qDfwBxhyArCcUPXWQM0xVxVDjEwSc7gok21DJymQxa2wueaRpDcCf45JxwsK24IxV9DiBYPi/P/00B7oCjJZ9BR4FpCZfOnlkXJuOjDTuUFhrVBu+Pne0FLF6/w+OGkuh3AGCXQgTEFb8t+UwgKKCMrRH25nHRJLMZUookAGcn2VRiH6/aog9PIaNYFEylNd2VhACISvHv/OVdrI4IlQHoOQBAfTaArK9RIaqWze5cEBFSAmvCo8s6AJED4Gk7i0UzP4AKUwCTkTMewCLG2GLGWDuAhwCcJ6U5D8ADzMHrAPoS0RAAYIy9CEC1gp4H4H739/0Azi+h/kZgjGHhmm3oKBSVcmnZ6iculo4I0QZZ5h64Aq1ZYZInD9qd7YXYXVCcqKKJ6wAUVRaJUCajjgXEOQCd7F4XDVR3OtOarbtKDgltgqtOH43jRg2ITce/ydcemIEP1+8IOYKZmqWqlMC73G/cLZ8tWWHHFwLuHxE+R5eQz4YXJMCxOa8EkR3Sp0UqhwkiIPU7LXnOgcbnH2UWq4I8jsQxl1ecDaAKBZEE8pxVhVXh/SB/L7muqu926k0vYOvOjlji3RVEQEMBiMFZlrv3kqaRsQdjbBUAuH8HqxIR0WVENIOIZqxbFz6MxQTfe/hdnHHLi3hs5ko1Gy2KCCRb/jgCLA5EWbHKyxLZVQ5ZcbVi807c+cIHkWWZcwDhdCoOQB64TVlnl8nrPViSOXd6VkDB93QEoNIQ5fnR6Zz6cv2Ed+pbQg5A3tmxom8nL3reJsW8VduwXrBbl6uTIYcDUDoWFlmC2JzmuOacgwKcXYHpRR4cWS/eUnx/rt7ahicN4wQBCg4gIjgjUL4ISF6YO4thkanOQKRVKlsXRHHx+h0hJzYZNRcBQa0rkatVyjnWRmCM3cUYG8cYGzdo0KCS8jj70CEAgE072pWKybAfQFAkFIWgCCj4wTjlb1YM1lJ2bbHESNIBiMhKZqBAWAaZzRDaO4teH8mLms4PIGoQ//qCI6IrXQYYY4YEIHjNN598XiaJS9RLsMRxdACcyJdOAH42dS7OvOVF71r1XfLZjNKs2AkumD6XlSEK2NcXi77nuK63MoqAe1H4+h/eMq6PvKCrNh3iOJRNkpNC3rQ99OZShR+AmiPa2R4U9+liMbW2F2LXga7AASwHMFy4HgZgZQlpZKzhYiL3rz4qU5k4cf+BAByRjG4XxSHL8YtMv8Pq2z0vKYGDitx2V5Si2h2qTNfiEG8G6hIA1TOBCGUyYasjwBFLdRSYJwLqJtW7QxMMLooDKFU+DQCjB/eMDKFQZNEetl46zQT0T/oyr5MoWhD9AFoUXF4SbNzhWw3J3yWTITRlSWsZprPmidslR4EQ3KQUGBP6S10e/xQ/fVx91jNH96bkxFKetyLhH71HT9z0ucNx5xeP8u6Z6ACiIG/a/vuxOQo/AOevTLBl8ZPu+7R3hsNYyOgKOoDpAEYT0UgiagJwAYDJUprJAC52rYEmANjCxTsRmAzgEvf3JQAeS1DvROAywg4NAQg6gnn/uc/0+c645vSgGSgTB4WoAwh3s8nOVUacrJoPJtWgyUtmoI4OIJime1MOHYWiZ8EiL2odhSJWb2kLxbaPIgAq8ZcpmnIZ7N2/u/Z5kekPBw+mk1j3ohvcjIuAEuyyAh7kjKGjUARReYROhvz5suT0haqanUWmVc43l1EnIoSibPJxpVN2m5xxC5S2+QkHZguW+5mjhmFwL19vUT4BCBMpeV7pQknLBCDqHI048WnNRUCMsU4AVwJ4EsBcAA8zxuYQ0eVEdLmbbCqAxQAWAbgbwDf5+0T0ZwCvATiAiJYT0VfcRzcAOIOIFgI4w72uCDIZQk4Sb4iI5QA0u8xcNhNYgLbs7MCt0xZ67/kiIJUSOPnkjBsLvoxWVVfJDFTiVr75sVEY1q8bOos+ByDXu7PAMOEX0/DSwvWB+/ks4dnvnYzffOHIULlN2dJFI9lM9NlTRcawZ++WiBQ8XfCauXknDQXBwRcfxpxdXFM2UxJBl7FIc3QgFwGpUCgw7YJcDvElosDuVDw+UT5M/kvHjcDlJ4+CrIu94pRRgWuZowTMFfDyIir6VfC+F8d4ucYHqk2bbKnEu/0v04PnV8iBFP1Q6QrdXKwOoLIUwCgUBGNsKpxFXrx3p/CbAbhC8+6FmvsbAJxmXNMykZcUnCLEcAlLNgRtquMUhOLEFJ2giox5xEY1mFTxS1QYPbgnFrphbE3jpajSifXMZXgYCj/daQftgdcXbwAAbHe9SlskVn1zq/rg9Hw2g30H9UT3pvBwKmcRymYo0oGKMWBov26x+cjfkPsC7FDE0jdBLuPI4//21nL834uL0as5F3tIiglWaswyiSIIAGOhoGUcqo2HKRwRUJAA8H6SRVE/PncMiAiPzVwRuL9D8k7u0ZwL7cx18nEZnUWG4f27YdlGx8NX7G7+uxxFvAyVTuejDcHooUXGMHvFFsyTwm2HRUBOf23ZGZ4/KgsmERVe/xvHE9ixo2ZKE0c+nlXU2FTuLqPI3DglWfUiZiK6cNIFRUxRaG3v9A6cD+UTCgURXBizGfLEGDzCZItLuI7fbwAuPnYfLVvNZc2qXXA5opEsUaQOoGdzzlPwR0HeeTEGHLVPP7ywYB3koziN6uW283+fnA/AaX9UsLdyQdDL8wsRSmDVxsO4TKLA4iSekyETAM4hy5yyfCBSt6awjur5+WaWfZ1Fhh7CBkMca3x+jRzYA9eff0jsgTUmUI1b2W+CMXW4kx2yEtgl0OsVEUqjREC/fmZhZDyxNNAwBCCfzeDp99coPwIf2DkFNY63AlKPNh4quCmbUQ5IVVkqiAt3XF1++dQCHHbtU7EcQJa4FzIL3ONt2e7GmOc7KgJF1pfnrZLtNuczeOKqEyPrrUMmggPo2z2Pr5+8L4b1644lN5wTmU9Yecfw6bFDsWRDK/4yfZm2Xyfs2x/Pff9joftyO5tyGRjS85KhI6SFYmWUwBlJB3Dbs4uwyeUAdeIV+fO3Cgvh4/9xgvId0XM8Cs6ZCuFFHwgSg4sm7IMhBmJBIFr8pOJc5fASRaYO8RF2DOPOp+2htFEbwVueWaB9lhYahgBsbm3His07vTM+RXgiIMXHiNsd7tlbLYLwOIBcRjlIVLvlvt2doGp5SV7PoTvZSVW2jAAByDr6EM5O8/o0uSKDrTtdDiDgbq8fqDmFDJajKZvBgXv2xkFDehvVXc5Xt5v78vEjI0UcARt2hfLu/COGoldzDvNWb1PGZLnzi2Px0GXHYuTA8IHechybfEo6gCjkc+r8K8cBhL/nzGWbAeij3MrjXBQBDe3brWxxhtjFKhEQh6l+LWpqq/purYIAmHx3TqBVZxTUyoeGo2EIQJSsMYoDYCw6aNRlJ+2rvF9kDLs6imjOZULKMUC9oH7l+JH48BdnB2Tp4i5FHoA6KB3BhPKyRNjVWQzEkMlmCP17OIdXr9ziEAbOATCwyJ0KX2xV/cdlqVGWEDpE6QDiFpNezUF7fbEfGXPEFc35DHZ1qmMfRZ1aJnMA+SxV9JSzfDajVaZHWgEpCOSAHuoDymUQSLuQ6mIIyWuhqAwlEr5ZiYRA7OKMhhsASrMykqHqO5nzKTKnn+LAx/4GJQGIf3/MkN649cKwgUUaaBgCEAVuBqgaOHFilwP27IVXJp2qeM9R+vRuySsHiWrnkHUXkrhDZeKgqrLsCSwjQ8AevR3P35WbHQLgi3YykTsVL52GAwCSBx5z6qRfWOO8okWHLVnfwX9xz2fVN1b5F3D7dfnb5bOZxKENkiCXDVrkHLhnL+/3tLlrcN8rSwAAv7/06MB73RT29vd/ebxRmY5pq7pNOpFTlA7AdA5E1knII+zYWHq+KphwT7IlnQ48hIocIBAwEwWfdcie+OThe8WmKwWWAMAXrShFQAYfeGjfsBiIMYbNO9vRt3teKcZoUuww+EAWdzSlEYBwnZsEEYJKYUlE2MOVna7c7Ci7+MKez+qtUPhzQK8DAJKFHubIZUgbRTNOb9tTIAAFxgIcIJ+0eTe8gmoSq+Zlv+5Nbr2CDystAsplMgGRljhOl2xoxfvuWcHyrlV1Klj/Hk3Yb3DP2DKJ9HofvQ4grBvxMxTClPu3tPjO6fsr8o//DaRDAESLojPH7IGD9wqLMIvMzIpJpzwHoDzMRkYlDQwsAQBw9d/DegGOUuWWRcawZWcn+nRrUooxWiKcw8TvXYppm2qzLZ8JLCNDwMCeDgewwuUARNFOFFvNiYNqt845gFKsGTIROoC4D8MPqQGCFiz8mtdt8rsr8cnfvBJ6X7VjHdDTJQDSRiGXjbZWKhf5LAWUwFnNwiwrLrvlw2a5RGaxj2QzUBF6HYDzt4fLedwheOYGREAuonbP40b0Cy/kQieT5jci6p0Eg4Q4WHddPE5JNOVxpUOhyPDG4g3KiKcm5uCVHFsNeyRkEpTijFEsAlta23HwXr2Vsj8Ve84XWT7wT9p/EIb2bfF2eKZQmaYFdACKCcIDjg3o0YQNblgCPjjzuUykYi3SE9g7u1VPAPbs3aKcHLkoHYA2Nwe9W3wCwCQOgH9OUwXcSfsPQv/uea+OIRFQpsIcQDbY/7pFQxZbqDgAAhnFDorybpbDYnPwb9VeKOK4UQMCnDFB4ZEd8REJYQ6BALx+9Wlo7ywGDtyRx4iOQOqw/x49sWDN9sC9UYOCyn++ORLhhI4PNqJ3S04Z+vvf7npdWbYJsaqkeNFyAAJ0/VyC9AJFxrB5Zwf6dMsrd8aqnX1W2kmfc+ieyQuG2g1eNgOVwW8NFkzoOHuaz1Cksipq8eNtiRIB3XrhkYEdF4fDAZQ2+EUdQKGI0KE/gLmZ5ANfHo9fXXCkx0XJ3FAuS0onOBk9m0vbb+UzFCTgmv6Wx5SqTkTQOo6Ji16GSLs4bdOcbcBFFR2FsHUMEXlE28j7l8JEPkPAnn1asPeA7oExLNOpJErgw4f3xakH7hG6L1t/qQgAU4iA9h6gD12igskmJA0nQ23eFcu5DqHT6Ju6q4to7yyitb2Avt3ySjl2VIA4nj6bUcd/iYMqFK7sCCaDDzIxBLSo3DX1W9AhSgTUqyWHE/cbGLqfjwgFEdcvPZuDIqAgB+CLgJIg63Fowfdy2YzyvN60kMsGnbJ0C7PMAaiCrhH0cusLx+/tiTocEZDCgz1L2LZLQwBITBd8l+D3u8mQVn35oNhHfR9IRgCypN749RI4SMAX/4lQiYCiYlepYEIAKikCsgTAACrljc6xhWOTy6I6SmAVB6DXAfDdjclAXnLDOVhywzmB/Ha2F5SWKnI5IngVuSWQ845bDynmUSmIEgFlJMsnv/yMVgmsEqGJ6CmFbZZDfgNhDiBgIqko14s5ExIBEXp3y4dfkFBqD8r9ryPGQ/p0wwVH+0F5lVE3Cfj9pWpLoHw248fr0YiAujflsE3jba5zzgLUOoAoqBY9OQCc6jegJpBPfeckZTm5jHqMyfUfJHAAP//UoQDUSmCdXxCHHA8pbhwDlgOoGnT9rLJ6OGRon8i8uNdkn+5Nag5AFSDOTcgJBg/ZYAKRJd7ZUQiVmQ8sIHoOYA9BBMQXzTgRkAmi5L3ZjMZ5LZdREk8i4KsnjowsT9yRF4tBay7PCkhq042fPSwyT95t8gKTyZCSA5h01oHBRbTELpT7v193NbHJZgg3fMZvg04sNWav3sqdalawusoQKcdJSz6jNGcEgvNH7luCT+QNJUCh7hKtnAJmoNIqptIB6EKLZzJ6zv+Z756EZ793MgDfSRMADhzimOE6G4tkB9XI/dLLgHOspI+hJQACdP1cigkjR59uaj8ApQ6Ai4Dcr+IEbTMrOxMgAMVQmSI7rzIr4xyCKALiDka5bKZsEVAUZN8HjiYNB/DNj42KDXQmKoGLUsA0zw9Amqx9uzfh5P2jDh1ScwAEKDmAcw4dgps+f7j0dnLI/X/ekUNx/fmHKL2UReiUwDpkBZ2LTgQUtRsVn8nlkCDTv83AqYkhLCoSRVw6bgAIf5/H/+MEEKlbruMAAGC/wb2w7yBHJNYUKJu8SsprQ5xERx7mvQz0QpU0MGhIApB0NytG+UyKvt3y3oIusuQmIiD5w180YR9tOeKi3tZeCLvHC8+3S0q8X19whKeE7SEMSB5mWHZEShsZUnM65XjYHjSkN04c7egVHB2AeDA35wCC3yBuovm6mfDi1rslTAB4GPJykcsGFbJN2QwumrCPMryyCP5cXLyiujNLPgdAREpHsKiD5MWFWLY0IvIV8eNG9MPEg6MNHFTjQdw0iSI+eYzI34dz6yrilcn4SoCovhE5OZ69SgeQNH5/DwMCUEkv84YkALKCh6MSHS3qAPoKu0TV5PWVwFz2HlwYe3czUzS2dRZCg11c7MTTp/p1z+O8I/zjm8XFgnMATdloT2AT3PHvY7XPshoOQNQB6EJu6NCcz+B7Zx4AIHxSm04HELdWe99FwQ2pvk0uQ5FE5ZJj9QRdRF7yxPb8M2KIMl9cxJ2z/EYgqqZgduvEAkr2zcWmyiofAgkOYNH9Ajie3iERkLBp6tvN19fI3cAvv3XaaMy7fqJ/XyXrJ9F/QT+/8tkwB6DSAcgE4bQDlUedezAhAFYHkAKO3Luv91tnjleJbu4rOIL16e4PWrUIiJuBBq/l5yqIA29ne1gHIC4WYkRUeQCLOx2u/M5lov0ATHBWRNhmIv0ZBpxxH7t3X3zvDMc71CT+ihieoViE0hNYVnIGJppiJ8e7P8QBgNQcgGRKKW4wltxwDn5y3iGx7eBtETkJvhjFLQyqMRbaLUvmlGKeSRX/Yt6ybJwI+NSRzkajOa/2m/jM2GHeb8ZUIiA1By2LNLnFWa/mXKAP1OEo/DEWdVRlXsFFqXQAosHIq5NOxY/PPTjwXI4sYGIabHUAKeDBr07wPBRNFC9pQSxLVN5FmYGK1ibicIlyCBFZ7p0dhdBEF80IN0ScPyvuijs8HUBlRUDZDGkimJIQkoESxRBrymUCE1UMmMZFETIHELcrFZXzwftqrtLhAMqfYo4ISOAAIsJvy+XLiOQAyPdoLhZZYvGVmFzeWBCAH597MGZdeyZa8lllX+87qAeO328AAJ0ISFyEg+IqEfzglf5y4DsVB5AROICIxTgQodctT3YwBBAYZ7kMhax85PkmlimeaSzCcgApoFtTFkNcz0Qt1a1AP2cy5B0O3y/AAUSFgtAvNDqIA6utQ6EDEAaw6JkcRQC46aYcioDjM2OHxShNzeDoAJx6nHPoEIwf0R+AszvjfSF2hcl8EA/iKTKGtdt82TXz0kQQAEUZPD+V+79KBCTrAErxJ3HqGSTAPK5THMGKMvflkP1DfKIZXafTDwqLNjIBDkAiAOSIfTinpKpbsxA63Ym1H3wub5q44l3eGHECINvuq7orl/E3CpwDUFlIiRuogAhIUgKLBIGIQlxFaL4JY3DiIaU5fpaDhgoFwQe7VgeQUjmfHzcMe/ZuwTtu/PSdrtlcH4EDUFmx8EnBdzQ5yQw0ar6L1gi7OoohtljcQX5u3HDc9eJiAOEBKcqL+WDWRQMVLVzKQYb8BefTY4fitQ824M0lG93zAIJE0RTNWX+XuX1XAd//49v+Q60OgPDfnxiD4j/nYMLIAaE8eQ1kT2sivU5HXOhKtSVzYjGJHEAwKun4Ef3xhWP29p5P/daJmLd6q7LPZDGIOE7E0BtOpMvgu09/5yScccuLAIBvn7Y/npm7NphXQAkc5gBEqLhZcewVGUPP5pxnTi0/BxzF++bWjtC82OK+M6BH0HtXJQLasGMXiJxAbz2acph3/URlv4kiIN7vHYWi105uHiv6u2QoTLRkEZAJZ20SkLJUNAwHAPgfLs754th9w5M/DiKlv/Gzh+O7Zx6AP3zlGABAq+s5KToaqeyF5WBwWckM1FRJvauzGNrpiRYdV591IF51Q1iHdQB+OzoEDqASB1d4i6OgA8gQeXXKCJ7ASSUp+ZwvIli+KXjOs14HAOw3uCf+8JVjlGOE57ezPexprfo2WY0VUJSsWYUwBxDUARy33wCcf6SvyB+zV298euwwr88CNZCqI1Y7QxTY3fJ+Gtq3G372qUMweg8/DLVKPyDmFeYAgmlV8aiaBA6AIbxRkxdTLl6VF+zNGg5ANX0+2tDqcwDNWbTks8q5KfY/5953dRQ9HcA1Zx8EIDifMhRWdjMGHD2in5Bv/MCu5MHwDUUA+GTUhQDgk/gbHxulfH7deQcr7wPA6z/Un2/P46KLMklVHWQroIxoPA1zl/D2QlFhGx2Un+rknSoroFw2U9LxgnEyZHFyePNGIAaihQZjyTxJ81nfiew3zy4KPNNZAcWbgTrPZft6nVI6K1sBMcfZbMq3/CMyuzdlcc6hQ3Dh+L0VOTg49/C9AgSGL0ZFjTmrV77GiU6HoAjI18J+ccI++PdjghZLqp1rJAcQoXyesK8v8vNM7BkL6erkCLrD+nXz6yqAEx9ZB6Da2a/btsv7firnTA5RBMS5vR88MguvLNoAwOfoZQKgwl8vPw4v/eAUfO3EkTg0xqEUSG5amgQNRQA4u6s67xPwN0e6heDiY0do81ZZgXC0ujtGcUCqjviTOYAiC7qBJRGDRFkBAXoiKE5s7gfgnEiVfKjE7W7ERYAJHACfwNlsJrAj5TDpBVH8Ii9GvFcjrYAU8IP0DcHw/oLLv+a1LAVjKDEAnx83PODA9f51E3H7v4/Ft07bT5nHA18ej77dgwsZJ1w6Tkaub+CedB0QMYoiIPh7D53sPNAHCHJpcpRMGeIc47ox0QqJsfCckjmAg/dyFs/F63cE7p/qml7K6VWft7Po6xqiYvOLojIx3xcWOIfa8zVFtAqiiOE/vH93XHPOGLNY/xXkABpSB6CfMMG/AHDh+OH45OFDMbh3OBqgKfjh2KJMUrU4yjoAeQ4lMcoIWQFJ5elkjyoroHyWtGfSRiGOa/A4BCaKgPwdXJZ8HUBSBapoJSKD92tfKaSCqSNYJkP49JHD8OtpC6PTyzqAiDboxqSKKPF+9cNa6zkQGeRxMVx8IqQX+kwMdaxaOPO5DJ7+zsmB3XcUB6CrW/emrLcIOiaZvPywtZ68cTv1wMG4+ekF2FfyiL7ji2PRukshplPU44sT9vbaZ2rpprLg88Oex3MASWE5gJTAB11zLhM4Vi+UTvhwv/j0YTh21ACMGhR/itK0752Mf3zzuNB9LgLi3rYfO2CQkgDI1j9MOsvWxP7dTwt86bgR3rWp16u4aHNdyOjBvUrSAcS9c8bBThjebk1Z8I1TgAPIiDbXiYvXnyXAZdv9gjtYUxEQISjbj3pL5LyimqAjlqoqcXEE3znrdq6qd/mtuy8+Cj+YeACG9fMtXsRQEMUiCzhuhetAaMlnA/GGgo5g0WcO8P7r173JC9MhKv0dEZDEAUgimkOG9sEbPzwNnx83PHC/OZdFP8XZxyqO6KfnH+p9V1MuVzVOOKcXFAEZZRcLk0NnSkXDEoCHLpuASWcdiM8eNcyLx8HX2lKPYBs1qCeO3Ltf6D4PntW/RxOevOok3PnFo5SDTTZ5jFOkRSFDhGs/6essZBEQnwxXnhIUPTQLSuALxw/H9GtOx5i9epdEAKLOVT31wMH4n88chteuPhW9WvIBJTC3esgKcVpKMaHUfUfercP6Bs394v0A1HnzvrzspH1x2LCgTDfIAejz1hEA1aLFy96zjxO4T6WUBqJ3oMP6dcc3Pxb89pkAB+DXV5WNyjEw6AjmvPyfHz9AXTe3oL7d88L3Dpqhyn2p2nnv0bvFeL7GzZ9yDB34/OoIWAGlxQFYEVCqaMpl0Ld7Ey4/eZR3/eAbS5VHMqaJXi05jwtQhZjmMtTvnLE/Lv/DWzhor96BXWOSUBVyG/IKM5olN5wTuheMG0NefUvTAQQr0bM553FD3ztzf+SzGQzp4+zCxcVG5ABEHYBpYDwO3XfkxIQvoGL9ouCL5tT1+KFrCTJi0hTvXsAPIKL+ehGQvj48cqsuPg9fGHVx9EN1EK1whPAZqnGnsgISFzzu1HTFKfvhilPC+o0AByCI/Pz2Mlx87D44YM9euMA9TcvkoPYo7L9HL7z54UbvWiYw5YQ85xyAuGnj3XH0iH6YvmQTgPBBMwBw0+cOD3GjIiq4/jcWAejo5Da7wZ3EdZ88GFedPhpfe+AtANGU+5VJp5Yc4EscwFEhmY8bNRCzrv04ABj7Acgo9ZxU3U60FE9gOa/H/+MEfOyXzwMI77aLCiWwaJbIDM2hpnzrBCzbuNN7XwXRCmj8iP54c4mzKMQtMILKInggScQ7gbDFEfXWEfeo3e3h7gK2Z+8W5XOlFZBUW96vZx2yJ8bu3TdIcDkBUOSt2lDwqu7RuzlgMqoCL6dv97x3hGI2Q179iszpkwmCSXbPMj34775oHA6/7ikAwEOXTcARw/sC8HVdpWxyOPj8UukAfvelo7FwzTY0ZbPKhf4zR/khMAb2bA6EagGSb3ySoKFEQG2uR65sApnLZjC4V0vAEuWTh++Fn3wybPY5tG+3QMx8E5zgnnYlH2R9xSlBc9NYKxQ4i1bIxV2BfaSj6UwXcJ0YxDQk7d+/eRy+5sbql1nqEQN7eIpXeXESlcBFQSTAVx9TMejBe/XxPCp1/SlO9HsvPdr7HcdhibvjjOGuWrQC+ts3wvqhcnDcfgPxyDeOw5dPGKl8rvKdkOt67bkH46AhvXHLvx0BIvJs51vyGZztxm86YXT4tDalH4D7sUzOsOUboN7d8p6+IJvxQ3+odr0mgdOi0Kd7Hvu6x14O7NnsbQR9f5fSl0M+v8SQLHyM9G7J46h9+uPQYX1i5+4IxZGSNVcCE9FEIppPRIuIaJLiORHRre7zWUQ0Nu5dIrqWiFYQ0Uz339npNEkPflSizhFH3HneeuGRuERQopaD331pHN790Zmh+//58QMD16o1NmAGmiHMvvbjeP1qvc8Bx6lSFMJyHblMxU9j9+6Hb5++P1ryGXz/zLD8l0/sMAfgl8N3UVkKeqYmZYVV/fm5o4bh0SuO966TcHOefLwYDFMg53DB0cMDznyAs9M9aEhv47I4OhSHEYk4ap9+WuJsIoMeN6I//vXtE73F8JpzDsL15x2Mk/cfhGP2HYAlN5yD/RW7eVW/8Z2qiUyeL7qDejZjnwHOoty/R5PAAYQ/dqlnKosoFv2NBkdnCgQgm5IS+MyDw+cTV1IJHNujRJQFcDuAMwAsBzCdiCYzxt4Xkp0FYLT77xgAdwA4xuDdWxhjv0ytNTHY5U4m3UDi/Zz22SfNuWzsASaAgRIS5geZD5NYzXKjeSZBz+Yc5l1/lvJZYHcvgCk4AMcu3X9P8BUzgmhVxSfRt04bjeFCrJckk140zxVFKfKB4Td85jDvZC7fosu4mAB2RRylGQdlKIiYzuvelMNFEf4ufj7hjHgfmxDVda6YY1CvZnzjY6NwxkF74JCh/g5ZpfAtlwMAEFA4c7QLQQ9LBW+zKAIqJbz8107cF8fuOxDn/uZl716pMaRMYDL6xwNYxBhbzBhrB/AQgPOkNOcBeIA5eB1AXyIaYvhu1cA5AN1AEkVA1cIj3zjW8yxUlRswA01Qr5CYq1Ka7YTwLK1CIiD//ogBPpsuOgZxmHYDJySi4i10nGMivYqfL//dvSmL72ssXcT8TSYxl0mLiOMAoqBqWhJT4qTwlfcGBGCbQwAG93JEMae4HOvVZx+IH587RhlsrnvMATgmEM2NOTgHUIoO4KfnH4ILx+/thYdQncKWBESEQyXldK1FQEMBLBOul7v3TNLEvXulKzK6l4j6qQonosuIaAYRzVi3bp1BdfXg5pi6gx9EEVC1cNQ+/TGwl7PriSs2SbVkAuAfZmKWSdxxg6WCaTgAUQn8g4kH4v4vj8dR+/QLOAYlxYCezfjJJw/GH77iH4Iuy6eTWVb5xIi/9/lxw5W71VLw6BXHh0IDtJfBAahQjaHdM8IrnoMTgEG9gtxT96YcLj1+ZKT5azkY6kYEbg44PPon3yXFFyfsg198+lCMGtQT/++cg/CbL4yNf8kAL//XKZ4fT9zJb+XAhKdS9Yo8HXVpot69A8D17vX1AG4C8OVQYsbuAnAXAIwbN64sWrjT4wDUHXro0D5YuHa78clbaeG+L43HwzOWeYNTREAHkGD2ymeNenGGDCbR2/99hjJcdRrgC7lMAC49fiS+/9d3sXf/7mjKZbww06V6AnPIepxyzlf1rYBYQDcRBU6Iv6fQh6jA67dXnxas3NKWylnMphZL5WK/wT3x3TP2x+fGDYtNe915h+DGJ+dhzF7J9SLl4M6LjsKrH6zHYMGQo93zeDfv69984Ugs3egHGSQifPXEfVOr57B+3fHDsw/C4N7NqekiVTBZ6ZYDEF3thgFYaZimSfcuY2wNv0lEdwN43LjWJaI9Rgfw808fiouPG+HZplcL+w3u6dmQyzhwj16YglUAkk3eEAfA3e0NiIiJlVGp0OkAPnvUMHz2qPDCIQaDSwPlLKiiDkCMmROFfDaj9LfQgRPqH507Bis3t+HMMb5S8KHLJsTG2BHBzSYnTTwQk/7+HoDS5NKmICJ867TRRmkPH94Xf/rqhIrVRYf+PZrwicP2CtzjHECcfm3kwB5YtcUxMZbzqASacpmQs17aMJkN0wGMJqKRRNQE4AIAk6U0kwFc7FoDTQCwhTG2KupdV0fA8SkAs8tsizG6awhASz6rlMPWEt88ZT/sN9gJQyFyAPdcPE57ghCg13PUWhfAly8TQgQAXz95FPYZ0B2nHTQ4FSpQzoZapQNI20uTiyF6t+Tx5RNGBji2CfsOwHH7hU0ydeDE5wIh0mjX0AR1LYhhz6Mw7bsnY85PJkamAaDk5LsqYjkAxlgnEV0J4EkAWQD3MsbmENHl7vM7AUwFcDaARQBaAVwa9a6b9Y1EdAScNWEJgK+n2K5IpKFMqhayGcJXThiJq//+HkYN9uMRnT4mbC4mQuZyuMzznMP0Z/NWBQnDbew3uCde+M9TAvfKUWSqOID/u+gojBoUr/NQ6QDSNtDwPEor6f5ZZzhsWJ/AgUdpw9QPwHTMPnbl8VgmiIe6MoyE3YyxqXAWefHencJvBuAK03fd+xclqmkKOGG/gXh50fpUlEnVxAVHD8dxowZ49tImkEUsLfkspl9zeuBc4lLx6BXHl7zz1YmAqgVVuR8/2OwoPv4mY6ws5XQUdCGsy8XkK4/H1PdWV0UJnDYmX3lCRfNv7+Tmq+novQb2bA6ZBndVNFQoiHsuGYdtrtt5PYGIEi3+OsgWF0nx+H+cgHw2gwMiIqnGgS9rtRJFlacEVugAUt6pH7xXb7ywYB36d09XD3PYsL44bFjfVPPcXcC9d5tKCHle72ioUBAt+WzZi2BXxykHlH9Iuw6HDO1T1uIPBM09a4Fy6I6okD79oMHo1ZJL3ULju2fsj79cNgGHdzFd1O6M/5p4IA4b1gfHKM6B3t3RUBzA7oxLjx+BRWu3475Lx8cnriF0oSCM3o141qMpix2asMgAsHf/7li6sbUsKxgxGujg3i14zw3YlyZy2QyOKeFM6mrhzi+OxQ7FYSv1jIOG9K64mKmrwhKA3QQ/Pld9XvGTV52E+Wu2Vbk2enx+3DA8PGN5STvxC8bvjSmzVuHzR4fNRV+ZdKp39KYKj3zjOCwssx/4iVS1tqSqJSYeUmMjAotUQZWMM5E2xo0bx2bMmFHraliUgUKRoa2jkEpcl2qjraOAm59egO+cvj+6aQIKWlh0RRDRW4yxcfL9+puFFnWNbIbqcvEHHB2SzmHPwqIe0VBKYAsLCwsLH5YAWFhYWDQoLAGwsLCwaFBYAmBhYWHRoLAEwMLCwqJBYQmAhYWFRYPCEgALCwuLBoUlABYWFhYNirryBCaidQA+KvH1gQDWp1idWsK2pWtid2nL7tIOwLaFYx/GWChSZF0RgHJARDNUrtD1CNuWrondpS27SzsA25Y4WBGQhYWFRYPCEgALCwuLBkUjEYC7al2BFGHb0jWxu7Rld2kHYNsSiYbRAVhYWFhYBNFIHICFhYWFhQBLACwsLCwaFA1BAIhoIhHNJ6JFRDSp1vWJAxHdS0RriWi2cK8/ET1NRAvdv/2EZ1e7bZtPROkfVFsiiGg4ET1HRHOJaA4Rfdu9X49taSGiN4noXbctP3Hv111bAICIskT0DhE97l7XazuWENF7RDSTiGa49+q1LX2J6G9ENM+dM8dWvC2Msd36H4AsgA8A7AugCcC7AMbUul4xdT4JwFgAs4V7NwKY5P6eBOB/3N9j3DY1AxjptjVb6za4dRsCYKz7uxeABW5967EtBKCn+zsP4A0AE+qxLW79vgvgQQCP1+v4cuu3BMBA6V69tuV+AF91fzcB6FvptjQCBzAewCLG2GLGWDuAhwCcV+M6RYIx9iKAjdLt8+AMELh/zxfuP8QY28UY+xDAIjhtrjkYY6sYY2+7v7cBmAtgKOqzLYwxtt29zLv/GOqwLUQ0DMA5AO4RbtddOyJQd20hot5wNn6/AwDGWDtjbDMq3JZGIABDASwTrpe79+oNezDGVgHOwgpgsHu/LtpHRCMAHAln51yXbXHFJjMBrAXwNGOsXtvyKwA/AFAU7tVjOwCHCD9FRG8R0WXuvXpsy74A1gG4zxXN3UNEPVDhtjQCASDFvd3J9rXLt4+IegJ4BMBVjLGtUUkV97pMWxhjBcbYEQCGARhPRIdEJO+SbSGiTwBYyxh7y/QVxb2at0PA8YyxsQDOAnAFEZ0UkbYrtyUHR+x7B2PsSAA74Ih8dEilLY1AAJYDGC5cDwOwskZ1KQdriGgIALh/17r3u3T7iCgPZ/H/E2Ps7+7tumwLh8uaPw9gIuqvLccD+CQRLYEjDj2ViP6I+msHAIAxttL9uxbAP+CIQeqxLcsBLHe5SgD4GxyCUNG2NAIBmA5gNBGNJKImABcAmFzjOpWCyQAucX9fAuAx4f4FRNRMRCMBjAbwZg3qFwIRERyZ5lzG2M3Co3psyyAi6uv+7gbgdADzUGdtYYxdzRgbxhgbAWcuPMsY+yLqrB0AQEQ9iKgX/w3gTACzUYdtYYytBrCMiA5wb50G4H1Uui211nxXSbt+NhwLlA8AXFPr+hjU988AVgHogEPpvwJgAIBpABa6f/sL6a9x2zYfwFm1rr9QrxPgsKWzAMx0/51dp205DMA7bltmA/iRe7/u2iLU72PwrYDqrh1w5Obvuv/m8Lldj21x63YEgBnuGHsUQL9Kt8WGgrCwsLBoUDSCCMjCwsLCQgFLACwsLCwaFJYAWFhYWDQoLAGwsLCwaFBYAmBhYWHRoLAEwMLCwqJBYQmAhYWFRYPi/wMA7/cWQFTIQwAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plt.plot(res._hat_matrix_diag)"
]
},
{
"cell_type": "code",
"execution_count": 21,
"id": "6ed8f2fc",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"[<matplotlib.lines.Line2D at 0x1d9d83df430>]"
]
},
"execution_count": 21,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAA5C0lEQVR4nO2df5QU13Xnv7d7ZjDYAxojBIhhBiFZRGL8I8xYgthrWZHkSF7ZikGyfiWxE8vI5yi70XpzEv9YE5acZJPddVbZszory4rXThZhLMDRj2OdlWTrh+WIMTNEmEEyCCHNMID4OSAsMDPTdfePqlfz6nVVdXV3VXdV9/2cw2G6urr6ve6u77vvvvvuJWaGIAiCkD1y9W6AIAiCUBki4IIgCBlFBFwQBCGjiIALgiBkFBFwQRCEjNJSyzc7//zzedGiRbV8S0EQhMwzODh4jJnnmMdrKuCLFi3CwMBALd9SEAQh8xDRsN9xcaEIgiBkFBFwQRCEjCICLgiCkFFEwAVBEDKKCLggCEJGEQEXBEHIKCLgQlkMDo/h/mf3YnB4rN5NEYSmp6SAE9F3iOgIEQ35PPenRMREdH4yzRPSxODwGO58aCu++dRu3PnQVhFxQagzUSzw7wK43jxIRAsBXAdgJOY2CSll677jGJ+0YDEwMWlh677j9W6SIDQ1JQWcmV8AcMLnqf8B4M8ASEWIJmH54tloa8khT0BrSw7LF8+ud5MEoampaCs9EX0awAFm3kFEMTdJSCu93R1Yf9dybN13HMsXz0Zvd0e9myQITU3ZAk5EMwB8HcAnIp6/GsBqAOjq6ir37YSU0dvdIcItCCmhkiiUiwFcBGAHEb0JoBPAdiKa53cyMz/IzH3M3DdnTlEyLUEQBKFCyrbAmXkngAvUY0fE+5j5WIztEgRBEEoQJYxwA4CXACwholEi+kLyzRIEQRBKUdICZ+bbSzy/KLbWCIIgCJGRnZiCIAgZRQRcEAQho4iAC4IgZBQRcEEQhIwiAi4IgpBRRMAFQRAyigi4IAhCRhEBFwRByCgi4EJdkMo+glA9FaWTFYRqUJV9xicttLXksP6u5ZLhUBAqQCxwoeZIZR9BiAcRcKHmSGUfQYiHhnahDA6PSfWYFCKVfQQhHhpWwMXPmm6kso8gVE/DulDEzyoIQqPTsAIuflZBEBqdhnWhiJ9VEIRGp2EFHBA/qyAIjU3DulCaDdnZKAjNR0kLnIi+A+BGAEeYucc59t8AfArAOIDXAfwhM59MsJ0CgsMiJeJGEJqTKBb4dwFcbxx7GkAPM38AwB4AX425XYKBEulvPrUbdz601WNpS8SNIDQnJQWcmV8AcMI49hQzTzoPtwLoTKBtgkaYSEvEjSA0J3EsYv4RgI1BTxLRagCrAaCrqyuGt2tOlEhPTFpFIi0RN4LQnBAzlz6JaBGAJ5QPXDv+dQB9AFZyhAv19fXxwMBAhU0VJDWAIDQnRDTIzH3m8YotcCL6HOzFzWuiiLdQPRIWKQiCTkUCTkTXA/hzAFcx85l4myQIgiBEoeQiJhFtAPASgCVENEpEXwDwvwC0A3iaiF4mogcSbqcgCIJgUNICZ+bbfQ7/QwJtEQRBEMpAdmLWCNkpKQhC3DR0LpS0IDslBUFIArHAa0Aj7pSUGYUg1B+xwGtA2CacLCIzCkFIByLgNaDRdkr6zSiy3idByCIi4DWikTbhNNqMQhCyigi4UDaNNqMQhKwiAi5URCPNKAQhq0gUiiAIQkYRARcEQcgoIuCCIAgZRQRcEAQho4iAC4IgZBQRcEEQhIwiAi4kjuRNEYRkaOo4cKkxmTySN0UQkqNpBVyEpTZI3hRBSI6mdaE0YorXNKLypuQJkjdFEGKmaS1wScgUTlzuJcmbIgjJQcwcfgLRdwDcCOAIM/c4x94LYCOARQDeBPBZZi65QtXX18cDAwNVNjk+0u4Dr1f7xL0kCOmCiAaZuc88HsWF8l0A1xvHvgLgx8z8PgA/dh5njt7uDtxz9SWpFCclot98ajfufGhrTSM4xL0kCNmgpIAz8wsAThiHbwLwPefv7wH43XibJdRTRMVvLQjZoFIf+FxmPgQAzHyIiC4IOpGIVgNYDQBdXV0Vvl3zUWsfvemuEb+1IKSfkj5wACCiRQCe0HzgJ5n5PO35MWYueZenzQeedmrlAxeftyCkmyAfeKUW+GEimu9Y3/MBHKmueemnHguKtSqaILHagpBNKhXwxwB8DsDfOP8/GluLUkijW6gSUikI2aSkgBPRBgAfB3A+EY0C+AvYwv0DIvoCgBEAtyTZyHrT6Baq6fMGgPuf3Sv+b0FIOSUFnJlvD3jqmpjbklqWL56NlhxhosDI56ghLVTlrmn02YYgNBJNu5W+bIi8/zcoEgMuCNlBBDwCW/cdx2TBAgMoFBpb1CQGXBCyQ9PmQimHZlrkkxhwQcgOkeLA4yLLceBpz5siCELjEncceNNRq5hsQRCEqIgPXBAEIaOIgAuCIGSUzAp4MxbKbcY+C4IQTCZ94M242aQZ+ywIQjiZtMCbcbNJ2vssswNBqD2ZtMCbKS5bkeY+y+xAEOpDJgW83ptNyo0JjyOGPM19bvRkX4KQVjIp4ED94rLLtTbjtE7j6nMlA1BYH9I8OxCERiazAl4vyrU202adVjKglOpDvWcHgtCsiICXSbnWZtqsU1OMt2wfLSm8UfogO1UFofZILpQKqIcPPC6UBT4xaSGfzwHMmLS4pDWepj4IQrMhuVBiJMvWpu7uOHDyLL7/85FI7p1K+yzCLwjJIQKeMGkMsdOr72zZPpqYeyeNfReERkIEPGHStoipk/TiY5r7LgiNQFUCTkT/AcBdABjATgB/yMy/jqNhjULaFjFNknQHpb3vgpB1Kl7EJKIFAF4EcDkznyWiHwD4ETN/N+g1jbKIWS7N7Adu5r4LQlwktYjZAmA6EU0AmAHgYJXXa0iyvOhZLc3cd0FImoqTWTHzAQD/HcAIgEMATjHzU+Z5RLSaiAaIaODo0aOVt1QQBEHwULGAE1EHgJsAXATgQgDvJqLfM89j5geZuY+Z++bMmVN5SwVBEAQP1aSTvRbAG8x8lJknAGwB8FvxNCu7pCmtapraIghC/FTjAx8BsJyIZgA4C+AaAM23QqmRprhnsy1rblyKsTPjspgoCA1ExQLOzP1EtAnAdgCTAP4VwINxNSyLpCnuWW/L+ISFNY8OweLSW+aFdCNRPYJOVVEozPwXAP4iprZknjTFPettISJYzKkYWITKSdMMT0gHshMzRtKUVlVvS8eMNqx7YldsA4tYgfUhTTM8IR2IgMdMmuKe9bYsmdcei+j6WYEARNBrQJpmeEI6EAFPgDRaqNUMLHp//PKJb94+KtP6GtDb3YE1Ny7Fk0OHcEPPfPmcBRHwuGk0P6VfNItuBTLQcNP6NA7AgN2udU/swvikhW1vnsCSee2pap9Qe0TAY6bR/JRmf8bOjHv8/AASTUlba9I8ADfab0uoHhFwjTgsr0bzU/r1x3THpGXhNg7SLJKN9tsSqkcE3CEuyytNkShxENYffcC75+pL6tjK+EizSDbab0uoHhFwhzgtrzRFooQRdcbh1580uxqqIe0imZXfllAbRMAd0mx5JUG1ApxmV0O1iEgKWUEE3CHtlldUolrV1Qpwsw14gpBGRMA14ra8ah2OVo5VXa0AN8qAJwhZRgQ8IQaHx3D7t7e6Arnhi8n7iEtZ1eaAYgpwuQOOuBriIa1x50L6EQFPiC3O7kTA3uiyZfto1VvYS93kYVZ1kHWurtWoi5JpRz53oRpEwMugHEvJLBVdWenoqfeNcpOHuTVKWeeNvCiZZqJ+7mKlC36IgDuUukHKtZRWLevEpoH9mCgwWvOEVcs6I7+XiZnb+75n9uDeay8NFHG/46V83rIoWR+ifO5ipQtBNIWABwmmOq7SrYbdIOVaqL3dHdiwekXR+1ZyM6qbfHzCggXgZ3uPYdubJ8q6kUstOkbdsCPCES9RFoNldiQE0fACHiSYD/ePuFVqckQoWAxG8A1SiYXqZw1XcjOqm/y+Z/bgZ3uPVXwjl1p0jGPDTlxi30yDRqnvRWZHQhANL+B+ggkAax4dwqRle6aZGfkcgZkDb5C4wubKuRlNEbv32kux7c0TNb2Ryxlw4prqi8vAi4RsCkE0vID7CebWfcdRsKaWFfM5wrqbekoW/Y0jbC7qzRgkYrW+kcsZcOKa6ovLoBgJ2RT8qErAieg8AA8B6IEdaPFHzPxSDO2KjSDRm9Zq+5RzjnjfcWVXnVvqJUjEan0jlzNoBIl9ue6QamYpaSPt7ROyDTFXHuBGRN8D8FNmfoiI2gDMYOaTQef39fXxwMBAxe8XJ/W6sfwKJPhZ/uo8JWJZcSOYn2ul7pAo30/aXS1pb5+QHYhokJn7zOMVW+BENBPAxwB8HgCYeRzAeKXXSxI/MajXlNQMCVQLqeYNnlW/p/m5VuoOifL9pN3Vkvb2CdmnGhfKYgBHAfwfIvoggEEAf8LM7+gnEdFqAKsBoKur9m6KKFZQkta4ee2OGW3IEQHMyOUIFnPgDa7+VguvWbz5k4ygSHt0RtrbJ2Sfil0oRNQHYCuAjzBzPxH9PYC3mfkbQa+phwvl/mf34ptP7YbFQJ6AL39iiaf4QJLTXD93iYo3zxHhro9ehO++9Gagm6RRpuC1HCDTRtrbJ2SD2F0oAEYBjDJzv/N4E4CvVHG9RChlBSU5zTWv/eTQIfcxgdE+vTXUTdIoU/Ak3VVpj85Ie/uEbFOxgDPzW0S0n4iWMPNuANcAeCW+psVDKV9yLaf4N/TML4rjDnOTxNk2sQQFofGoNgrlQ7DDCNsA7APwh8w8FnR+mqJQdGo5xS8nSmNweAxbto+CYedWqbRtca4DyEAgCLUnCRcKmPllAEUXzRpxTHODhM28dtQoDVN09WRY5RIlT7h6r5Yc4Za+hVjpM2A0ik9eEBqFXL0bkEUGh8dw/7N7MTg85j6+86Gt+OZTu3HnQ1vd41FQbpI8wXWTDA6P4b5n9kyFG07aGQjLuW6p99DxhDYWGA/3j/j2IygtgSAI9aHht9LHjZ8VWslio26x6z56ALjzoa04N2GBARAAi4EXXys/A6Ei6jqAes+gpF4SFicI6SJzAl5vH6yfWJcrbGEui/uf3YvxSVtIcwC6Zs/A8PEzRaIadaeifk5Yjpf1dy3H5u2j2DQ4ikLBvx9Jbi7SU/uWykkjCIJNpgS83j7YweExHDh5Fi35nEfkyhU202Wxvn8Em7ePYv1dy4sGg9UfuxjrntjlGRyiLkqW81kpgV+1rDO0H0mExelttUMs7Vw14mMXhHAyJeDlpjaNKqjl5t1oyRFuu6ILK52Fxfuf3Yvli2d7NgiFsXzxbLTkCOOFqQig8Qm7P/dcfYlrDROAJfPaiwYHZaWHfQ5JbmGPG72tQLALRxAEL5kQcH16HcVVUY71GfVcXWQKFuPC86YDmPJX58vIatjb3YFb+hZiff+IeyyXI09/VFFkZZnrg0OUrH9ZyujnVhxyPt8c/BdbhfpR79+I4E/qBTxq9j6dcqzPsHNLCeLWfcfdhb9Ji7Hm0SEsmdce6Qe+clknNjsinSNb/PVNPWHt93PZ+A1E1eQdryV6f8QHXjlJiWwafiOCP6kXcFPMxs6Ml3RVRLU+g3za6rkogpjPkVvZx2Iuy1URJLBR2h8lnvyeqy8pO6Pf5u2jiVpaUePlhfJIUmQbJaVDI5J6Aa+0FmUp61P94M9NWMgRcM1lc3H3VRf7WsHjkxbWPb4LSxfM8uyI7O3uwLqbejwpYf3aV65oVRLtUWmIn/66fI6waXAUk4XKRKCUBSiWXHIkKbISPppeUi/gQWJWSixKWXS6+6PAwE9+eQR3X3Wx+7zpl90xego7Rk9h08B+bFi9wr32HVd2Ycm89sC2VCpaQe0PGwwqCfHTX3fw5Fls+PlIRSIQpZ9iySVHkiKbZPioUB2pF3CgWMz8xAJA2RZrmPtD/Wjve2YPfvraMfd1EwXGusd3oWfBLDd2O2ywiNNFUUokK3VDqNcNDo9h8/bRikQgijiLJZccSYusuLjSSSYE3MQUiy3bR90FwahWbhT3R2+3XQm+/40TGJ+0ANghbsoaf2RwFBu+6B+D/cDzr+PI27/GioguiigLUGF5U8ISZkWllAiEXTeq314sueQQkW0+MingplgwUNHUvJT7A7BvirWfWopv/PNOaGHbAIITQ936rX+Bo/fYMXoKX/rYYrRPbw10UZTKSBgWCRNWNKKSzINhrptSs54o4lwLkZGQN6FZyKSAm2FnQwdP+UaSRL1WqZt87Mw4/JLuBiWGUuKt2HXobfzTF64sclF0zGjD/c/uxcGTZyNlJPSLhDE39ehFI0plHizHJ+/nCtriM+upt2DKQqnQTGRSwIGpwgd+uyPjvmHNSI2PL7kAc9qn+b6XvcsSHhG/oWe+2+b1dy3Hlu2jOHr6HNY+NoRJi9GSI98BKEpoYJSiETpR3TBhn0FrSw6E4lmPun45+VniQL+mLJQKzURmBRzw3x0Z5WatRERWLusEOf+Hvaa3uwMb7/4t1wd+64e7inZnbt4+6kbAAMBkgXHblQux4LzpnjaZoqksdv0cP9dFmFuolBsmKB+4+T6qH3rbSrlYkrCO/VxIslAqNAuZFvBSC2d+Ql2uiJjnr4xQWKG3uwPf/gP/Ohdq0NFdMhaAngtnFQm96Sry82379THMleEn+LobRuUDV1v4w6Jc9OtEcbHEaR2rfh8w3E9jZ8ZloVRoGjIt4GELZ0FCrYvIuQlbaMIiQaoVHfOafomscmT72YP66Ofr3rrvOHa/dRprHh1CweKi7H3lLFSqgTAoH3hQqlfzOqVcLHGFEXpmDPkcWnKEgsWe7JAi3EIzkGkBB4KtzSDh1QWUATwysB+AbQErC1flJrnjyq6qRCcocsNcEG3J+yekCvN1d8xo80TGqGyGlbgq1EDolw88aqrXIBfL+KQFIkLHjLbYwgg9rrOChWsum4uzEwXc0DM/88ItETRCOVQt4ESUBzAA4AAz31h9k+KhY0YbckQA2CO8vd0d+PiSC/DUK4cB2BtzNvSPIO9YcQx7U4+emKpS0fEbRADb561zc29nSeE127F5+6gnrJEIoYufUeLie7uL84Hrlj8QnurVHEzX3LjUnSGse2KX+3lWK0yeReV8Ds/tPoJJi7HtzRORk4mlkTRE0MgAki3iqIn5JwBejeE6sTE4PIZ1T+xCwWLkiLDmxqUe18Jzu494zmcAlsUgmjpmWewKrnltvR5mGEpo9FqUyxfPRkt+6o3a8oRVyzqL6mD61Zzs7e5wo1DIeK9rLptbtPgZVAMzDP099GvlnDdUqV7VgmrY5zB2ZhwWs0f040ANZl/+xBLc3NuJSYsbok5n0IBfK6qp7SrUh6oscCLqBPBvAfwVgC/H0qIY0BcKmRljZ8Y9i15q+zwBrjC1tuTw+RWL8NCLb8CyGG2t/hEaIIqc7Km3uwNrblyKJ4cOudP7weExkPbeaz/dA8BbBzMXQXhXLuvEI4NTUSB6Hpc4dzzqfVg6fybap7cGLqiaJJ2fQ32eWyrc/p826p1qQEIws0e1LpT7APwZgPagE4hoNYDVANDVVbrYQRz4+Ytvf/AlTBQY+Tx5Fr3M/OLXLZ0XvFGmYC/xRa0Yo2YC45OWO723N/qw6wcfOzPuGXByAD5yyfm499pLAQBf++FO3/DF3u4ObPhisEjHtZBn9qGcaBK/gSTuKXqcg1W9qXdf6j2ACOVTsYAT0Y0AjjDzIBF9POg8Zn4QwIMA0NfX57ehMXb8/MUq6mOywLjkgvfgiove60kNq1wjyxfPDqx+k3cscL8dn37CtEWL9y4ViaEyH+aI3I0/tz/4kttuv7wruhVqxodXUvTYDz+xLudG1weSpHy8jRR1Us++1HsAEcqnGgv8IwA+TUSfBPAuADOJ6P8y8+/F07Tq0G+ELdtHPc/tPfIrjI6dwSonpjvK4qGqUdk+rQW7Dr3tiXgIijZ5ZGC/a2nn894CyOp66j3MBb9Vyzodi98myNJ9uH/ETcilNuEsdSJqwkq9RRVTP7Gu9EZvtil6FhcEG2kwbAYqFnBm/iqArwKAY4H/aVrE20T5i1VGQcC7SKQvHp6bsLMbmj9itSlFhdLpEQ+B0Saar11FmujXOzdhYeO2/Vh3U0/Rgh8DaM1PxYsHbVRa8+iQ+z5qE46eJjeo1Jve5vEy3SDquBl3X0qskpqiRxXKWgpqGiJKhMYn83HgUVD+YjPOWW3/1re1q9hw3eesi506x4wtD3KLqGOrtB2cW/cV19Jcd1NP0fmrnLqZx06fw/nt09zX6wuyBcvrlVJhkDmC216/Um8dM9q05+3H5vWDNuyYlCNWUVMSRCXqe9daUJtttiHUh1gEnJmfA/BcHNdKCr84Z79t7YDtJ9dvOCXQZtV0PTeJn5UaVvPSLCYRtgVcCc+W7aOedLEt+Rxa84TJgh0CmcsRLGdx1o2oCch1PnZmHHaUvN0ftRO0EqGLIlaVpCSIQlShrERQq7HYZUGwcrLoeqoXTWGBB6HfZLkcuT5nC16LVHcjqO3kfqF0+uJnqa3sfsUk/CxdU3j0dLGFgoXbrujChU4SLHW+GVHTMaPNE1Osjk1rLRaYqK6VoM8xSKyiCmi5N29UoSxXUKu12Bt9QTApkRXXU3k0lYCH5dc+cPIsvu8UW/DLTWKKq19uknKiLaIUkxgcHsPBk2c9YY9muli/8ELzbz1vCJgxabGbuU8PoQTCXStBRBGroCyIUbIVlhoMoxaSKEdQ43CBNOqCYJIiK66n8mgqAdd9zyp3iLKaH3j+dRARcppFXGkJMb/MfH7XCVsIBKaEFwAWn/9u/NFHF0cSfrPPeluAKR++Em9lnfd2dwS6VkpRSqxMAdX7F5at0O881a+oPvqobdQRF0gwSYqsfO7l0VQC3jGjzfV3KzfJ4PCYJ946nyN8fsUid8EzaNdlmEVnxo4/MrAfkwXbgvaroQkUWzWrlnV6Fk73Hn0Hax8vnU/EL/uhnjcEzK4175fDe/ni2b6ulTjQ2+03g/G7eU2xCKp/KhuEakfSO2zlc49Owwl40I08ODyGJ4cOTVmXNLULUo+3tizGQy++4Sa2Amxf8H3P7MG9114aakHrx9WPcMf+k27iLLUQqQuO8qmbea0ZABEBHB4Lrsq0EeDGf5viZlq+5iKu7u++5+pLanIDlRNfHlT/dDzESo9LxEVAiklaZOVzj04mBTxMpIN8qOq4cg3okRl6vHVOy0oI2DHcFgMvvnYM29484es39kP9CL/2w52e46y189cTlvserS3evNarlnXi2OlzrvgD3oyDqr/m7MGyirf6mzeE+nv3W6d9/d3lWPiVEiQCfm3Vz/Nrs/hNa4+IbDrInICHLaAE3cj68RxN5RpR4n5L30IcOX0OF7RPm7JiJywQARfNeQ9eP/Ir12+uR45Eqfi+alknNg3sx0SB0epkHlS+eAUDmCxYuF2LKOnt7sDdV12Mn+w+4qafJbJzbANw+2XOHlQoIeUIB0+exeDwWOCNVq6/O+7Fq0p811v3HS9qc5r9phISJyRJ5gR8s09+ET2CIue4HfQb2bzB9URRfn7ukePv4MGf7oPFwMiJM2jN25YxEcFiDrT0dIHTi0JsWL3CE87XMaPNEwcOADki34iSW/sWYn3/CACgYAEb+kewxSl3tnzxbM/soSU/ZYFPFhgbfu5fGk0R1d+tROig4eZR7otaCpRfm9PqN611SJwMFs1HpgR8cHgMmwZHp/KL5MgVHZU1z2LbCtVzgPv5gW//9lSEB+AVJHsDzNTx6y6fiw8uPM+N/Y5S8d0sCgF4/bR3ffQi19eu8pX43XRLL5zleawPXPdcfQk2rF7h+sCPnD6HpzWXS1CIo36TlxI+M51uSz7n7mQ9fXYCt37rJd+SbtVQSdhgGqf0tXTtSPx0c5IpAd+67zgmC1N+41v6Fvq6TwiMXQdPeTL06Tf413+40yPeygetFvb07ekM4LndR3D3VRejt7t0xfecY6UDU0UhTDfOxKSF9umt2Hj3ipIW066DpzyP9bYC8PT/6OlzoecG3eRhN7re7oLFuPWKhVhw3nR0zGjz5mGZKC1QUTMklhKioDanzQINqgqVBLIO0JxkSsBNV4i+HdsMl3tkYL+7YcUUAXPr/Ac6Z2HNp6Ys9mmtOXeBEbDdEese34WlC2Zh1bJOz45Lnd5ubYelVhTCr+1BOy91BofHMHTAK+BmW00LWW2tz+cJt/YtDMzpopJ2qeNhG3BanF2q+Ry5KXjvf3avO1ABxQusfn2JYiH6xeqbbqpyFrDrhV4VKm/MCJMgzesAQnJkSsDDpvz6c/quSj9rxFxY1AVRXceNA5+0YAHYMXoKO0ZPYdPAfmxYvQKAv/AFbbQpx0+rQgOVf17RZrRVlWFTgmdbyF1YoC2E6ixfPBst+ZwbjbNxYH/oQOfiWJF6zTkl7Hr0ThhRLUQzVn/H/pPuQqyKulHf24bVK3wHp6Cwz7goJ486Y6oqVJLUch0gbTOdZiYzAq7/aPwsYHMXY1iZrd7uDndhMcjPqiJI7ntmD3762jH3uYkCR9pMYu5w1K9bqp+mf15FztzQM9+T08QMjcznCAQE3li93R24ubcTG/pHbMEvMAoozq6o90O5rezzvaGJt/QtxMPOtXR3kR+mJR9kIY6dGfdkUnz6lcN44bWj7qCqBozxAmOzlvZXWaBKxFXYZ9yWeClLX4/vr7VFXM06QDkpedM002l2MiHgUW6aoBwnYWIW5Ye38L0z0JqfSnTVmifPZhK/Ld/l1s7U+dbzr3vEGwBacnaVHn2Tjr5TM0fA+xfMwiuH3sb6/hF8f9t+/KVPEQfAnn2owS2v5UYhIndnqv5Zfn7FokA/7kon3W1kkfKx5E2UECsXlj64mK/SHysL9L5n9uDF145FLntXLmEzCfOzC9szkCYrthxRFl97usiEgJf60fg9r1dWL4XfzWT6lj9x+Vyc3z7NzeutXDDKmvRM4Z3amUBwgQg/Hu4fwdOvHi46fkvfQoydGS/aqalbeHNnvgs7Rm1/ecFifMOniAPgvzFGrwS0UhsYxicsT6TM51cs8swq1LW2bB8tWlcwCbLkTXq77epE33DaBHijjfRCzmZK2t7uDtx77aWeZF9+g0pSaWLNNYahg6fw1595v+/7p8mKLUeUxdeeLjIh4KV+NNX8qIJuJjP64oMLz8M9V1/i+qddwXKsyY4ZbUVb3wH/AhFB7Vjz6JD5ck/ubL+CD0qI1AYfRSHEpaH7jQ+cPOupBHTs9Dk3lj6XI/c5lWLAbxOTciep+HQArqirRc9yvqOxM+Ng54Mwo43CCjmrvoXNvpJME6uvCzCATYOjnrqriris2Lis+HK+m1r62oXSZELA9YVFv8l3NT+qoKgHP5+tuvn1Cj6FwlSCJbM6jsIUU78bzwxfJAKuu2yuG74I+BeI0Pv6g20jUN4XPVWA+X6eOpr5qS38+RzhuT1H3Vj6uz56Eb770psYn7AAgptiIGiXq8q8uEkrX6cWfcv5jsKijcJcX6XWSdTnnFSaWHNdIGimEYcVG6cVX+79U42vXYiXTAi4QtWl9NtdWOmPyi9DoYvhs9UjC4CpOGvlEzfJ56hoV2jQjad2GI5P2MUl/AoRh/Wxt7sDG+/+LXeQUxa/n19Wd09MTFr4QOcsLF0wCwRggxO9Q2C0T2/1FFtWScDCdrkSptLWAvair77wGeU7qmRArqZIc5xEWReIw4qN2xctopxNMiPgSS2e6FEPeiEHP5+tJ9bcqQCvrEO96DFgW8BrP1W8iGX6SVUkRRw3td9NaH5uG7eNFG1U2nngFHYfPo01Ny71Teeq3Cg5ePPIqPc0d7nqBaRb88ERJ+X2JYwt20fdhc+wSkJJuwCiXr9awRRftABUIeBEtBDAPwKYB9t4fZCZ/z6uhpkk9YMNum7QxhuztBrgX3JNXUcP+1PXjeInjQuzHxfMfBeAU55zlLgH1eU088j4LYzqxzZ8cXmRDzxJBofHsHHbiKc/YZWEkrY2a2HNii9aAKqzwCcB/Edm3k5E7QAGiehpZn4lprZ5UNEJTw4dwg0982P7wQbdCGHHAf8qMTq73zrtm5tb+UlVgqrJQrw5rf0qypsW8vO7jzi+fSCXm8pt4rc7tBKhSFrAzD7a6wfec5LeOJMGxO0hVCzgzHwIwCHn79NE9CqABQASEXC1NXl80sK2N0/4hsiFvTZKTPjg8Fhg/hRPOx7f5cmIqBYx1bEc2dkFld9Yd5UMDo95cpbEmdM6LNeJfq21n+5xB8Io5dnKEYpqIiMGh8dKWu5+fVy+2JuVsa1GLoU0xXILzUksPnAiWgTgNwH0+zy3GsBqAOjqKt5YEpVSGyiqCRtTwlFqW7lZQAGwFyqVSLuLoQwAjBwBBYbrKulxco3rucBVTms98VE+R9ix/yS+/sOdgeGHQZEspSrKmwPh+ruWB0ZslIO5/b/cWYS5A1WPXtEJivnXszKWCtmMg6i/KxF4IUmqFnAieg+AzQDuZea3zeeZ+UEADwJAX19fqf0egQT5qtWNr46bNSdLWbZ+oYFBuTS27vMWUADsGOWlF87ybGQh2Fbgx943x62oM1mw8OTQoaIolrbWnJumtmDZoj9psfu6Hwzsx/cNIQsSD7Oi/OmzE0WfYxKLwX6foZqZRBWwrfuOe6JXxgv+cexBv4NauxOi/q5KbbkXcReqoSoBJ6JW2OK9npm3xNMkf4J8sSr6A/DWnFSUWvz0Cw1UuTReev24J5zPb6q+0tlMoyJZCMBH32dHaux+67QrxBYDs989JbAAcO3lc/Glqy72tMEYH9zcK1EGJb3CDmDnNb9u6bySn0e1YuL3GZbKCGmyfPFs5J1Migq/hci0LN6pxGAqJUHQ7yrKlvt678QUsks1USgE4B8AvMrMfxdfk4Lxs7LMHNi6/ilhCstJ4QkNzOfQ1TEde4++Y1e1sbxFGXq7O4qm6gBw8ORZTz3LpfNn4r5n9mB6a95T/uv4O94SZh9aeJ77+lyOYJnq7dMns81m5Ixe6ccvwZTfomYUMQkTeb/wSgZCM0Ka9HZ7qw/pIZ1+58YteGFpagMHC7Vt1tw+i+hb7iWfiFAN1VjgHwHw+wB2EtHLzrGvMfOPqm5VRAaHx/Dc7iPu4xan5qR6Lqrv+9+8b45bD3PtY0Oec0wR1MXDky8ln8OtVyzEmXOTeOCFfe7r8znAsoCWlhxu6Jnv5ulQfu77ntmDyQIH5ndqycHtkyIsQiYoH7n5erWgqqejVTspw/LC+H2Wfu1Rn205YZ9lJ8eKiSjFsP3qr046i9R+aQvCZgoSwy3ERTVRKC8Cvjvba4a6ieA05LN9CwEA9z+711O/0S+h1ODwGG5zcksDtlAuvXBWkY/bFEEz1aqbL6VggQA8tuOg5/WWcusyY8m8dk+u8adfOexa1z5GHC6Z82787c0fLGshcOzMONbd1BM449DP1f3WBHtBVl+EVDOXA9pnOT5hrw/c0DPfXXxV76UvhlYafhjlNXH7j4Ms4jBLOYoIB80U0uIGErJPZnZi+mHeRD0XzvJkEMw7bgm/hFJbto96xHrSgpvNLwfbYr65t9MTzua3LV1/fwY8Pm5gyv2hrLR7rr7EXbALW9FtyVOgeJuZEtVCql/ceRBFfmsCPr7kAjzz6mFXqM18KZMFhgV7feCnrx1z3UEExFYT05zhlDsbKIdSubvDRFoXYRUGqo6X209BqJRMCrh+Y+uWjMcithg9C2bhF6OnfKe5QeJJAN7fOQs9C2a5Pm4VG26G6Zk7FwF7YFBhgpfNb8feI79yfeNqwfDl/SdLpl/9rJaBz8TTjgLj4f4R5HPkm2wqiOWLvfU7wcCc9mmuYJHznJpdXHPZXOw5fBrDx89MzRq0/+NeqIuSJbLaTH6lcneXspTVY1mQFOpF5gTc78bWp+76Ytrcme9Ca/5tj4Aq9LJqgFeMdo6ews4Dp/DIwH63MEOOCDd+YL4nTK9jRluRJaW7SH751mnXN95z4Sw31tx00wDwVKFpM3z5poComHElvmqAUoRVvFH4+ctXLut0o2pUaKNa3H3O2b1phkuqBdmoC3VRCxubvnn1+rj8x2b7xs6M+8bDl7KUZUFSqCeZEXB14+u+bfOGURaTEtBnXj3sCqi5s09FlKhrqhSggJ3YBewtzGAxe/zbagOO3jZ9a/dkwesbVxt4wix/vWiEvoh2bsIekNbd1AMAbnbAfI5ABFiOa0Oh/x0mmGH1OwG4z6kao3rbW5yiyUsvnOXrbw8KV4yysGyWijPjvYOsYrOsXthAEddAIAuSQj3JhICbPt+WvDd/h44poOOTdpEC/SZ+uH8ETw4dwtL5M9E+vRVLL5yFaa25UIEF7IXGlhzBsux82R0z2tzc2gWLXT+weVMzUBQnrXZp6tdWRSMUW/dN5SqftBjf+OedgLNF334N47YPd2HkxBlP3c6CswkGgLvJqSVv+8r9BrIgi1GPVtmipQoA7EHjwvOm+5ZtU6xc1unZGXn/s3tLWqu6RatqgZobqvzabEYEgdnN5R6UmjeOhURZkBTqSSYE3PRt33rFwsDK64Dj39Xiqp965TAe7h/BHVd24eH+EXzthzsBwBW9d7XaPtBdB0+5Lg5dbFVa8LZWu0akqkyz9rEhFHjKfaEKQtxz9SW+vnHljrjq0jl4fs9RWJqoqwHB7IcnrpvhVqoB7Hwryk/f/8aJohSu+ianiQJjQ/+IWzWnHKFRIvU3T76KbW+O2W2B/0Ybv7QEqo1RrFXznKiV5U1XhvqUzFh+s19xCK4sSAr1IhMCbt7UpVKU9nZ3YOn8mW5UCQA8OXQId1zZhSeHDhWdr3ygf/WZ97uxyJsGR10rX1/gUpV3GCjyCec037O6qR/uH8HGbSO4aPa70daSw60f7sLYmXH8+NXDbgQHke2iWffELgDwuCRcPzXzVLHkyamiD+pz8EvhapZZi7rAGcT2kZOex0MHT3keB22pV+8XxVqt1KI1N2QVCpY7w7E4uLycIGSZTAi47tv2Czw3/byDw2OYa+S9vqFnvvu/7m4AvIt+Smj0epPqxrfFeL8rTiodq1rk1AVVna+sfcWrh4aw9tM9nmgPNSCMT1j4hrOoqPK6mH5qoNi3q/pvJnHSF2oZxdV0SqFnBySgqGSc+V34bak389ZEEeZKLFpT+FWxZlXDM6zPkpdEyCqZEHCFX0k1v3AwFQ/dmidcPn8mbv1wl+sDvePKLjy3+4ibo4Rgxz+bcbymiKiiw7qI3frhLqxc1umK3JJ57R4x8LX2C+wJP9SjPfRokvFJC996/nU8+Ad9nt1/5oYZv4VBdW7HjDbc3LcQBHgWGwF40ubq19IHCj07YN5Q63wORVXhTStYj6P3Szqm9ylqhEoY+nfW290RKVVuLfKSyAAhJEVmBDzqbjmV8c9igCzGJ5bOK1rAuvuqi/HCa0fdcMPn9hx1I1ZM0QmzQJdeOAvAVFX2TVrYYVuL7S83rX3ln9bFZsm8dmzePooN/SOec5965TA++8C/4LwZbXhu9xHfxFD6Que5CQvrHt+FV9867Q4IObIXXtE3JbhRto2vWtbpyQ5oRj7+9m/M9fUpm75/NTCaScf+fNMOjIyd9ez6VAOv2pxUbVrYKJZ80mGAjZa4SgajdJEZAQ9L4KQfXzp/Jl56/TgARj6fw4GTZzE4PFYUxaCE5uDJs24h3/FJy13oW3PjUqx9fJcrOi0528etRJxg+6o9AuCEHSpfc/v0Vvz1Z96PjdtGMK0lh0vmtnsGB2Ulj50Zx7HT53wjYH7uLBoqTJHRizIz4PH7A3A3+6zvH8EjA/txc99CX8HSC1Io8W9tybn9z8EbnnhB+zTf70mPXNGF62Pvm+M5b+/Rdzx90gdetTnJr3h13CQdBthIceKNNhg1ApkR8KDFLf24mVfbYsb3f+4feaELjS5eSsCeHDrktUAt4LrLL8BPfnmkKFGU68/OEdiJFNFLlN1xZZcr2IDtG/+G4Y7x8+2bqHN27D/pDkp6UeZSjBcYx06fK0qDOjg8hk2Do5pv395ItEpzD82c1uJJ0qVmH0GYwnV++zS0aal49T61aom+zO8hacFLOgywkeLEG2kwahQyI+BAeHIgPc7YzUlSmLKGg4oL9HbbtTY3bhvBroOnwDwlKGZo3t1XXYy7nfzd+nXW37V8apels8Hm8ysWuaXWerQ8JS1OWGBQzpQg8jkCM6PAtmvlud1HsGH1Clcg1I2lQwDmzpyGt972ptw106CquHn1mlu0bfzq//uf3evZeVmq5qRf5JAaEB4Z2O9uRNJdJcqVpEcA1ULw4goD9HMvNFKceCMNRo1CpgS8FGZeahChULAtzY3bRlCwbCHWS3UNDo9h7eO73M0uv/0bF+D89mlYMq89sLq6n+9XTy9qWYxvv/iGx93ihtUZoYd+qI0+qkDEtc4OzYc1H/lEYSo5lj4DGXJi2ScL9izg319zqdu/1pYczm+f5rZzosDYvH0Uq5Z1em5Mc3FSfbbTWv0rIgUNjEEzppU+ET7qOT0CqJIkUaVIyocb5l5olDjxRhqMGgVivzymCdHX18cDAwOJvodyieiRF8/tPuJuQAHsSBQlEjv2n3QjUoApS7fcOofqBlahgZNRfBoB5Ai4/YouXHjedNdHfvrshMeF0ZInbFy9AgA8vvSOGW1Y+9gQJgrsDlb6ObsOnsLGbSNQ3qG2gIiQoM/WjBpJyieaxLWTuqa5lpIn4MufWBJLrVFBAAAiGmTmPvN4Q1ngChXxoKIbzA0ox06fc29k0/esh/GZPj5dAFTct4pwMX3xpo/bjwUd07F0/kz8ZPcRt5QYYapUGwDP++nMa5+GB55/3U0ypcg7lrsKSVRWun4toqk5QaEwtXu0lJiZlqTpE/UrBqE+t3KttiT8rXFfs5wUD/VAIkYan4YTcL+wQsvI1DenfZon34YSPT0/ico2CPgn0rKY8Z/+eSdGjr+D9umtRaGBADwinoNtNSvXBQAcffvXuPu233T96npxhN7uDnz9hzu1XY3sWawcPflrjJ78dVH/CzwVs60LiSfHCOxdnfpiayWYLiu/ivSVWr1J+Fvjvma5KR5qSa0iRmSQqC8NJ+DmTaqiG8Ynprafq8UydY7aKq9nJSTA3YijfMj5PNnWK0/lJnnghX3IEYpEa+zMOP7SqYyjC7OK9zYtZPPH/3D/iCcDYF4LYSzFNZfNxQcXnue5qczPJahOaDk3pFoAfnLoEKa35t1iEFGr2oS9VxL+1rivWW6Kh1pSi4gRCSuMTlIDXUMJeFARY78deX43svKfj09YsAD8bO8x/Mvrx13hnCwwrljUgW3DY54SaPpNApRO8F+qVqS769N5DwJw+fyZ2HngVNG5Oiok7+6rLq5IEMu9IQeHxzzRNX4uhCCrN8p7VbP4F7a4GtcNlOZFvVpEjEhYYTSSHOiqEnAiuh7A3wPIA3iImf8mllZVgN+Wej2CIYo4qBvyvmf24Gd7j9lRIMYi74l3xotSweo5Rkr9qMNueiU6B06enaqUA9v6vvXDXdh92NlY5EwRLLajam7pswtGlKqDGSRe+vuWc0NGcSEE9Teq/7wSwtYq4ibqgFBrV0OpwSWO9khYYTSSHOgqFnAiygO4H8B1AEYBbCOix5j5lVhaVib6h6TXcyx3xOvt7sC9116K/jdOuKGFzIyCo5v6DkICcN3lxe6KSordmvmsVby4LkBqJqEiHgA7ZHFBibzcYfi9r18FIz+iuhD8+hvFf14p+m/B4uB0srWiXq6GsAE7jvakeQaSJpIc6KqxwK8AsJeZ9wEAEX0fwE0A6iLg+oek13OseMTjqaiQ/3zT+/Hk0CG8+NoxT6a9aa3F7opKf9Qea7Zg4TYnjNC0ZpWPfXMJN0xUorxvENXcwPpr9RC8OCyU5Yu99T4tox5qrUmbqyHO9jRKjHuSJDnQVSPgCwDs1x6PArjSPImIVgNYDQBdXclMYwH/LfWVCpy+Kadg2dkD7732Umx784RrMYYlW6rkR22O0mGJnOL8QZTzvkFtqfbmj3NAUtc1633Wc3qfNldD2trTDCQ10FW8kYeIbgHwO8x8l/P49wFcwcz/Lug1tdjIo6jGx6dvymk1okuSnC7WKyQrDaFgSbQhDf1KY1vS2B4hnKCNPNUI+AoAa5n5d5zHXwUAZv4vQa+ppYBXi/zABUFIC0nsxNwG4H1EdBGAAwBuA3BHFddLFeLbEwQh7VQs4Mw8SUR/DOD/wQ4j/A4z74qtZYIgCEIoVcWBM/OPAPwoprYIgiAIZZCrdwMEQRCEyhABFwRByCgi4IIgCBlFBFwQBCGj1LQiDxEdBTBcszec4nwAx+rwvrWkGfoINEc/m6GPQHP0M64+djPzHPNgTQW8XhDRgF8QfCPRDH0EmqOfzdBHoDn6mXQfxYUiCIKQUUTABUEQMkqzCPiD9W5ADWiGPgLN0c9m6CPQHP1MtI9N4QMXBEFoRJrFAhcEQWg4RMAFQRAySuYEnIiuJ6LdRLSXiL7i8zwR0f90nv8FES0r9VoiWktEB4joZeffJ2vVnyCq7Od3iOgIEQ0Zr3kvET1NRK85/9c1X25CfWyY75KIFhLRs0T0KhHtIqI/0V7TEN9liT420nf5LiL6ORHtcPr5n7XXVP5dMnNm/sFOW/s6gMUA2gDsAHC5cc4nATwJu2zlcgD9pV4LYC2AP613/+Lop/PcxwAsAzBkvOa/AviK8/dXAPxtA/axYb5LAPMBLHP+bgewR/vNNsR3WaKPjfRdEoD3OH+3AugHsLza7zJrFrhbSJmZxwGoQso6NwH4R7bZCuA8Ipof8bVpoZp+gplfAHDC57o3Afie8/f3APxuEo2PSFJ9TBsV95OZDzHzdgBg5tMAXoVdi1a9JvPfZYk+po1q+snM/CvnnFbnH2uvqei7zJqA+xVSNr/soHNKvfaPnSnPd+o9HUV1/QxjLjMfAgDn/wuqbGc1JNVHoAG/SyJaBOA3YVtuQAN+lz59BBrouySiPBG9DOAIgKeZuervMmsCTj7HzDjIoHPCXvu/AVwM4EMADgH4ZoXti4tq+pkVkupjw32XRPQeAJsB3MvMb8fYtrhIqo8N9V0yc4GZPwSgE8AVRNRTbYOyJuCjABZqjzsBHIx4TuBrmfmw8+FaAL4Ne6pUT6rpZxiHlQvC+f9Ile2shkT62GjfJRG1wha29cy8RTunYb7LoD422nepYOaTAJ4DcL1zqOLvMmsC7hZSJqI22IWUHzPOeQzAHzirwcsBnHKmJYGvVR+ew2cADKG+VNPPMB4D8Dnn788BeDTORpdJIn1spO+SiAjAPwB4lZn/zuc1mf8uw/rYYN/lHCI6DwCIaDqAawH8UntNZd9lkqu2SfyDvcq7B/Zq8NedY18C8CWeWu2933l+J4C+sNc6x//JOfcXzoc5P+P93AB7yjkB2yL4gnN8NoAfA3jN+f+9DdjHhvkuAXwU9vT7FwBedv59spG+yxJ9bKTv8gMA/tXpyxCANdo1K/4uZSu9IAhCRsmaC0UQBEFwEAEXBEHIKCLggiAIGUUEXBAEIaOIgAuCIGQUEXBBEISMIgIuCIKQUf4/C2jy0SyUHZoAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plt.plot(res._hat_matrix_diag, res.resid_distance, \".\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d7e3a79d",
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.10"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment