{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Homework 9 Prinicipal Component Regression and Partial Least Squares\n", "\n", "## Part 1 - Data preparation\n", "\n", "We start this homework with probably the most important step in data science: data preparation.\n", "This week, we are going to investigate a data set that contains baseball data from the (North American) Major League during 1986 and 1987." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Task**: Import the Hitters data set, available on the class web page and drop all rows containing missing values." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import pandas as pd\n", "import numpy as np\n", "df = pd.read_csv('./datasets/Hitters.csv').dropna(axis=0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Task**: Identify the three variables containing categorical variables and store their labels in a list `dummy_vars`." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "scrolled": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " AtBat Hits HmRun Runs RBI Walks Years CAtBat CHits CHmRun CRuns \\\n", "1 315 81 7 24 38 39 14 3449 835 69 321 \n", "2 479 130 18 66 72 76 3 1624 457 63 224 \n", "3 496 141 20 65 78 37 11 5628 1575 225 828 \n", "4 321 87 10 39 42 30 2 396 101 12 48 \n", "5 594 169 4 74 51 35 11 4408 1133 19 501 \n", "\n", " CRBI CWalks League Division PutOuts Assists Errors Salary NewLeague \n", "1 414 375 N W 632 43 10 475.0 N \n", "2 266 263 A W 880 82 14 480.0 A \n", "3 838 354 N E 200 11 3 500.0 N \n", "4 46 33 N E 805 40 4 91.5 N \n", "5 336 194 A W 282 421 25 750.0 A \n", "\n", "Int64Index: 263 entries, 1 to 321\n", "Data columns (total 20 columns):\n", "AtBat 263 non-null int64\n", "Hits 263 non-null int64\n", "HmRun 263 non-null int64\n", "Runs 263 non-null int64\n", "RBI 263 non-null int64\n", "Walks 263 non-null int64\n", "Years 263 non-null int64\n", "CAtBat 263 non-null int64\n", "CHits 263 non-null int64\n", "CHmRun 263 non-null int64\n", "CRuns 263 non-null int64\n", "CRBI 263 non-null int64\n", "CWalks 263 non-null int64\n", "League 263 non-null object\n", "Division 263 non-null object\n", "PutOuts 263 non-null int64\n", "Assists 263 non-null int64\n", "Errors 263 non-null int64\n", "Salary 263 non-null float64\n", "NewLeague 263 non-null object\n", "dtypes: float64(1), int64(16), object(3)\n", "memory usage: 43.1+ KB\n", "None\n" ] } ], "source": [ "print(df.head())\n", "print(df.info())\n", "dummy_vars = ['League', 'Division', 'NewLeague']" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The pandas function `get_dummies` converts categorical variables into 0-1-dummy variables. This has been discussed in Chapter 3 (slide 112).\n", "\n", "**Task**: Convert the three categorical variables in the dataset into dummy variables and store them in a new `DataFrame` called `df_dummy`.\n", "Take a look at the new `DataFrame` using the method `head`." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " League_A League_N Division_E Division_W NewLeague_A NewLeague_N\n", "1 0 1 0 1 0 1\n", "2 1 0 0 1 1 0\n", "3 0 1 1 0 0 1\n", "4 0 1 1 0 0 1\n", "5 1 0 0 1 1 0\n" ] } ], "source": [ "df_dummy=pd.get_dummies(df.loc[:,dummy_vars])\n", "print(df_dummy.head())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Once you have done this, you should see that there are only two categories in each of the dummy variables.\n", "Thus we should only include one of each into our final data frame.\n", "\n", "**Task**: If you did everything right so far, the following code should execute without errors." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "real_vars = ['AtBat', 'Hits', 'HmRun', 'Runs', 'RBI', 'Walks', 'Years', 'CAtBat',\\\n", " 'CHits', 'CHmRun', 'CRuns', 'CRBI', 'CWalks','PutOuts', 'Assists', 'Errors']\n", "dfX = pd.concat([df.loc[:,real_vars], df_dummy.loc[:,['League_A', 'Division_E', 'NewLeague_A']]], axis=1)\n", "dfy = df.Salary" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Part 2 - Applying PCR" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We have studied intensively PCA during the lab.\n", "We will now combine this knowledge to perform a PCR.\n", "\n", "According to [Wikipedia](https://en.wikipedia.org/wiki/Principal_component_regression), the PCR method may be broadly divided into three major steps:\n", "\n", "1. Perform PCA on the observed data matrix for the explanatory variables to obtain the principal components, and then (usually) select a subset, based on some appropriate criteria, of the principal components so obtained for further use.\n", "2. Now regress the observed vector of outcomes on the selected principal components as covariates, using ordinary least squares regression (linear regression) to get a vector of estimated regression coefficients (with dimension equal to the number of selected principal components).\n", "3. Now transform this vector back to the scale of the actual covariates, using the selected PCA loadings (the eigenvectors corresponding to the selected principal components) to get the final PCR estimator (with dimension equal to the total number of covariates) for estimating the regression coefficients characterizing the original model." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Task**: Scale the data using `StandardScaler` from `sklearn.preprocessing` and perform a (full) principal component analysis." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/LOCAL/Software/DataScience2018/miniconda3/envs/DS2018/lib/python3.6/site-packages/sklearn/preprocessing/data.py:617: DataConversionWarning: Data with input dtype uint8, int64 were all converted to float64 by StandardScaler.\n", " return self.partial_fit(X, y)\n", "/LOCAL/Software/DataScience2018/miniconda3/envs/DS2018/lib/python3.6/site-packages/sklearn/base.py:462: DataConversionWarning: Data with input dtype uint8, int64 were all converted to float64 by StandardScaler.\n", " return self.fit(X, **fit_params).transform(X)\n" ] } ], "source": [ "from sklearn.preprocessing import StandardScaler\n", "from sklearn.decomposition import PCA\n", "\n", "X = StandardScaler().fit_transform(dfX)\n", "y = dfy.values\n", "pca = PCA().fit(X)\n", "pc = pca.transform(X)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now you should implement a loop over the number of principal components in your model.\n", "We want to measure the quality by a cross-validated mean squared error using 10 folds.\n", "\n", "**Task**:\n", "Implement a loop over the number of principal components in your model.\n", "Use the function `LinearRegression` as an estimator in `cross_val_score`.\n", "As data, you should choose the first $j$ principal components.\n", "Store the means of the mean squared errors in a list called `mse`.\n", "You can use an appropriate `scoring` option." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "123924.13089980236\n", "125219.8186020469\n", "125936.4971858609\n", "123338.38059497124\n", "120267.60729131897\n", "118508.15848939035\n", "118789.33335711621\n", "119588.53890740126\n", "120446.21349148313\n", "122072.10474467053\n", "122818.16267921469\n", "123533.98633524493\n", "123501.98662363982\n", "119791.38051645881\n", "120171.99674286325\n", "116609.51312630426\n", "116359.68557326551\n", "115083.91154069177\n", "116599.0136738026\n", "[123924.13089980236, 125219.8186020469, 125936.4971858609, 123338.38059497124, 120267.60729131897, 118508.15848939035, 118789.33335711621, 119588.53890740126, 120446.21349148313, 122072.10474467053, 122818.16267921469, 123533.98633524493, 123501.98662363982, 119791.38051645881, 120171.99674286325, 116609.51312630426, 116359.68557326551, 115083.91154069177, 116599.0136738026]\n" ] } ], "source": [ "from sklearn.linear_model import LinearRegression\n", "from sklearn.model_selection import cross_val_score\n", "N = pca.n_components_\n", "\n", "mse = []\n", "mseminus = []\n", "mseplus= []\n", "\n", "for i in range(N):\n", " lr = LinearRegression()\n", " cvs = -1*cross_val_score(lr, pc[:,:i+1], y, cv=10, scoring='neg_mean_squared_error')\n", " mse.append(cvs.mean())\n", " mseminus.append(cvs.mean() - cvs.std())\n", " mseplus.append(cvs.mean() + cvs.std())\n", " print(cvs.mean())\n", "print(mse)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Task**: Determine the number of components for which the MSE is smallest." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Minimum MSE for 18 principal components\n" ] } ], "source": [ "print('Minimum MSE for %d principal components' % (np.argmin(mse)+1))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You should observe that the MSE is minimized by taking all but one principal components into consideration.\n", "This corresponds to no dimensionality reduction at all, and simply performs a linear regression using all of the variables.\n", "But we also observe that the values do not change very much, even using only one predictor yields a good fit." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Task**: Plot the MSE against the number of components in your model." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAD8CAYAAACLrvgBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAIABJREFUeJzt3Xd8VFXe+PHPNx1CIKQCoYQUeidCKCpSFMsK6rp2eXbZB3Ut67ruqo9u+215Vt3VtbsorlgesAOWlaYuSg+9k1ATSkISCCWkzvn9MTc4xoRMkpm5ycz3/XrlNTPnnnPnOzcD39xzzj1XjDEopZRS7giyOwCllFKthyYNpZRSbtOkoZRSym2aNJRSSrlNk4ZSSim3adJQSinlNk0aSiml3KZJQymllNs0aSillHJbiN0BeFpcXJxJTk62OwyllGpV1q1bV2iMiW+ont8ljeTkZLKysuwOQymlWhUROeBOPe2eUkop5TZNGkoppdymSUMppZTbGkwaIvKaiBSIyFaXsidFZKeIbBaRj0Qk2mXbIBFZKSLbRGSLiERY5cOt1zki8qyIiFUeIyKLRSTbeuxolYtVL8d6n2Ge//hKKaUaw50zjdeBybXKFgMDjDGDgN3AIwAiEgK8BdxpjOkPjAMqrTYvATOAdOunZp8PA0uNMenAUus1wOUudWdY7ZVSStmowaRhjFkGFNcqW2SMqbJergK6Ws8vBTYbYzZZ9YqMMdUi0hlob4xZaZx3fXoDmGq1mQLMtp7PrlX+hnFaBURb+1FKKWUTT4xp/AT4t/W8F2BEZKGIrBeRX1vlSUCeS5s8qwwg0RhzBMB6THBpk1tPG6WUUjZo1nUaIvIoUAW87bK/scAFQCmwVETWASfraN7QfWbF3TYiMgNnFxbdu3dvOPAWbMPB45wpr2ZsepzdoSil1Pc0+UxDRKYBVwG3mG9vNJ4H/McYU2iMKQU+A4ZZ5V1dmncFDlvP82u6nazHApd9daunzXcYY2YaYzKMMRnx8Q1e0NgiORyGF7/K4bqXVvDfb2Rxpryq4UZKKeVjTUoaIjIZeAi42koONRYCg0SkrTUofjGw3ep2OiUimdasqduB+VabBcA06/m0WuW3W7OoMoGSmm4sf1NytpIZb67jic93MbR7R85WVrN4e77dYSml1Pe4M+V2DrAS6C0ieSIyHXgeiAIWi8hGEXkZwBhzHHgKWAtsBNYbYz61dnUX8CqQA+zh23GQvwKTRCQbmGS9BudZyl6r/ivAz5r5WVukbYdLuPr5b/hqVwG/vaof790xiqToNszbeMju0JRS6nsaHNMwxtxUR/Gs89R/C+e029rlWcCAOsqLgAl1lBvg7obia83ey8rlsXlbiW4bytwZmWQkxwAwZUgX/rlsL4Wny4lrF25zlEop9S29ItwGZZXVPPzBZn71/maG9+jIp/ddeC5hAFwzNIlqh+GTTXUO4SillG00afjYwaJSrntpBXPX5nL3Jam8OX3k984m0hOj6Ne5PR9t1KShlGpZNGn40NId+Vz13NccLC7l1dsz+NVlfQgOqmtmMUwd2oVNuSfYV3jGx1EqpVT9NGn4QLXD8OTCnUyfnUW3mLZ8eu+FTOyXeN42Vw9OQgTmbdABcaVUy6FJw8sKT5dz+2ureeHLPdyQ0Y0P7hpN99i2Dbbr1CGCUSmxzN94iG8vg1FKKXtp0vCidQeKuerZb8jaf5wnrhvE4z8cRERosNvtpw5NYn9RKRtzT3gxSqWUcp8mDS8wxvCv5fu44Z+rCAsJ4oO7RvOjC7o13LCWyQM6ERYSxHwdEFdKtRCaNDzsTHkV987ZwB8+3s643vF8fO9YBiR1aNK+2keEMqlvIh9vOkxltcPDkSqlVONp0vCgfYVnmPLCcj7bcoRfT+7NzNsy6NAmtFn7nDKkC0VnKvgmp9BDUSqlVNNp0vCQsxXVzHgji6LT5bz105H8bFwaQfVMp22Mcb0TiG4bqrOolFItgiYND/nzZ9vJLjjNszcNZXSq55Y1DwsJ4oqBnVm0LV9XvlVK2U6Thgcs3p7PW6sOMuOiFC5M9/zS7NcMTeJsZTWLth/1+L6VUqoxNGk0U/7JMn79/ib6d2nPg5f29sp7DO/e0bny7QadRaWUspcmjWZwOAy/fHcTZyureebGoYSFeOdwBgUJU4d24evsYxw7Ve6V91BKKXdo0miGWd/s45ucQn73g/6kJbTz6ntNHZKEw8Anm/VsQyllH00aTbT1UAlPLNzJZf0TubEJF+41VnpiFP27tGeeXuinlLKRJo0mKK2o4r65G4iNDOev1w7CeQdb75s6JElXvlVK2UqTRhP88ZPt7Cs8w1M3DKZjZJjP3vcHg7voyrdKKVtp0mikz7ceYc6aXO64KNWj12O4o1OHCEanxjJPV75VStlEk0YjHCk5y0MfbGFQ1w48MKmXLTFMGZLEAV35VillE00abqp2GB54ZxOV1Q6vTq9tyOQBnQgPCdIuKqWULTRpuGnmsr2s3FvE76/uT8+4SNviaB8RysS+iXyy+YiufKuU8jlNGm7YlHuCvy/axZUDO3P98K52h8PUoUnOlW+zdeVbpZRvadJowJnyKn4+dwMJUeH85ZqBPpteez4X94p3rny7UbuolFK+pUmjAb9fsI0DxaU8fcMQOrRt3r0xPCUsJIgrdeVbpZQNNGmcxyebD/PeujzuHpfGyJRYu8P5jqm68q1SygaaNOpx6MRZHvlwC0O6RfPziel2h/M9w7t3pGvHNnykK98qpXxIk0Ydqh2GX8zdiMNheObGIYQGt7zDFBQkTBnShW905VullA81+L+hiLwmIgUistWl7EkR2Skim0XkIxGJrtWmu4icFpEHXcomi8guEckRkYddynuKyGoRyRaRd0QkzCoPt17nWNuTPfGB3fHSVzms2V/M/5sygB6x9k2vbUjNyrcfb9KzDaWUb7jzJ/TrwORaZYuBAcaYQcBu4JFa258G/l3zQkSCgReAy4F+wE0i0s/a/DjwtDEmHTgOTLfKpwPHjTFp1v4ed/MzNcv6g8d5ekk2Vw/uwrXDknzxlk1Ws/LtfJ1FpZTykQaThjFmGVBcq2yRMaZm2s4q4NzFCyIyFdgLbHNpMgLIMcbsNcZUAHOBKeKcvzoeeN+qNxuYaj2fYr3G2j5BvDzf9VRZJffP3Uin9hH86ZoBLWJ6bUOuGZrEprwS9h47bXcoSqkA4InO+p9gnVWISCTwEPCHWnWSgFyX13lWWSxwwiUB1ZR/p421vcSq7zW/m7+NvOOlPHPjENpHtIzptQ05t/Kt3mdDKeUDzUoaIvIoUAW8bRX9AWdXU+0/e+v6k92cp/x8beqKY4aIZIlI1rFjxxoOvA7zNx7iww2HuHd8OhnJMU3ahx0S20cwJjWOeRt05VullPc1OWmIyDTgKuAW8+3/ViOBJ0RkP3A/8D8icg/OMwjX29t1BQ4DhUC0iITUKse1jbW9A7W6yWoYY2YaYzKMMRnx8fFN+jzx7cK5YmAn7h2f1qT2dpoypAsHi0vZoCvfKqW8rElJQ0Qm4+yGutoYU1pTboy50BiTbIxJBv4B/MUY8zywFki3ZkqFATcCC6xk8yXwQ2sX04D51vMF1mus7V8YL/4pPTotjhdvGU5IC5xe2xBd+VYp5SvuTLmdA6wEeotInohMB54HooDFIrJRRF4+3z6sMYl7gIXADuBdY0zNQPlDwAMikoNzzGKWVT4LiLXKHwAeRtUpKiKUif105VullPeJv/WDZ2RkmKysLLvD8LnF2/P57zeyeO2/MhjfJ9HucJRSrYyIrDPGZDRUr/X1xag6nVv5VpcVUUp5kSYNP3Fu5dvtRzmtK98qpbwkpOEqqrW4ZmgSb68+yKJtR7l2mP03i1Itz2dbjjBnzUGCg4SIkGAiQoOICA0mIjSY8NAgwmvKQoKtcpcyl9dJ0W3oGBlm98dRNtCk4UeG96hZ+faQJg31HRVVDv733zv41/L99IhtS4c2oZRVVlNW6bAeqymrclBR5d5Eirh2Yax6ZEKrnG2omkeThh8REaYOSeLFr3IoOFVGQlSE3SGpFuBIyVnufns96w+e4CdjevLIFX3qXbnZ4TBUVDu+m1Cqqik/99zBmn1FvPDlHrYfOcmgrtF17kf5L00afmbq0C48/2UOH286wvSxPe0OR9lseU4h983ZQFllNc/fPJSrBnU5b/2gICEiyNk1VZ++naJ44cs9rNpbpEkjAOm5pZ9JS3CufKvLpQc2h8Pw/BfZ3DZrNTGRYcy/Z2yDCcNdCe0jSImPZNXeOhdoUH5Ok4Yfmty/ExtzT1BwqszuUJQNSkor+ekbWfxt0W5+MLgL8+4eQ1pCO4++R2ZKLGv3FVOlF5MGHE0afmhiP+fFfV/sKLA5EuVrW/JKuPK5r/k6+xh/nNKff9wwhMhwz/dCZ6bEcqq8iu1HTnp836pl06Thh/p0iiIpug1LNGkEDGMMc9Yc5LqXV+BwGN69YxS3jUr22j1hMns6V4JetbfIK/tXLZcmDT8kIkzom8A3Occoq6y2OxzlZWcrqnnwvc088uEWRvaM4ZP7LmRo945efU8d1whcmjT81MS+iZRVOlieU2h3KMqL9hee4ZoXl/Phhjzum5DO6z8eQYyPLrrTcY3ApEnDT41MiSEyLFi7qPzYwm1H+cFz33D0ZBn/+q8LeGBSL4KDfHeL4lE6rhGQNGn4qfCQYC7uHc8XO/NxOPxrJeNAV1Xt4H8/28Edb66jZ3wkn9w7lnG9E3wex8gUHdcIRJo0/NiEPonknyxn6+ESu0NRHlJwsoybX13NP5ft5dbM7rx35yi6dmxrSywJURGk6rhGwNErwv3YJX0SCBJYsqNAr9z1A2v3F/Ozt9dzqqySp28YzDVD7V9fLDMllgUbD1NV7dB1qAKE/pb9WExkGMN7dGTJ9ny7Q1HNYIzhrVUHuGnmKtqFhzD/7rEtImGAXq8RiDRp+LkJfRPZfuQkh0+ctTsU1QTlVdU88uEWHpu3lQvT45h39xh6d4qyO6xzdFwj8GjS8HMT+zqvDl+6U2dRtTb5J8u4ceYq5q7N5Z5L0nh12gV0aBNqd1jfoeMagUeThp9LjY8kObatdlG1MusOHOeq575h19FTvHjLMB68rLdPp9M2hl6vEVg0afg559XhiazcU8QZvQ1sqzB3zUFunLmSNqHBfPSzMVwxsLPdIZ2XjmsEFk0aAWBi30Qqqh18nX3M7lDUeVRUOXhs3hYe/nCLc1bSPS1r/KI+Oq4RWDRpBICM5I60jwjRq8NbsIJTZdzy6ireWnWQOy9O5fUfjyC6beu4B3fNuMbKPZo0AoFepxEAQoODGNc7gS93FlDtMC22bzxQbcw9wZ1vruPE2QqevWkoVw/2zM2SfCkzJZb5er1GQNDfboCY2C+RojMVbMw9bncoysV7Wbn86J8rCQkWPrxrTKtMGOBMGqfLq9h2WMc1/J0mjQBxca94QoJEu6haiMpqB79fsI1fvb+ZC5I78vE9Y+nXpb3dYTWZjmsEDk0aAaJDm1BG9Ixh6Q6demu3otPl3Prqal5fsZ+fju3J7B+PoKOPljP3lm+v19Ck4e8aTBoi8pqIFIjIVpeyJ0Vkp4hsFpGPRCTaKp8kIutEZIv1ON6lzXCrPEdEnhXrlmIiEiMii0Uk23rsaJWLVS/Hep9hnv/4gWVC30R255/mYFGp3aEErK2HSrj6+eVszD3B0zcM5rGr+vnNGEBmSixr9x/X6zX8nDvf1teBybXKFgMDjDGDgN3AI1Z5IfADY8xAYBrwpkubl4AZQLr1U7PPh4Glxph0YKn1GuByl7ozrPaqGSb2dS6fvUTPNmzx0YY8rntpBcYYPrhrdItZP8pTdFwjMDSYNIwxy4DiWmWLjDE1V4qtArpa5RuMMYet8m1AhIiEi0hnoL0xZqUxxgBvAFOtelOA2dbz2bXK3zBOq4Boaz+qiXrERpKe0E6Thg0WbTvKL97ZxJBu0Sy4dywDkjrYHZLH6bhGYPDEefFPgH/XUX4dsMEYUw4kAXku2/KsMoBEY8wRAOux5m4ySUBuPW1UE03om8iafcWcLKu0O5SAUV5VzZ8+3UGvxHa89dORxLULtzskr9BxjcDQrKQhIo8CVcDbtcr7A48Dd9QU1dG8odvJud1GRGaISJaIZB07plc9n8+kfglUOQz/2aXHyVdeX76fg8Wl/OaqfoT6yfhFfXRcw/81+RssItOAq4BbrC6nmvKuwEfA7caYPVZxHlYXlqUrUNONlV/T7WQ9Fri06VZPm+8wxsw0xmQYYzLi4+Ob+pECwpBuHYmJDNMuKh85dqqc577IYWLfBC5M9//vpo5r+L8mJQ0RmQw8BFxtjCl1KY8GPgUeMcYsrym3up1OiUimNWvqdmC+tXkBzkFzrEfX8tutWVSZQElNN5ZquuAg4ZLeCXy165j+NegDTy3eRVllNf9zRV+7Q/EJHdfwf+5MuZ0DrAR6i0ieiEwHngeigMUislFEXraq3wOkAb+xyjeKSM0YxV3Aq0AOsIdvx0H+CkwSkWxgkvUa4DNgr1X/FeBnzfuoqsakfgmUnK0k64BeHe5N2w6XMHdtLtNGJ5MS387ucHxCxzX8X4NrTxljbqqjeFY9df8E/KmebVnAgDrKi4AJdZQb4O6G4lONd2F6PGHBQSzZnk9mSqzd4fglYwx//GQ70W1CuW98ut3h+JSuQ+Xf9DcagCLDQxiVGsuSHfm4DEcpD1q0PZ9Ve4t5YFIvOrRtWXfb8zYd1/BvmjQC1MS+CewvKmXPsTN2h+J3yquq+ctnO0hPaMdNI7rbHY7P6biGf9OkEaDG19w7XGdRedzsFfs5UOScYhuI3TMJURGkJbTTpOGnAu8brQBIim5Dv87tWaqr3npU4elynluaw/g+CVzUy/+n2NYnMyVGr9fwU5o0AtjEvglkHSjm+JkKu0PxG39ftJuzldU8emVgTLGtj45r+C9NGgFsQt9EHAa+3KVnG56w/fBJ3ll7kNtHJZMaIFNs6zOyp3NW3krtovI7mjQC2MCkDiREhevV4R5QM8W2fZtQfj4hsKbY1iU+KlzHNfyUJo0AFhQkTOibwLLdhVRUad9zcyzens/KvUUBOcW2PpkpMazdV6zjGn5Gk0aAm9g3kdPlVazep38RNlV5VTV/tqbY3hyAU2zrk5kSy5mKarbquIZf0aQR4MakxRER6rw6XDVNzRRbf7oLnyfUjGtoF5V/0W94gIsIDWZsWhxLdhTo1eFNUDPF9pLe8VwcwFNs66LjGv5Jk4ZiYt9EDp04y678U3aH0uo8tbhmim0/u0NpkXRcw/9o0lCM72PdO1y7qBplx5GTzF1zkNtG9SAtIbCn2NZHxzX8jyYNRUL7CAZ3i2aJXh3uNmMMf/pUp9g2RMc1/I8mDQXAxD4JbMw9QcGpMrtDaRWW7ChgeU4Rv5jYi+i2YXaH02LpuIb/0aShAOfV4QBf7tSzjYZUVDn486fbSUtox80jdYptQ3Rcw79o0lAA9O0cRVJ0G+2icsPsFfvZX1TKY1f2JVSn2DZIxzX8i37jFQAizqvDv84+Rllltd3htFhFp8t5dmk243rHM653QsMNlI5r+BlNGuqcCX0TKat0sGJPod2htFhPLd5NaWU1j+kUW7fpuIZ/0aShzslMiSEyLJjF27WLqi47j55kzpqD3JapU2wbS8c1/IcmDXVOeEgwF/WK54udeu/w2mpWsY2KCOX+iTrFtrF0XMN/aNJQ3zGxbyL5J8vZekj/cbtaem6KbbpOsW0CHdfwH5o01Hdc0ieBIIHFeo+Nc8oqnavYpsZHcktmD7vDaZV0XMN/aNJQ3xETGcaw7h1ZqkkDgJKzldw+aw37i87wux/01ym2zaDjGv5B/wWo75nYL5Fth09ypOSs3aHYquBUGTfOXMWG3OM8d9NQLtJVbJulZlxjy6ESu0NRzaBJQ33PxL7WAoYBfKFfbnEp17+8kv2FZ3h12gVcNaiL3SG1et+OaxTbHIlqDk0a6ntS49vRMy6SD9blBeQsql1HT3HdSys4UVrJ2/89Uu+T4SE6ruEfNGmo7xERZlyUwsbcEwF3trHuwHF+9M+VALx7xyiGde9oc0T+JTMlhqz9xVQ2c1yjrLKaBZsOk7Vfz1p8rcGkISKviUiBiGx1KXtSRHaKyGYR+UhEol22PSIiOSKyS0QucymfbJXliMjDLuU9RWS1iGSLyDsiEmaVh1uvc6ztyZ760Kph1w/vSs+4SP62cBfVjsA42/jP7mPc+upqOrYN5YO7RtO7U5TdIfmdc9drNHFc40DRGf7y2Q5G/e9S7puzgQfe3eThCFVD3DnTeB2YXKtsMTDAGDMI2A08AiAi/YAbgf5WmxdFJFhEgoEXgMuBfsBNVl2Ax4GnjTHpwHFgulU+HThujEkDnrbqKR8JCQ7igUm92JV/igWbDtkdjtd9vOkwP529luS4SN67czTdYtraHZJfykxp/LhGVbWDhduOctus1Vz85FfM+mYfmSmx3HhBNw4Wl5JbXOqtcFUdGkwaxphlQHGtskXGmCrr5Sqgq/V8CjDXGFNujNkH5AAjrJ8cY8xeY0wFMBeYIiICjAfet9rPBqa67Gu29fx9YIJVX/nIlQM7069ze55avJuKKv+dJvnWqgPcN3cDQ7t1ZO6MTOKjwu0OyW/FtQsn3c1xjfyTZTyzJJsLn/iSO95cR3b+aX4xsRcrHh7PS7cO5ydjewKwco+OkfhSiAf28RPgHet5Es4kUiPPKgPIrVU+EogFTrgkINf6STVtjDFVIlJi1f/eanoiMgOYAdC9u97fwFOCgoRfTe7Nj/+1lnfWHuS2Ucl2h+RRxhhe+DKHvy3azYQ+CbxwyzAiQoPtDsvvZabE8uH6PCqrHd+77sUYw4o9Rby16gCLtudT7TBcmB7H76/uz4Q+CYS41E9PaEdcu3CW7ynkRxd08/XHCFjNShoi8ihQBbxdU1RHNUPdZzTmPPXPt6/vFxozE5gJkJGRERgd8D4yrlc8I3rG8MzSHK4b3pW2YZ74O8N+Dofhz5/tYNY3+7hmaBJP/HCQXrjnI5kpsby56gBbD5Uw1JpocKK0gvfX5fF/qw+yt/AMHduGMn1sT24e0Z3kuMg69yMijE6NZcWeIowxaEeEbzT5fwARmQZcBUww387LzANcU35X4LD1vK7yQiBaREKssw3X+jX7yhOREKADtbrJlPeJCA9N7s11L63kX8v3c/claXaH1GyV1Q4e+mAzH64/xI/HJPObK/sRFKT/4fjKyJQYAFZaXVRvrz7Ix5sOU17lYFj3aJ760WCuGNjZrbO+MWmxLNh0mJyC06Qn6sQFX2hS0hCRycBDwMXGGNdRqAXA/4nIU0AXIB1Yg/OsIV1EegKHcA6W32yMMSLyJfBDnOMc04D5LvuaBqy0tn9hAvGigRZgeI8YJvRJ4J//2cOtI3vQoW2o3SE1WVllNff833qW7CjggUm9uHd8mv6F6mM14xr/WJzNE9W7aBsWzHXDu3LryB7069K+UfsanRoHwPKcQk0aPuLOlNs5OP/j7i0ieSIyHXgeiAIWi8hGEXkZwBizDXgX2A58DtxtjKm2ziLuARYCO4B3rbrgTD4PiEgOzjGLWVb5LCDWKn8AODdNV/neg5f15lR5FS8v22N3KE12sqyS219bw9KdBfxxSn/um5CuCcMmN1zQjb6do/jjlP6s/p8J/OWagY1OGADdYtrSPaYty3Uw3GfE3/54z8jIMFlZWXaH4Zd+PncDC7cdZdmvLiGhfYTd4TRK4elypr22hl1HT/H3Hw1mypCkhhupVuHhDzbz6ZYjbPjNpO8MlKvGEZF1xpiMhurpEVZue2BSL6qqDc99kWN3KI1Ss47UnmOneXVahiYMPzM6LY5TZVVs0xs8+YQmDeW2HrGR3HBBN+asOcjBotZxQdX2wye5/uWVFJ0u5+2fjmRc7wS7Q1IeNsq6YHC53tveJzRpqEa5b0I6IcHC00t22x1Kg77cVcD1L68A4N07RzG8R4zNESlviI8Kp3diFCtydFzDFzRpqEZJbB/BtNHJzNt4iJ1HW253wJurDvDT2Vkkx0Uy7+4x9OnU+EFW1XqMTotl7f5iyquq7Q7F72nSUI1218WptAsP4W8LW97ZhsNh+POn2/nNvK2M6xXPu3eMolOH1jVorxpvTGoc5VUO1h84YXcofk+Thmq06LZh3HFRCkt25LPuwHG7wznnbEU1d729jle+3se0UT2YeXsGkeH+cQW7Or8RKTEECazQcQ2v06ShmuTHY3oS1y6MJz7f2SJu1OS8NetKFm3P53c/6McfpgwgWK/yDhjtI0IZ1DWa5TmaNLxNk4ZqksjwEO4dn87qfcUsy7b3H+ru/FNc88IKduefZuZtGfx4TE9b41H2GJMWy6a8Ek6XVzVcWTWZJg3VZDeN6E7Xjm14cuFOHDbdqGl5TiHXvbSCimoH794xikn9Em2JQ9lvTGoc1Q7Dmn06i8qbNGmoJgsLCeIXE3ux9dBJ/r31qM/f/921uUx7bQ1dOrRh3t1jGNi1g89jUC3HsB4dCQsJYrlOvfUqTRqqWaYOTaJXYjv+vngXVc2877O7HA7DE5/v5NcfbGZUaizv3zWKpOg2Pnlv1XJFhAaT0aMjK3QdKq/SpKGaJThI+OWlvdl77AwfrM/z+vuVVVZz39wNvPjVHm4a0Y3X/usCoiJa76q7yrPGpMWx48hJik6X2x2K39KkoZrt0n6JDOkWzT+WZFNW6b2Lq4pOl3PLq6v5ZPMRHrm8D3+5ZqDeOEl9x6hU55IiK924naxqGv0Xp5pNRPj1Zb05UlLGW6sOeOU99hw7zbUvrWDroRJevGUYd1ycqsuaq+8ZlNSBqPAQHdfwIk0ayiNGp8UxNi2OF7/a4/Epj6v3FnHtiys4XVbFnBmZXDGws0f3r/xHSHAQI1NiWKkX+XmNJg3lMb+6rDfFZyp49eu9HtlfWWU1c9Yc5LZZa4hrF8ZHPxvDMOue0krVZ3RqHPuLSjl04qzdofglXWNBeczgbtFM7t+JV5bt5bbMHsS2C2/0PgpPl/PFzgKW7sjn6+xCSiuqyUyJ4Z+3ZrTq28wq3xlCwml2AAARXElEQVSdZi2VnlPIjzK62RyN/9GkoTzqwct6sWj7UV78ag+/uapfg/WNMWQXnGbJjnyWbM9nQ+4JjIHOHSK4blhXJvRNYGxanN6RTbmtd2IUce3CWLmnSJOGF2jSUB6VlhDFtcO68uaqA0wf25MudVw/UVntYO2+YhbvyGfpjgIOFjtv6DQwqQP3T+jFxH4J9OvcXge6VZOICKNS41ieU4gxRr9HHqZJQ3nc/RPTWbDxMM8syebxHw4CoKS0kq92F7BkRwFf7SrgVFkVYSFBjE2L446LU5jQJ1GXMFceMzo1lo83HWbPsdOkJUTZHY5f0aShPK5rx7bcktmd2Sv20yW6Dav2FrFmfzHVDkNcuzAuH9CJiX0TGZseR9sw/QoqzxuTGgfA8pwiTRoepv9ilVfcfUka767N5eklu+mdGMUdF6UwsV8iQ7pGE6RLlisv6x7blq4d27BiTyHTRifbHY5f0aShvCKuXTif3HchwSJ0j21rdzgqAI1OjeXzrUepdhi9t4oH6ZQU5TU94yI1YSjbjEmL42RZFdsOl9gdil/RpKGU8ks161DpqreepUlDKeWXEqIi6JXYTm8B62GaNJRSfmt0ahxr9xdTXuW91ZcDTYNJQ0ReE5ECEdnqUna9iGwTEYeIZLiUh4rIbBHZIiI7ROQRl22TRWSXiOSIyMMu5T1FZLWIZIvIOyISZpWHW69zrO3JnvrQSqnAMDo1lrJKBxsOnrA7FL/hzpnG68DkWmVbgWuBZbXKrwfCjTEDgeHAHSKSLCLBwAvA5UA/4CYRqVlj4nHgaWNMOnAcmG6VTweOG2PSgKetekop5baRKbEEiY5reFKDScMYswworlW2wxizq67qQKSIhABtgArgJDACyDHG7DXGVABzgSnivL5/PPC+1X42MNV6PsV6jbV9guh6AEqpRujQJpSBSR1YoeMaHuPpMY33gTPAEeAg8DdjTDGQBOS61MuzymKBE8aYqlrluLaxtpdY9ZVSym2j0+LYmHuCMx6+z0ug8nTSGAFUA12AnsAvRSQFqOsMwZynnAa2fYeIzBCRLBHJOnbsWOOjVkr5rTGpcVQ5DGv2FzdcWTXI00njZuBzY0ylMaYAWA5k4DyDcF2juCtwGCgEoq3uLNdyXNtY2ztQq5ushjFmpjEmwxiTER8f7+GPpJRqzYb36EhYcJB2UXmIp5PGQWC8OEUCmcBOYC2Qbs2UCgNuBBYYYwzwJfBDq/00YL71fIH1Gmv7F1Z9pZRyW5uwYIb1iNb7hnuIO1Nu5wArgd4ikici00XkGhHJA0YBn4rIQqv6C0A7nLOr1gL/MsZstsYk7gEWAjuAd40x26w2DwEPiEgOzjGLWVb5LCDWKn8AODdNVymlGmNMahzbj5yk+EyF3aG0euJvf7xnZGSYrKwsu8NQSrUg6w4c57qXVvDiLcO4YmBnu8NpkURknTEmo6F6ekW4UsrvDeragciwYF1SxAM0aSil/F5ocBAjU2L1Ij8P0KShlAoIo1Nj2Vd4hsMnztodSqumSUMpFRBGW7eA1bON5tGkoZQKCH06RRETGabXazSTJg2lVEAIChJGpcayfE8h/jZr1Jc0aSilAsaY1DjyT5azt/CM3aG0Wpo0lFIBY3TNLWC1i6rJNGkopQJGj9i2JEW30SVFmkGThlIqYIgIo1NjWbm3CIdDxzWaQpOGUiqgjE6LpeRsJduPnLQ7lFZJk4ZSKqDUXK+hS4o0jSYNpVRASWwfQVpCO5brRX5NoklDKRVwRqfGsnZfMRVVDrtDaXU0aSilAs7o1DjOVlazMfeE3aG0Opo0lFIBZ1RKLEGi4xpNoUlDKRVwOrQNZUBSB1bquEajadJQSgWkUamxbMg9TmlFld2htCqaNJRSAWlMahyV1YY1+4rtDqVV0aShlApIFyTHEBosvLHyAOsOFFNV3XpnUhljeOHLHPJPlnn9vUK8/g5KKdUCtQkL5tbMHry+Yj9f7CwgKjyEkSmxXJgex9j0OFLiIhERu8N0y6dbjvDkwl20Cw9h2uhkr76X+Nu68hkZGSYrK8vuMJRSrcSJ0gpW7Cni6+xCvsk5Rm6x83awXTpEMDY9jrHp8YxJjSW2XbjNkdbtVFklE/7+HxLahzP/7rEEBzUt0YnIOmNMRkP19ExDKRXQotuGccXAzlwxsDMAB4rOOBNIdiGfbz3Ku1l5APTr3P7cWcgFyTFEhAbbGfY5Ty3ezbHT5bxye0aTE0ZjaNJQSikXPWIj6REbya2ZPah2GLYcKuGb7GN8nV3Ia8v38c9lewkLCWJEcgxj0uK4MD2O/l3a29KVtfVQCbNX7OeWkd0Z3C3aJ++p3VNKKeWmM+VVrNlXfK4ra3f+aQCuH96VJ68f7NNYHA7DtS+tIO94KUt/OY4ObUKbtT/tnlJKKQ+LDA/hkj4JXNInAYD8k2W89NUeXl+xn1GpsVw7rKvPYpm7NpeNuSd4+obBzU4YjaFTbpVSqokS20fw2JV9GZEcw2/mbWWfj+49Xni6nMc/30lmSgxThyT55D1raNJQSqlmCAkO4h83DiEkOIj75mzwycq5//vZTkorqvjT1AE+H0tpMGmIyGsiUiAiW13KrheRbSLiEJGMWvUHichKa/sWEYmwyodbr3NE5FmxPqmIxIjIYhHJth47WuVi1csRkc0iMsyzH10ppTyjS3QbHr9uEFsOlfDkwp1efa/Ve4v4YH0e/31hCmkJUV59r7q4c6bxOjC5VtlW4FpgmWuhiIQAbwF3GmP6A+OASmvzS8AMIN36qdnnw8BSY0w6sNR6DXC5S90ZVnullGqRJg/oxK2Z3Xnl6318tavAK+9RUeXgsXlb6dqxDfeOT/fKezSkwaRhjFkGFNcq22GM2VVH9UuBzcaYTVa9ImNMtYh0BtobY1Ya53StN4CpVpspwGzr+exa5W8Yp1VAtLUfpZRqkR67sh+9E6N48L1NFJzy/JIes77ZR3bBaf5wdX/ahNlznYinxzR6AUZEForIehH5tVWeBOS51MuzygASjTFHAKzHBJc2ufW0+Q4RmSEiWSKSdezYMQ99FKWUapyI0GCeu3kop8qq+OW7m3A4PHdJQ97xUp5dms2l/RKZ0DfRY/ttLE8njRBgLHCL9XiNiEwA6hqpaehout3GGDPTGJNhjMmIj49vTLxKKeVRvRKj+O0P+vF1diGvfL3XY/v9w8fbAfjd1f09ts+m8HTSyAP+Y4wpNMaUAp8Bw6xy1wnMXYHD1vP8mm4n67HAZV/d6mmjlFIt1s0junP5gE48uXAXmzxwS9kl2/NZvD2fn09MJym6jQcibDpPJ42FwCARaWsNil8MbLe6nU6JSKY1a+p2YL7VZgEwzXo+rVb57dYsqkygpKYbSymlWjIR4a/XDiIhKpx752zgVFllw43qUVpRxe8WbKNXYjumj+3pwSibxp0pt3OAlUBvEckTkekico2I5AGjgE9FZCGAMeY48BSwFtgIrDfGfGrt6i7gVSAH2AP82yr/KzBJRLKBSdZrcJ6l7LXqvwL8rLkfVimlfKVD21CeuWkoecdLeWzeVpq6ZNNzX+Rw6MRZ/jR1IKHB9l9ap2tPKaWUFz2zJJunl+zm79cP5rrhjVtmJDv/FJc/8zVThybxNy+vbeXu2lP2py2llPJj94xPY0TPGH4zfyt7j512u50xhsfmbSUyPIRHLu/jxQgbR5OGUkp5UXCQ8MyNQwgLCeK+uRsor6p2q92H6w+xel8xD1/ep0XdAEqThlJKeVnnDs5lRrYeOsmTn9d1XfR3lZRW8pfPdjC0ezQ3ZHRrsL4vadJQSikfuKx/J27L7MGr3+zjywaWGXli4U6Ol1bwp6kDCPLB3fgaQ5OGUkr5yKNX9qVPpygefHcTBSfrXmZkY+4J/m/NQf5rdE/6d+ng4wgbpklDKaV8JCI0mOduGsqZiioeqGOZkapqB49+tIWEqHAeuLSXTVGenyYNpZTyofTEKH57VX++ySlkZq1lRt5cdYBth0/y26v60y68Zd5YVZOGUkr52E0junH5gE78beEuNhw8DjhvHfv3Rbu5qFc8VwzsZHOE9dOkoZRSPlazzEhi+wjum7uBk2WV/PGT7VRUO/h/V/f3+d34GkOThlJK2aBD21CeuXEIh46f5bZXV/PJ5iP8bFwqyXGRdod2Xpo0lFLKJhnJMdw/sReb8kpIjm3LnRen2h1Sg1rmSItSSgWIuy9Jw2EME/smEhFqz934GkOThlJK2Sg4SLh/YsucXlsX7Z5SSinlNk0aSiml3KZJQymllNs0aSillHKbJg2llFJu06ShlFLKbZo0lFJKuU2ThlJKKbeJMabhWq2IiBwDDtgdRwPigEK7g3CDxulZrSVOaD2xapye08MYE99QJb9LGq2BiGQZYzLsjqMhGqdntZY4ofXEqnH6nnZPKaWUcpsmDaWUUm7TpGGPmXYH4CaN07NaS5zQemLVOH1MxzSUUkq5Tc80lFJKuU2ThheISDcR+VJEdojINhH5eR11xolIiYhstH5+a0esViz7RWSLFUdWHdtFRJ4VkRwR2Swiw2yIsbfLsdooIidF5P5adWw7piLymogUiMhWl7IYEVksItnWY8d62k6z6mSLyDQb4nxSRHZav9uPRCS6nrbn/Z74IM7fi8ghl9/vFfW0nSwiu6zv68M2xPmOS4z7RWRjPW19djw9yhijPx7+AToDw6znUcBuoF+tOuOAT+yO1YplPxB3nu1XAP8GBMgEVtscbzBwFOe88hZxTIGLgGHAVpeyJ4CHrecPA4/X0S4G2Gs9drSed/RxnJcCIdbzx+uK053viQ/i/D3woBvfjT1AChAGbKr9b8/bcdba/nfgt3YfT0/+6JmGFxhjjhhj1lvPTwE7gCR7o2qWKcAbxmkVEC0inW2MZwKwxxjTYi7iNMYsA4prFU8BZlvPZwNT62h6GbDYGFNsjDkOLAYm+zJOY8wiY0yV9XIV0NVb7++ueo6nO0YAOcaYvcaYCmAuzt+DV5wvThER4EfAHG+9vx00aXiZiCQDQ4HVdWweJSKbROTfItLfp4F9lwEWicg6EZlRx/YkINfldR72JsEbqf8fYks5pgCJxpgj4PxDAkioo05LO7Y/wXlWWZeGvie+cI/VjfZaPd19Lel4XgjkG2Oy69neEo5no2nS8CIRaQd8ANxvjDlZa/N6nN0rg4HngHm+js/FGGPMMOBy4G4RuajWdqmjjS3T7kQkDLgaeK+OzS3pmLqrJR3bR4Eq4O16qjT0PfG2l4BUYAhwBGfXT20t5ngCN3H+swy7j2eTaNLwEhEJxZkw3jbGfFh7uzHmpDHmtPX8MyBUROJ8HGZNLIetxwLgI5yn+K7ygG4ur7sCh30T3fdcDqw3xuTX3tCSjqklv6Ybz3osqKNOizi21gD8VcAtxupwr82N74lXGWPyjTHVxhgH8Eo9799SjmcIcC3wTn117D6eTaVJwwusvsxZwA5jzFP11Olk1UNERuD8XRT5LspzcUSKSFTNc5yDoltrVVsA3G7NosoESmq6XWxQ719vLeWYulgA1MyGmgbMr6POQuBSEelodbdcapX5jIhMBh4CrjbGlNZTx53viVfVGke7pp73Xwuki0hP66z0Rpy/B1+bCOw0xuTVtbElHM8ms3sk3h9/gLE4T4k3AxutnyuAO4E7rTr3ANtwzu5YBYy2KdYUK4ZNVjyPWuWusQrwAs5ZKVuADJtibYszCXRwKWsRxxRnIjsCVOL8a3c6EAssBbKtxxirbgbwqkvbnwA51s+PbYgzB+c4QM139WWrbhfgs/N9T3wc55vW928zzkTQuXac1usrcM5Y3GNHnFb56zXfS5e6th1PT/7oFeFKKaXcpt1TSiml3KZJQymllNs0aSillHKbJg2llFJu06ShlFLKbZo0lFJKuU2ThlJKKbdp0lBKKeW2/w/Af8q8k18JugAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "import numpy as np\n", "%matplotlib inline\n", "nrange = 1+np.asarray(range(N))\n", "plt.plot(nrange,mse)\n", "#plt.fill_between(nrange, mseminus, mseplus, alpha=0.2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Task**: Plot the percentage of variance explained by the first $j$ principal components against the number of principal components $j$.\n", "You should use the attribute `explained_variance_ratio_` from your `PCA()`." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAH3xJREFUeJzt3Xt81Hed7/HXZ3IPSYBAoIGES1vaQu2VlFKv9dbSHm2tun2A9a7luEfU3a1H6/Fst6fu0fV+dLfq0lq17lqsdq3o4tau21qtcgmVIpdCA+GShkuaAEkIuUzmc/6YgU7DhAwwyW/mN+/n4zGP+V2+M/PJL8Obb77znd/P3B0REQmXSNAFiIhI5incRURCSOEuIhJCCncRkRBSuIuIhJDCXUQkhBTuIiIhpHAXEQkhhbuISAgVBvXCkydP9lmzZgX18iIiOWn9+vUvunvNSO0CC/dZs2bR2NgY1MuLiOQkM9udTjsNy4iIhJDCXUQkhBTuIiIhpHAXEQkhhbuISAiNGO5m9oCZHTSzTcPsNzP7ppk1mdlGM7sy82WKiMjpSKfn/n1g0Sn23wDMSdyWAt8++7JERORsjDjP3d2fMrNZp2hyM/Cgx6/Xt9rMJphZrbvvy1CNIpIn3J1ozOmPxuiLxug/fhscZGDQGYzFb9GYE3MnOpi4jzmDsRiDMU7cR2OxE+0HY86gO7GY44B7/LViTmI9frnRmHt8X6JNLLHdE9tjQ65KGn+247UP3feyH+xl+944dyqX1U/IyDEbTia+xDQd2Ju03pLYdlK4m9lS4r17ZsyYkYGXFpGxFIs5PQODHO2L0t0Xpbs3+tJy3/HlwZO29fQPJkI6ObBjiRAffCnIB2MnhWRYmL20PKWqNCfC3VJsS/nrcfflwHKAhoaGkP4KRbKTu9PdF6WrN0pn70D8/tjAy5aP7+s8lrjvjdLVO8DRvihH+wY52h9NK3wjBuNKCqksKWRcSSFlxQUUF0QoLYpQVVpIcWGE4sL4tuLCCCWJW3Fh5MS2+PaCE8tFEaMgYhQWGBEzCiMRIhEojEQoOL4vkthXEF8vMDuxryCxzyweWi8txzdEDMwMIx7EkUQaH2+TvD9ZcmibpYrDYGQi3FuA+qT1OqA1A88rIqfg7nT1RTnY2UdbVx9t3X0c7OylrbuPts74ent3/4nw7uodOGlYYajSogiVpUVUlRZSVVbE+LIi6iaUUZEI6YrSQipKCuLLidu4IfcVJYWUFkWyKujyUSbCfSWwzMxWAFcDRzTeLnLm3J32o/20Hj5GW1cfB7sS4d3Vx8Gu3qQg76MvGjvp8cUFEWoqS6ipLKF2fCkXnVNJZSKsq0qLTixXlhZSVVp0YrmytJCSwoIAfmIZDSOGu5k9BFwLTDazFuDvgCIAd/8OsAq4EWgCeoAPjFaxImHg7nQei7L3UA8th3rY23EssXyMvR3x+2MDgyc9bkJ5ETUVJUypKmH+jIlMqSqlpiIe4lMqj9+XUlVWqF6zpDVbZskI+x34aMYqEgmB7r7oieAeGuAtHT109UVf1r6ypJC66nJmTx7Ha+bUUF9dxrQJZUytKqWmsoTJFcXqVctpCeyUvyK5zN1p6+5jT3sPu9t72N3Rw572o+xq72FPRw8dR/tf1r6sqID66jLqJpazYNZE6iaWn1ivn1jO+PKigH4SCSuFu8gwBgZjtB4+9rLw3p0I7z0dPfT0vzR0EjGoHV/GzEnlXH/xVOqr46Edvy+jelyxhkpkTCncJe/1R2PsaOtm675OtrR2su1AF7vbe3jh8DEGk6aXFBdGmFFdzszqcq45bxIzq8uZOWkcMyeVUzexnOJCnapJsofCXfLKoaP98RBP3Lbu66LpYBcDg/EQLy6MMGdKBZfWjeetl9Uys3ocMyaVM3NSOVMrS4lE1PuW3KBwl1AajDm72o+ydV/niR751n1d7O/sPdGmprKEubVVvPaCycyrrWJubRXnTh5HYYF64JL7FO4SCvuP9LKmuZ01zR1sbu1k+/6uE9MJCyLG+TUVLDy3mrmJEJ9bW0VNZUnAVYuMHoW75KTjYb56Zzurd3bQ/OJRID6l8OLpVSxeUM/c2irm1VZx/pQKSos0jVDyi8JdcsKwYV5ayNWzq7nt6hksPHcSc2urKNC4uIjCXbKTwlzk7CjcJStEB2M8tvkAv29qU5iLZIDCXQI1GHN+8Wwr3/jN8zS/eFRhLpIhCncJRCzm/GrTfr7+n9tpOtjN3Noqlr9nPm+cO1VhLpIBCncZU+7O41sO8LXHt/Pc/i7mTKngW7ddyaKLz9EXhEQySOEuY8LdeXJ7G19/fDsbW44we/I4vrH4ct5y6TT11EVGgcJdRpW783RTO197fBvP7DlMfXUZX37npdxyxXR9E1RkFCncZdSsbe7gq7/exprmDmrHl/L5Wy7hnfPrdIItkTGgcJeMe2bPIb7++HZ+9/yL1FSWcPdb57F4wQx9S1RkDCncJWP+3HKErz2+jSe2tVE9rpjP3jiXdy+cSVmxQl1krCnc5awd7Ozl7l9sZtWf9zO+rIhPLbqQ910zi3ElenuJBEX/+uSs/PvGfXz20T9zrH+QT7xxDh96zWyqSnXJOJGgKdzljBzpGeCulZv4+YZWLqsbz1dvvZzzp1QEXZaIJCjc5bQ9tb2NT/10Iy929/HXb7qAj77+PE1rFMkyCndJW09/lC+seo4frt7N+VMquO+9DVxSNz7oskQkhbTC3cwWAd8ACoD73f0fhuyfCTwA1AAdwLvdvSXDtUqA1u8+xB0Pb2B3Rw8ffvVsPnn9hZraKJLFRgx3MysA7gXeDLQA68xspbtvSWr2FeBBd/+Bmb0B+ALwntEoWMZWfzTGN36znW8/uYPa8WX86MMLuea8SUGXJSIjSKfnvgBocvedAGa2ArgZSA73ecBfJ5afAB7NZJESjOf2d/I3P36WLfs6+Yv5ddz11nlUaiaMSE5IJ9ynA3uT1luAq4e0eRZ4B/Ghm1uASjOb5O7tyY3MbCmwFGDGjBlnWrOMssGYc//vdvLVX2+nqqyQ+97bwJvnTQ26LBE5DemEe6pT9vmQ9U8C/2Rm7weeAl4Aoic9yH05sBygoaFh6HNIFtjT3sMdP9nAul2HuP7iqXz+lkuYVFESdFkicprSCfcWoD5pvQ5oTW7g7q3A2wHMrAJ4h7sfyVSRMvrcnRXr9vK5X26hwIyv3XoZt1wxHTOdjlckF6UT7uuAOWY2m3iPfDHwruQGZjYZ6HD3GPAZ4jNnJEcc7Ozl049s5Iltbbzq/El86Z2XMX1CWdBlichZGDHc3T1qZsuAx4hPhXzA3Teb2T1Ao7uvBK4FvmBmTnxY5qOjWLNk0LpdHXzkh+vp7oty91vn8d5rZumKSCIhYO7BDH03NDR4Y2NjIK8tcSvW7uFvf76J+onlLH/vfM6fUhl0SSIyAjNb7+4NI7XTN1TzUHQwxt//+1a+/4ddvGbOZP5pyZWML9cUR5EwUbjnmSM9Ayx76Bl+9/yLfOjVs/nMDRfpvDAiIaRwzyNNB7u5/cFGWg718KV3XMqtV9WP/CARyUkK9zzx5LaDfOyhP1FSGOGh2xfSMKs66JJEZBQp3EPO3fnu75v5/KqtXHROFfe9r0HTHEXygMI9xPqig3z2Z5v46foWbnjFOXz11ssoL9avXCQf6F96SB3s6uUjP1zPM3sO81dvmsPH3zBH89dF8ojCPYQ2vXCE2x9s5HDPAN+67UpuvKQ26JJEZIwp3EPm3zfu446fbKC6vJif/uU1XDxNV0oSyUcK95CIxZz/95vn+eZvnmf+zIl8593zqanU2RxF8pXCPQSO9kW54+Fn+Y/N+7m1oY7Pve0VlBTqEngi+UzhnuNaDvXw4R80sv1AF3/7lnl88FWzdJpeEVG457Idbd0sWb6aYwODfO8DC3jdBTVBlyQiWULhnqOaDnaz5L7VuDuP/OUruWCqzugoIi9RuOegl4IdHrp9IXMU7CIyhMI9xzQd7GLx8jUArFh6tc7BLiIp6VyvOeT5Awp2EUmPwj1HPH+giyX3rcYMVixdqGAXkVNSuOeA7SeC3Xjo9oWcP6Ui6JJEJMsp3LPc9gNdLFm+mohZoseuYBeRkSncs9i2/fFgL4gYDy1dyHk1CnYRSY9my2Spbfu7eNd9qyksiA/FnKtgF5HToHDPQs/t7+Rd962hqMBYsfQaZk8eF3RJIpJj0hqWMbNFZrbNzJrM7M4U+2eY2RNm9icz22hmN2a+1PywdV882IsLIgp2ETljI4a7mRUA9wI3APOAJWY2b0iz/w087O5XAIuBb2W60HwQD/bVlBRGWLF0oYJdRM5YOj33BUCTu+90935gBXDzkDYOVCWWxwOtmSsxP2xpjQd7aVEBK5YuZJaCXUTOQjrhPh3Ym7TektiW7G7g3WbWAqwCPpbqicxsqZk1mlljW1vbGZQbTltaO7nt/tWUJYJ95iQFu4icnXTCPdXJwX3I+hLg++5eB9wI/NDMTnpud1/u7g3u3lBTo9PTAmxuPcK7TgT7NQp2EcmIdMK9BahPWq/j5GGXDwEPA7j7H4FSYHImCgyzza1HuO3+NYwrLmTF0muYMak86JJEJCTSCfd1wBwzm21mxcQ/MF05pM0e4I0AZjaXeLhr3OUUjvZFWfrgesqLCnjo9oUKdhHJqBHD3d2jwDLgMWAr8Vkxm83sHjO7KdHsDuB2M3sWeAh4v7sPHbqRJF/59TZajxzjH991hYJdRDIurS8xufsq4h+UJm+7K2l5C/CqzJYWXs/sOcT3/7CL9yycyfyZ1UGXIyIhpHPLjLH+aIw7H9lIbVUpn1p0UdDliEhI6fQDY+xbTzax/UA333v/VVSU6PCLyOhQz30MbT/Qxb1PNHHz5dN4/UVTgi5HREJM4T5GBmPOpx/ZSEVJIXe9ZejZG0REMkvhPkYe/OMu/rTnMH/31ouZVFESdDkiEnIK9zHQcqiHLz+2jWsvrOHmy6cFXY6I5AGF+yhzd/7XzzYB8H9vuQSzVGdzEBHJLIX7KHt0wws8tb2NT11/IdMnlAVdjojkCYX7KGrv7uOeX2zhyhkTeM81s4IuR0TyiMJ9FP2fX2zhaN8gX3zHpRRENBwjImNH4T5K/uu5A6x8tpWPvv585kytDLocEckzCvdR0NU7wGd/tokLp1byl9eeF3Q5IpKH9P33UfCl/9jG/s5evnXblRQX6v9PERl7Sp4MW7ergx+u3s0HXjmbK2ZMDLocEclTCvcM6h0Y5NOPbGT6hDLuuO6CoMsRkTymYZkM+qf/amJn21Ee/OACxumMjyISIPXcM2Trvk6+89sdvP3K6bz2Al38W0SCpXDPgMGYc+cjGxlfVsTf/jed8VFEgqexgwz43tPNPNtyhH9ccgUTxxUHXY6IiHruZ2tPew9f+fU23jR3Cm+5tDbockREAIX7WYmf8fHPFEYifO5tr9AZH0Ukayjcz8JP1rfw+6YXufOGi6gdrzM+ikj2ULifoYNdvfz9L7ewYFY171owI+hyREReJq1wN7NFZrbNzJrM7M4U+79uZhsSt+1mdjjzpWaXf/jVc/RGY3zhHZcQ0RkfRSTLjDhbxswKgHuBNwMtwDozW+nuW463cfe/Tmr/MeCKUag1a+xs6+bRP73Ah149m/NqKoIuR0TkJOn03BcATe6+0937gRXAzadovwR4KBPFZat7n9hBcWGEpa/VGR9FJDulE+7Tgb1J6y2JbScxs5nAbOC/htm/1Mwazayxra3tdGvNCrvbj/Lohhe47eqZ1FSWBF2OiEhK6YR7qgFlH6btYuCn7j6Yaqe7L3f3BndvqKnJza/o3/tEE4UR47+/9tygSxERGVY64d4C1Cet1wGtw7RdTIiHZPZ29PBvz7zAkgUzmFJVGnQ5IiLDSifc1wFzzGy2mRUTD/CVQxuZ2YXAROCPmS0xe3zrySYiZnzkdRprF5HsNmK4u3sUWAY8BmwFHnb3zWZ2j5ndlNR0CbDC3YcbsslpLYd6+On6FhYvqOec8eq1i0h2S+vEYe6+Clg1ZNtdQ9bvzlxZ2efbT+4AUK9dRHKCvqGahtbDx3i4cS9/0VDPtAk6zYCIZD+Fexr++bc7cIf/ca167SKSGxTuIzjQ2ctD6/byzvl11E0sD7ocEZG0KNxH8J3f7mAw5nz09ecHXYqISNoU7qdwsLOXH63Zw9uvmE59tXrtIpI7FO6nsPypnURjzrI3qNcuIrlF4T6MF7v7+Jc1u7n58mnMnDQu6HJERE6Lwn0Y9z21k/5ojGUaaxeRHKRwT6G9u48H/7ibmy6bxrk6X7uI5CCFewr3/76Z3uigxtpFJGcp3Ic4dLSfB/+wi7dcOo3zp1QGXY6IyBlRuA/xwNPNHO0f5GPqtYtIDlO4JznSM8D3n97FjZecwwVT1WsXkdylcE/ywNPNdPVF+dgb5gRdiojIWVG4Jxw5NsADTzdz/cVTmVtbFXQ5IiJnReGe8IM/7KKrN8rH36heu4jkPoU70NU7wHd/38yb5k7l4mnjgy5HROSsKdyBB/+4myPHBviEeu0iEhJ5H+7dfVHu+91O3nDRFC6pU69dRMIh78P9X1bv5nDPgMbaRSRU8jrce/qj3PfUTl53QQ2X108IuhwRkYzJ63D/19V7aD/ar167iIRO3ob7sf5B/vmpHbxmzmTmz5wYdDkiIhmVVrib2SIz22ZmTWZ25zBtbjWzLWa22cx+lNkyM+9Ha/fwYrd67SISToUjNTCzAuBe4M1AC7DOzFa6+5akNnOAzwCvcvdDZjZltArOhN6BQb7z2x1cc+4krppVHXQ5IiIZl07PfQHQ5O473b0fWAHcPKTN7cC97n4IwN0PZrbMzPrJ+hbauvr4xJvUaxeRcEon3KcDe5PWWxLbkl0AXGBmT5vZajNblOqJzGypmTWaWWNbW9uZVZwBv9l6gPNqxrHw3EmB1SAiMprSCXdLsc2HrBcCc4BrgSXA/WZ20txCd1/u7g3u3lBTU3O6tWbEYMxp3HVIwS4ioZZOuLcA9UnrdUBrijY/d/cBd28GthEP+6yzdV8n3X1RFszWWLuIhFc64b4OmGNms82sGFgMrBzS5lHg9QBmNpn4MM3OTBaaKWuaOwAU7iISaiOGu7tHgWXAY8BW4GF332xm95jZTYlmjwHtZrYFeAL4n+7ePlpFn411zR3UV5dRO74s6FJEREbNiFMhAdx9FbBqyLa7kpYd+JvELWu5O2t3dfD6C7N6pqaIyFnLq2+o7mjrpuNoP1drSEZEQi6vwl3j7SKSL/Iq3Nc2dzClsoSZk8qDLkVEZFTlTbi7O2ubO1gwuxqzVFP3RUTCI2/CveXQMfYd6dWQjIjkhbwJ97UabxeRPJJX4T6+rIgLplQGXYqIyKjLn3Df1cFVs6qJRDTeLiLhlxfhfrCzl+YXj2p+u4jkjbwI97W7NN4uIvklP8K9uYPy4gIunlYVdCkiImMib8J9/syJFBbkxY8rIhL+cD/c08+2A10s0LVSRSSPhD7cG3cdwl3j7SKSX0If7mt3dVBcEOGy+pOu+iciElqhD/c1zR1cXj+B0qKCoEsRERkzoQ73o31RNr1whKtmTwy6FBGRMRXqcH9mzyEGY86C2ZOCLkVEZEyFOtzXNXcQMZg/Uz13EckvoQ73Nc0dvGL6eCpK0rpUrIhIaIQ23Puig/xp72HNbxeRvBTacN/YcoT+aEzz20UkL4U23I9fnOMq9dxFJA+lFe5mtsjMtplZk5ndmWL/+82szcw2JG4fznypp2dNcwcXTK1g4rjioEsRERlzI37SaGYFwL3Am4EWYJ2ZrXT3LUOa/tjdl41CjactOhjjmd2HeNsV04IuRUQkEOn03BcATe6+0937gRXAzaNb1tnZuq+L7r6o5reLSN5KJ9ynA3uT1lsS24Z6h5ltNLOfmll9qicys6Vm1mhmjW1tbWdQbnrWNLcDaKaMiOStdMI91UVHfcj6L4BZ7n4p8J/AD1I9kbsvd/cGd2+oqak5vUpPw9rmDmZOKuec8aWj9hoiItksnXBvAZJ74nVAa3IDd293977E6n3A/MyUd/piMWdd4mLYIiL5Kp1wXwfMMbPZZlYMLAZWJjcws9qk1ZuArZkr8fQ0tXVzqGdA89tFJK+NOFvG3aNmtgx4DCgAHnD3zWZ2D9Do7iuBj5vZTUAU6ADeP4o1n9KaxPz2qxXuIpLH0jrpiruvAlYN2XZX0vJngM9ktrQzs665g6lVJcyoLg+6FBGRwITqG6ruztrmDhbMnoRZqs+BRUTyQ6jCfW/HMfZ39mq8XUTyXqjCXfPbRUTiQhXua5s7mFBexJwpFUGXIiISqHCFe2J+eySi8XYRyW+hCfcDnb3sbu/RFEgREUIU7sfP364PU0VEQhbu44oLmFdbFXQpIiKBC1W4z59VTWFBaH4kEZEzFookPHS0n20Hulgwa2LQpYiIZIVQhPu6XcfH23VxDhERCEm4r23uoLgwwqV144MuRUQkK4Qi3Nft6uDy+gmUFhUEXYqISFbI+XDv7ouyqbVT89tFRJLkfLg/s/sQgzHX/HYRkSQ5H+5rmzsoiBhXztBMGRGR40IR7q+YVsW4krSuOyIikhdyOtx7BwbZsPewhmRERIbI6XB/du9h+gdjmt8uIjJETof78S8vXaVvpoqIvExOh/ua5g4uOqeSCeXFQZciIpJVcjbco4Mx1u8+xFW6pJ6IyElyNtw3t3bS0z+oD1NFRFJIK9zNbJGZbTOzJjO78xTt3mlmbmYNmSsxNV2cQ0RkeCOGu5kVAPcCNwDzgCVmNi9Fu0rg48CaTBeZyprmDmZNKmdqVelYvJyISE5Jp+e+AGhy953u3g+sAG5O0e5zwJeA3gzWl1Is5qzb1aFeu4jIMNIJ9+nA3qT1lsS2E8zsCqDe3X+ZwdqG9fzBbo4cG9D8dhGRYaQT7pZim5/YaRYBvg7cMeITmS01s0Yza2xra0u/yiHWNrcDsEAzZUREUkon3FuA+qT1OqA1ab0SeAXwpJntAhYCK1N9qOruy929wd0bampqzrjoNc0dnFNVSn112Rk/h4hImKUT7uuAOWY228yKgcXAyuM73f2Iu09291nuPgtYDdzk7o2jUbC7s7Y5Pt5uluqPChERGTHc3T0KLAMeA7YCD7v7ZjO7x8xuGu0Ch9rd3sPBrj59mCoicgppnSfX3VcBq4Zsu2uYtteefVnDOz6/XVdeEhEZXs59Q3VCeRFvnjeV86dUBF2KiEjWyrkrXFx38Tlcd/E5QZchIpLVcq7nLiIiI1O4i4iEkMJdRCSEFO4iIiGkcBcRCSGFu4hICCncRURCSOEuIhJC5u4jtxqNFzZrA3YH8uLpmwy8GHQRaVCdmZUrdULu1Ko6M2emu494Wt3Awj0XmFmju4/69WDPlurMrFypE3KnVtU59jQsIyISQgp3EZEQUrif2vKgC0iT6sysXKkTcqdW1TnGNOYuIhJC6rmLiIRQXoe7mdWb2RNmttXMNpvZJ1K0udbMjpjZhsQt5RWoxoKZ7TKzPyfqOOkatRb3TTNrMrONZnZlADVemHSsNphZp5n91ZA2gR1TM3vAzA6a2aakbdVm9riZPZ+4nzjMY9+XaPO8mb0vgDq/bGbPJX63PzOzCcM89pTvkzGo824zeyHp93vjMI9dZGbbEu/XOwOo88dJNe4ysw3DPHbMjmdGuXve3oBa4MrEciWwHZg3pM21wC+DrjVRyy5g8in23wj8CjBgIbAm4HoLgP3E5+VmxTEFXgtcCWxK2vYl4M7E8p3AF1M8rhrYmbifmFieOMZ1XgcUJpa/mKrOdN4nY1Dn3cAn03hv7ADOBYqBZ4f+2xvtOofs/ypwV9DHM5O3vO65u/s+d38msdxF/ALg04Ot6qzcDDzocauBCWZWG2A9bwR2uHvWfFnN3Z8COoZsvhn4QWL5B8DbUjz0euBxd+9w90PA48CisazT3X/t8QvWA6wG6kbr9dM1zPFMxwKgyd13uns/sIL472FUnKpOMzPgVuCh0Xr9IOR1uCczs1nAFcCaFLuvMbNnzexXZnbxmBb2cg782szWm9nSFPunA3uT1lsI9j+rxQz/DyZbjinAVHffB/H/8IEpKdpk27H9IPG/0lIZ6X0yFpYlho8eGGaYK5uO52uAA+7+/DD7s+F4njaFO2BmFcAjwF+5e+eQ3c8QH1a4DPhH4NGxri/Jq9z9SuAG4KNm9toh+y3FYwKZDmVmxcBNwE9S7M6mY5qubDq2nwWiwL8O02Sk98lo+zZwHnA5sI/4kMdQWXM8gSWcutce9PE8I3kf7mZWRDzY/9Xd/23ofnfvdPfuxPIqoMjMJo9xmcdraU3cHwR+RvxP22QtQH3Seh3QOjbVneQG4Bl3PzB0RzYd04QDx4evEvcHU7TJimOb+CD3LcBtnhgQHiqN98mocvcD7j7o7jHgvmFeP1uOZyHwduDHw7UJ+nieqbwO98RY23eBre7+tWHanJNoh5ktIH7M2seuyhN1jDOzyuPLxD9c2zSk2UrgvYlZMwuBI8eHGwIwbG8oW45pkpXA8dkv7wN+nqLNY8B1ZjYxMcxwXWLbmDGzRcCngZvcvWeYNum8T0bVkM95bhnm9dcBc8xsduKvvMXEfw9j7U3Ac+7ekmpnNhzPMxb0J7pB3oBXE/9TcCOwIXG7EfgI8JFEm2XAZuKf5q8GXhlQrecmang2Uc9nE9uTazXgXuKzEP4MNARUaznxsB6ftC0rjinx/3D2AQPEe48fAiYBvwGeT9xXJ9o2APcnPfaDQFPi9oEA6mwiPk59/L36nUTbacCqU71PxrjOHybefxuJB3bt0DoT6zcSn6G2I4g6E9u/f/x9mdQ2sOOZyZu+oSoiEkJ5PSwjIhJWCncRkRBSuIuIhJDCXUQkhBTuIiIhpHAXEQkhhbuISAgp3EVEQuj/A+/LwWel8O9ZAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plt.plot(nrange,pca.explained_variance_ratio_.cumsum())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Part 3 - Partial Least Squares" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We want to apply partial least squares regression to this data set.\n", "The function `PLSRegresssion` is provided by sklearn in the module `cross_decomposition`.\n", "\n", "**Task**: Implement a loop over the number of components in your PLS regression model.\n", "Use 10-fold cross-validation and store the means of the mean squared errors in a list called `mse`.\n", "You can use an appropriate `scoring` option." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "122318.40415029772\n", "120869.20847675091\n", "120342.57274839829\n", "119462.09544296702\n", "122189.40601362588\n", "120883.69908106592\n", "119227.66421886161\n", "115685.97332435506\n", "115984.97109227064\n", "115717.44911488841\n", "115944.5092787043\n", "115798.45521840893\n", "116644.281406398\n", "114824.45034020235\n", "115712.18085901762\n", "116097.3041963903\n", "115564.97080748554\n", "115676.89512608883\n", "116599.01367380246\n" ] } ], "source": [ "from sklearn.cross_decomposition import PLSRegression\n", "\n", "mse = []\n", "mseminus = []\n", "mseplus= []\n", "for i in range(N):\n", " pls = PLSRegression(n_components=i+1)\n", " cvs = -1*cross_val_score(pls, X, y, cv=10, scoring='neg_mean_squared_error')\n", " mse.append(cvs.mean())\n", " mseminus.append(cvs.mean() - cvs.std())\n", " mseplus.append(cvs.mean() + cvs.std())\n", " print(cvs.mean())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Task**: Determine the number of components, for which the MSE is minimized." ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "scrolled": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Minimum MSE for 14 components\n" ] } ], "source": [ "print('Minimum MSE for %d components' % (np.argmin(mse)+1))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Task**: Plot the MSE against the number of components in your model." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAD8CAYAAACLrvgBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvqOYd8AAAIABJREFUeJzt3XtclNed+PHPl6siICAgF1FQUbyjotHEpEmsjblUcjPG5lfdxNamzaVttt0k20u623Y33ew23dw3t0bTxEvTJJqYm7FJTIwaUVHxjgqKCKiAgIjczu+PebAjARlwZp4Z+L5fr3nNzHnO88x3xpHvPOec5xwxxqCUUkq5IsDuAJRSSvkPTRpKKaVcpklDKaWUyzRpKKWUcpkmDaWUUi7TpKGUUsplmjSUUkq5TJOGUkopl2nSUEop5bIguwNwt9jYWJOammp3GEop5Vc2b958whgT11G9bpc0UlNTycnJsTsMpZTyKyJS6Eo9bZ5SSinlMk0aSimlXKZJQymllMs0aSillHKZJg2llFIu06ShlFLKZZo0lFJKuUyThmV7USV/+GAPuvytUkq1T5OGZduRSp799ABbDlfYHYpSSvksTRqWWyYOoG/vYF764pDdoSillM/SpGEJCwli7uSBfJBXwpHyWrvDsV1jUzONTc12h6GU8jGaNJzMv3QQASK88mWB3aHY7o4XN/LjZbl2h6GU8jGaNJwk9u3N9WMTWbbpCNV1DXaHY5uqugY2FZSzavsxdhafsjscpZQP6TBpiMjLIlImInlOZY+JyB4R2S4ib4lIlFU+Q0Q2i8gO6/5qp30mWuX5IvKEiIhVHiMiq0Vkv3UfbZWLVS/fep0J7n/7X7dgWho1ZxtZnlPkjZfzSVsPV9JsDSJ76u/59gajlPIprpxpvALMbFW2GhhtjBkL7AMetspPAN82xowB5gOvOu3zLLAQSLduLcd8CFhjjEkH1ljPAa51qrvQ2t/jxg6IYlJqNH9ed4im5p45/HbToXICA4QF09J4P6+EvSXVdoeklPIRHSYNY8xaoLxV2UfGmEbr6QZggFW+1RhTbJXvBHqJSKiIJAKRxpj1xnEhxGLgRqteNrDIeryoVfli47ABiLKO43ELpqVRVHGGj3aWeOPlfM6mgnJGJUVy39VD6RMSyJN/3293SEopH+GOPo27gPfbKL8F2GqMOQskA87tPUVWGUB/Y8wxAOs+3ipPBo60s49HzRiZQEpM7x45/La+sZncI5VkDYohKiyEeZemsmrHMfLLauwOTSnlAy4qaYjIL4BG4LVW5aOAPwA/aClqY/eO2n5c3kdEFopIjojkHD9+vIPDdiwwQLjz0jRyCivYdqTyoo/nT/KKT3G2sZlJqdEAfG9aGr2CAnn6E+3bUEpdRNIQkfnADcAdxmnuDREZALwFzDPGHLCKi7CasCwDgJZmrNKWZifrvsxpn5R29jmPMeZ5Y0yWMSYrLq7DJW5dctukFCJCg3rc2UZOgaMlMis1BoB+4aF8d+ogVuQe5dCJ03aGppTyAV1KGiIyE3gQmGWMqXUqjwJWAQ8bY9a1lFvNTtUiMsUaNTUPWGFtXomj0xzr3rl8njWKagpwqqUZyxvCQ4OYMymFVTuOUVx5xlsva7tNBRWkxfYhLiL0XNn3Lk8jODCAZ/RsQ6kez5Uht0uA9cBwESkSkQXAU0AEsFpEckXkOav6vcBQ4FdWea6ItPRR/BB4EcgHDvCPfpBHgRkish+YYT0HeA84aNV/AfjRxb3Vzpt/aSrGGBatL/D2S9uiudmQU1BO1qDo88rjI3rxnUsG8ubWo3q1vFI9XFBHFYwxc9sofqmdur8DftfOthxgdBvlJ4HpbZQb4J6O4vOklJgwrh2dyJKNh7n/6nT6hHb4cfm1gydqqKhtYJLVNOXs7m8M4bWNh3nm03z+8+axNkSnlPIFekV4B+6alkZVXSN/29L9L/bbVOCY4TcrNfpr2/pH9mJOVgpvbC7iaA9qrlNKnU+TRgcmDoomMyWKl784RHM3v9hvU0E5seEhpMX2aXP73VcOAeC5Tw+0uV0p1f1p0nDBgmlpFJysZc2eso4r+7GcggqyBsVgzfDyNclRvbl1YgrLNh2h5FSdl6NTSvkCTRouuHZ0Akl9e/HSFwftDsVjSqvqOFxe22bTlLMfXTmEZmN47jM921CqJ9Kk4YKgwAD+6bJUNhwsJ+9o95z1dZN1fUZbneDOUmLCuGl8Mku+OkxZtZ5tKNXTaNJw0ZxJAwkLCeTlbnqxX05BBb2DAxmZFNlh3XuuGkpDUzMvrO2+Z15KqbZp0nBR397B3JaVwjvbiymr6n6/sDcVlDN+YBTBgR1/JVJj+3BjZjJ/2XCYkzVnvRCdUspXaNLohDsvS6Wx2bB4faHdobhVdV0Du49Vddg05exHVw2lrrGJF7vpmZdSqm2aNDphUL8+zBjRn9c2FlLX0GR3OG7TsuhSZ5LG0PhwbhibxOIvC6g4Xe/B6JRSvkSTRictmJZGRW0Db245ancobpNT4Fh0KXNgVKf2u+/qoZyub+LldXq2oVRPoUmjkyanxTA6OZKXvjjYbS72+6qgnJGJkYR3cpqUYf0juHZ0Aq+sK+DUmZ67prpSPYkmjU4ScSyDeuD4aT7bf/Frd9jt3KJLHVyf0Z57rx5K9dlGXllX4N7AlFI+SZNGF1w/Jon4iNBuMfx2Z/Ep6hqaO9Wf4WxUUl9mjOzPS18cpLpOzzaU6u40aXRBSFAA8y9N5fP9J9hbUm13OBcl5wKTFLrq/qvTqapr7HajypRSX6dJo4vuuGQgvYID/H5qkU0F5aT2CyM+oleXjzFmQF+uGh7Hi58f5PTZRjdGp5TyNZo0uigqLIRbJgzg7dxiTvjpBW7GGHIKK84t7Xox7pueTkVtA3/ZoGcbSnVnmjQuwl3T0qhvbPbbP5QHjp+m/HQ9ky6iaarFhIHRXJ4ey/NrD3Kmvvtcw6KUOp8ry72+LCJlIpLnVPaYiOwRke0i8pa1Njgi0k9EPhGRGhF5qtVxJorIDhHJF5EnrLXCEZEYEVktIvut+2irXKx6+dbrTHDvW794Q+LCuTojnr9s8M+L/XKsSQrdcaYBcP/0dE6erue1jf6ZRJVSHXPlTOMVYGarstXAaGPMWGAf8LBVXgf8CvhZG8d5FlgIpFu3lmM+BKwxxqQDa6znANc61V1o7e9zFkxL40RNPStzi+0OpdM2FVTQr08Ig9tZdKmzJqXGMHVwP/5v7UG/TKJKqY51mDSMMWuB8lZlHxljWno8NwADrPLTxpgvcCSPc0QkEYg0xqy31v5eDNxobc4GFlmPF7UqX2wcNgBR1nF8yqVD+pGREMHL6w7heGv+I6ewnKzU6HYXXeqK+6enc7z6LMs2HXHbMZVSvsMdfRp3Ae93UCcZcF5ku8gqA+hvjDkGYN3HO+1zpJ19fIaIcNe0NPaUVLMu/6Td4bisrKqOwpO1Xb4+oz1TBscwKTWaZz89wNlGPdtQqru5qKQhIr8AGoHXOqraRllHP8td3kdEFopIjojkHD/u/au0Z41LIjY8xK+G3246d32Ge5OGiHD/9HRKqup4Y3NRxzsopfxKl5OGiMwHbgDuMB23yxRhNWFZBgAtnQClLc1O1n2Z0z4p7exzHmPM88aYLGNMVlxcXOfeiBv0Cg7ku1NS+WTvcfLL/ONiv00F5fQKDmCUC4sudda0obGMHxjFM58coL6x2e3HV0rZp0tJQ0RmAg8Cs4wxtR3Vt5qdqkVkijVqah6wwtq8EphvPZ7fqnyeNYpqCnCqpRnLF90xZSAhQQG87CdzMOUUljM+JdqlRZc6S0S4/+p0jlae4a2terahVHfiypDbJcB6YLiIFInIAuApIAJYLSK5IvKcU/0C4I/AP1n1R1qbfgi8COQDB/hHP8ijwAwR2Q/MsJ4DvAcctOq/APzoYt6op8WGh3LLhGSWfnWYl77w7U7xmrON7Cqucsv1Ge25cngcGQkROrWIUt1Mh3NhG2PmtlH80gXqp7ZTngOMbqP8JDC9jXID3NNRfL7kl9eP5GRNPb99dxd7S6r47Y2jCQ0KtDusr9l6uMKx6FKae/sznIkIcycP5JGVO8k7eorRyX099lpKKe/RK8LdqE9oEM/9v4ncf/VQlucUcccLGzle7XtTjGwqqCBAYPxAz51pANyYmUxIUADLc3T4rVLdhSYNNwsIEB741nCe+s548opPkf3UF+QdPWV3WOfZdKickUmdX3Sps/qGBTNzVAJvbz2qF/sp1U1o0vCQG8Ym8cbdl2KAW5/7klXbfaMPv6Gpma1HKsga5LmmKWdzJqVQVdfIhztLvPJ6SinP0qThQaOT+7Ly3mmMSurLPa9v4Y+r99m+ROzO4qqLWnSps6YO7kdKTG+9QlypbkKThofFRYTy+vcvYfbEATyxZj8/fG2zrWtOtExS6MmRU84CAoTZE1P48sBJDp/scHS2UsrHadLwgtCgQP7r1rH86oaRrN5Vyi3PfsmRcnv+gG4qKGdQvzDiI7u+6FJn3TpxAAECf92sZxtK+TtNGl4iIiyYlsYrd06muPIM2U+vY+NB785VZYwhp8B7/RktkqJ6c8WwON7YXESTzc1zSqmLo0nDy64YFsfb91xGVFgwd7y4kdc3Hvbaax88cZqTblp0qbPmZKVw7FQda/d7f24wpZT7aNKwweC4cN6+5zKmpcfyr2/t4JEVeTQ0eX6OJncvutQZ00f0p1+fEJZ9pU1USvkzTRo2iewVzEvzJ7HwisEsWl/I/Je/ouJ0vUdfc1NBBTF9QhgS555FlzojJCiAm8Yn8/HuUr9dU10ppUnDVoEBwr9eN4L/mT2OnIIKbnxmHftKPTdLbk5BOVmD3LvoUmfMmZRCY7PhrS1HbXl9pdTF06ThA26ZOIClP5hCbX0TNz/zJX/fU+r21yirrqPAA4sudUZ6/wjGD4xiWc4Rn57QUSnVPk0aPmLCwGjeuXcaabF9uPsvW9hTUuXW428+t+iS9zvBnc3JSiG/rIYthyttjUMp1TWaNHxIQt9e/PnOSUT2Cua+17dypt598zV9dW7RJXtnm71hXBJhIYEs1yvElfJLmjR8TGx4KI/PGcf+shp+u2qX246bU1BBZkoUIUH2/pOHhwZx/ZhE3t1ebOuV8UqprtGk4YMuT4/j7m8M4fWNh3l/x8VPdFhztpGdxads7c9wdvvkFE7XN/nMJI5KKddp0vBR//ytYYxLieLBv22nqOLiphzJPVzpWHTJR5LGhIHRDInrwzJdZ0Mpv6NJw0cFBwbw5O3jaTbwk6W5NF7ExX+bCsqtRZei3Bhh14kIcyalsLmwgvwyzw0xVkq5nytrhL8sImUikudU9piI7BGR7SLylohEOW17WETyRWSviFzjVD7TKssXkYecytNEZKOI7BeRZSISYpWHWs/zre2p7nrT/mJgvzB+f9NocgoreGLN/i4fJ6ewnBGJkUT0CnZjdBfn5gkDCAoQlucU2R2KUqoTXDnTeAWY2apsNTDaGDMW2Ac8DCAiI4HbgVHWPs+ISKCIBAJPA9cCI4G5Vl2APwCPG2PSgQpggVW+AKgwxgwFHrfq9TjZmcncOnEAT36Sz/oDnZ/gsKGpmS2FlT7TNNUiNjyU6SPieXNLkVemUFFKuUeHScMYsxYob1X2kTGmZejLBmCA9TgbWGqMOWuMOQTkA5OtW74x5qAxph5YCmSL49Lkq4E3rP0XATc6HWuR9fgNYLrYdSmzzf5t1ijS+vXhp8tyOz3VyK7iKs40NNl+fUZb5kxK4URNPWt2l9kdilLKRe7o07gLeN96nAw4924WWWXtlfcDKp0SUEv5eceytp+y6n+NiCwUkRwRyTl+vPvNotonNIgn5o6n/HQ9P39je6eupt50btEl3zrTALgiPY7+kaEs1w5xpfzGRSUNEfkF0Ai81lLURjXThfILHevrhcY8b4zJMsZkxcXFXThoPzU6uS8PXpvBx7tLeXVDocv75RRUMDAmjP5eXHTJVUGBAdw6cQCf7i2j5FSd3eEopVzQ5aQhIvOBG4A7zD9++hYBKU7VBgDFFyg/AUSJSFCr8vOOZW3vS6tmsp7mrstSuTojnt+t2s2u4o6nGTHGkFNY7pNNUy1uy0qh2cAbuqqfUn6hS0lDRGYCDwKzjDHOFxGsBG63Rj6lAenAV8AmIN0aKRWCo7N8pZVsPgFutfafD6xwOtZ86/GtwN9ND5/lTkR47Nax9O0dzH1LtlBbf+Erqg+dOM2JmnqfbJpqMahfH6YO7sfynCKadVU/pXyeK0NulwDrgeEiUiQiC4CngAhgtYjkishzAMaYncByYBfwAXCPMabJ6pO4F/gQ2A0st+qCI/k8ICL5OPosXrLKXwL6WeUPAOeG6fZk/cJD+dOcTA6eOM1v373wNCM51iSFdqzU1xlzJqVwuLyWDYe8u/ytUqrzgjqqYIyZ20bxS22UtdT/PfD7NsrfA95ro/wgjtFVrcvrgNkdxdcTXTY0lh9+YwjPfHqAaUPjuH5sYpv1NhWUEx0WzJC4cC9H2DkzRycQsSKI5ZuOcOmQWLvDUUpdgF4R7qd+OmMY4wdG8dCb2zlS3vY0IzmFFWSlxti26JKregUHcmNmMu/nlXDqTIPd4SilLkCThp8KDgzgidvHg4EfL936tWlGjlef5dCJ0z7fNNVizqQUzjY2szJXV/VTypdp0vBjKTFh/MfNY9hyuJI/fXz+NCObCx0DzbJ8uBPc2ejkvoxMjNRJDJXycZo0/Ny3xyVxW9YAnv40ny8PnDhX/tWhCkKDAhht86JLnTFnUgp5R6vYWXzK7lCUUu3QpNEN/GbWKNJiHdOMlFvTjOQUlvvEokudcWNmMiFBAbqqn1I+zH/+oqh2hYUE8eTc8VScbuDnf91mLbpU5dPXZ7Slb1gwM0cl8HZuMXUN7lvqVinlPpo0uolRSX15+LoM1uwp42fLt9HUbJiU5l9JAxxNVKfONPDhzhK7Q1FKtUGTRjfyT5emMj0jng92lhAgMMFHFl3qjKmD+5ES01snMVTKR2nS6EZEhMdmjyM+IpRRSX19atElVwUECLdNTGFd/kkOn7y4ZW6VUu6nSaObiekTwps/upRn7phgdyhddmvWAAIE/qqTGCrlczRpdEMDosNIiQmzO4wuS+zbmyuGxfHG5iKadBJDpXyKJg3lk+ZkpXDsVB1r93e/RbWU8meaNJRPmj6iP/36hOg1G0r5GE0ayieFBAVw0/hkPt5dysmas3aHo5SyaNJQPmvOpBQamgzvbCvuuLJSyis0aSifld4/gqS+vcg9Uml3KEopiyYN5dMyEiPZU1JtdxhKKYsry72+LCJlIpLnVDZbRHaKSLOIZDmVh4jIn0Vkh4hsE5ErnbZNtMrzReQJsVYGEpEYEVktIvut+2irXKx6+SKyXUT898ID1WUZCRHkl9VQ39jccWWllMe5cqbxCjCzVVkecDOwtlX59wGMMWOAGcD/iEjLazwLLATSrVvLMR8C1hhj0oE1/GMt8Gud6i609lc9zIjESBqbDfllNXaHopTChaRhjFkLlLcq222M2dtG9ZE4/vBjjCkDKoEsEUkEIo0x640xBlgM3Gjtkw0ssh4valW+2DhsAKKs46geZERiBAB7SqpsjkQpBe7v09gGZItIkIikAROBFCAZKHKqV2SVAfQ3xhwDsO7jrfJk4Eg7+6geIrVfH0KCArRfQykfEeTm470MjABygELgS6ARkDbqdjQ/hMv7iMhCHE1YDBw40NVYlR8ICgxgeP8Idh/TMw2lfIFbzzSMMY3GmJ8aYzKNMdlAFLAfx1nCAKeqA4CWwfelLc1O1n2ZVV6E4yylrX1av+7zxpgsY0xWXFyc+96Q8gkZCRHsPqZnGkr5ArcmDREJE5E+1uMZQKMxZpfV7FQtIlOsUVPzgBXWbiuB+dbj+a3K51mjqKYAp1qasVTPkpEYyYmasxyv1ivDlbJbh81TIrIEuBKIFZEi4BEcHeNPAnHAKhHJNcZcg6M/4kMRaQaOAt91OtQPcYzE6g28b90AHgWWi8gC4DAw2yp/D7gOyAdqgTu7/C6VX3PuDI+L0DNJpezUYdIwxsxtZ9NbbdQtAIa3c5wcYHQb5SeB6W2UG+CejuJT3V9GQiQAe45Vc3m6Jg2l7KRXhCufF9MnhP6RoezWYbdK2U6ThvILIxIjtTNcKR+gSUP5hYyESPLLqmlo0ulElLKTJg3lF0YkRtDQZDh4/LTdoSjVo2nSUH5hRKKjM1wv8lPKXpo0lF9Ii+1DSGCAdoYrZTNNGsovBAcGMDQ+nD3aGa6UrTRpKL/hGEGlZxpK2UmThvIbIxIjKKs+y8kanU5EKbto0lB+49yV4TpNulK20aSh/EbLHFTaRKWUfTRpKL/RLzyUuIhQPdNQykaaNJRfcaytoWcaStlFk4byKyMTI9lfWkOjTieilC00aSi/kpEYQX1TM4dO6HQiStlBk4byKy0jqHZpE5VSttCkofzKkLhwggNFO8OVskmHSUNEXhaRMhHJcyqbLSI7RaRZRLKcyoNFZJGI7BCR3SLysNO2mSKyV0TyReQhp/I0EdkoIvtFZJmIhFjlodbzfGt7qrvetPJfIUEBDIkLZ4+eaShlC1fONF4BZrYqywNuBta2Kp8NhBpjxgATgR+ISKqIBAJPA9cCI4G5IjLS2ucPwOPGmHSgAlhglS8AKowxQ4HHrXpK6YJMStmow6RhjFkLlLcq222M2dtWdaCPiAQBvYF6oAqYDOQbYw4aY+qBpUC2iAhwNfCGtf8i4Ebrcbb1HGv7dKu+6uFGJEZQUlVHxel6u0NRqsdxd5/GG8Bp4BhwGPhvY0w5kAwccapXZJX1AyqNMY2tynHex9p+yqqvejidTkQp+7g7aUwGmoAkIA34ZxEZDLR1hmAuUE4H284jIgtFJEdEco4fP975qJVfydDpRJSyjbuTxneAD4wxDcaYMmAdkIXjDCLFqd4AoBg4AURZzVnO5TjvY23vS6tmshbGmOeNMVnGmKy4uDg3vyXla+IjehEbHsIeXZBJKa9zd9I4DFwtDn2AKcAeYBOQbo2UCgFuB1YaYwzwCXCrtf98YIX1eKX1HGv73636SpGREKnNU0rZwJUht0uA9cBwESkSkQUicpOIFAFTgVUi8qFV/WkgHMfoqk3An40x260+iXuBD4HdwHJjzE5rnweBB0QkH0efxUtW+UtAP6v8AeDcMF2lMhIi2FtSrdOJKOVlQR1VMMbMbWfTW23UrcEx7Lat47wHvNdG+UEcfSGty+vaO5ZSIxIjOdvYTMHJWobGh9sdjlI9hl4RrvxSS2e49mso5V2aNJRfGhofTlCA6AgqpbxMk4byS6FBgdZ0ItoZrpQ3adJQfisjMUJHUCnlZZo0lN/KSIjkaOUZTtU22B2KUj2GJg3lt0ZoZ7hSXqdJQ/mtEYk6B5VS3qZJQ/mt+IhQosOCdQSVUl6kSUP5LRFxrK2hZxpKeY0mDeXXMhIi2VdSTVOzTkumlDdo0lB+LSMxgjMNTRSePG13KEr1CJo0lF8bqZ3hSnmVJg3l14bGhxMgsEc7w5XyCk0ayq/1Cg5kcFw4u3Q6EaW8QpOG8nsjEiP1Aj+lvESThvJ7GQkRFFWcoapOpxNRytM0aSi/1zKdyF7tDFfK4zRpKL93bjoR7QxXyuNcWSP8ZREpE5E8p7LZIrJTRJpFJMup/A4RyXW6NYtIprVtoojsEJF8EXlCRMQqjxGR1SKy37qPtsrFqpcvIttFZIL7377qDhIie9G3d7BeGa6UF7hypvEKMLNVWR5wM7DWudAY85oxJtMYkwl8FygwxuRam58FFgLp1q3lmA8Ba4wx6cAa6znAtU51F1r7K/U1IkJGQoTOQaWUF3SYNIwxa4HyVmW7jTF7O9h1LrAEQEQSgUhjzHpjjAEWAzda9bKBRdbjRa3KFxuHDUCUdRylvmZEYiR7S6pp1ulElPIoT/ZpzMFKGkAyUOS0rcgqA+hvjDkGYN3HO+1zpJ19lDrPiMQIauubOFxea3coSnVrHkkaInIJUGuMaekHkTaqdfST0OV9RGShiOSISM7x48c7EanqLjISWqYT0SYqpTzJU2cat/OPswxwnCUMcHo+ACi2Hpe2NDtZ92VO+6S0s895jDHPG2OyjDFZcXFxbghf+Zth/SMIENitV4Yr5VFuTxoiEgDMBpa2lFnNTtUiMsUaNTUPWGFtXgnMtx7Pb1U+zxpFNQU41dKMpVRrvUMCSY3to53hSnmYK0NulwDrgeEiUiQiC0TkJhEpAqYCq0TkQ6ddrgCKjDEHWx3qh8CLQD5wAHjfKn8UmCEi+4EZ1nOA94CDVv0XgB915Q2qnmNEQqTOdquUhwV1VMEYM7edTW+1U/9TYEob5TnA6DbKTwLT2yg3wD0dxadUixGJEazacYyas42Eh3b41VZKdYFeEa66jZbO8L3aGa6Ux2jSUN1GhjUHlXaGK+U5mjRUt5Ec1ZuIXkE67FYpD9KkoboNEWFEQqSeaXjYW1uLePCN7Ti6HVVPo0lDdSsZiRE6nYgH5R6p5F/e2M6ynCPsOHrK7nCU5WxjEwsX57C5sLzjyhdJk4bqVkYkRlJztpGjlWfsDqXbqayt557XthAf0YuQwABW5LZ5ra2yweOr9/PRrlIqaz2/EJkmDdWtZCQ4OsN36UV+btXcbPjn5dsoq67jmTsmcOXwON7ZVkyTntHZ7qtD5fzf2gPMnTyQ6SP6e/z1NGmobmV4QgQisEf7Ndzqhc8PsmZPGb+8fiTjUqLIzkymrPosGw6etDu0Hq26roEHlucyMCaMX14/wiuvqUlDdSthIUGk9uujI6jcaFNBOf/14V6uH5PIvKmDAJg+Ip4+IYGsyD1qc3Q927+/s4viyjP88bZM+njpglZNGqrb0QWZ3OdkzVnufX0LKdG9efSWMVgLbtIrOJBrRifwfl4JdQ1NNkfZM32QV8JfNxdxz1VDmTgo2muvq0lDdTsZCZEUltdy+myj3aH4taZmw0+W5VJR28Azd0wkolfweduzM5Oprmvk0726HIG3lVXX8a9v7WBMcl/un57u1dfWpKG6nRGJERgD+0q1X+NiPP3S9lisAAAV80lEQVRJPp/vP8G/zxrFyKTIr22/bEg/YsNDWLlNm6i8yRjDQ3/bwemzjTw+ZxzBgd79M65JQ3U7IxIdf+D0Ir+uW5d/gsc/3sfN45OZMymlzTpBgQHcMDaJj3eXUV3n+aGeymHJV0f4+54yHr42g6HxEV5/fU0aqttJjupNeKhOJ9JVZVV1/HjpVobGhfO7m0af68doy6zMJOobm/kgr8SLEfZcBSdO89t3d3F5eizzpqbaEoMmDdXtBAQIGQkROuy2CxqbmrlvyVZOn23imTsmEBZy4RE541OiSInpzcpteqGfpzU2NfPT5bkEBwqP3TqOgID2k7knadJQ3VJGYgS7S6p0fqROevzjfWw8VM5/3Dya9P4dN32ICNnjklmXf4Ky6jovRNhzPfvpAbYeruR3N40hoW8v2+LQpKG6pYyESKrrdDqRzvhkbxlPf3KAuZNTuGn8AJf3y85MotnAqu26GrOnbC+q5H/X7GfWuCRmjUuyNRZXlnt9WUTKRCTPqWy2iOwUkWYRyWpVf6yIrLe27xCRXlb5ROt5vog8Ya0VjojEiMhqEdlv3Udb5WLVyxeR7SIywb1vXXVnLZ3h2kTlmuLKM/x0WS4jEiN55NujOrVvev8IRiZG+uVcVF8dKme/j4+yO1PfxE+X5RIbHspvs7+2+KnXuXKm8Qows1VZHnAzsNa5UESCgL8AdxtjRgFXAi3DKp4FFgLp1q3lmA8Ba4wx6cAa6znAtU51F1r7K+WS4QktCzJpZ3hH6hubuef1LTQ2GZ65YwK9ggM7fYzszCRyj1RSePK0ByJ0v6KKWn7wag63/d96vv3UF7y73XcT3h8+2MOB46f579nj6BsW3PEOHtZh0jDGrAXKW5XtNsbsbaP6t4DtxphtVr2TxpgmEUkEIo0x6621vxcDN1r7ZAOLrMeLWpUvNg4bgCjrOEp1KDw0iIExYewp8e1fkb7gvz7Yw9bDlfzhlrGkxfbp0jG+bTWZ+PrZxtnGJp7+JJ9v/vEz1u47wQMzhjE6qS/3vr6Vx1fv87kp9dfuO84rXxZw52WpTEuPtTscwP19GsMAIyIfisgWEfkXqzwZKHKqV2SVAfQ3xhwDsO7jnfY50s4+SnVohNUZfrEqTtfz5pYijpTXuiEq3/JBXgkvfnGIf7o0levHdv03WVJUbyanxfB27lGfHXzw2b7jzPzT5zz24V6uGh7Px//8De6fns5r37+EWycO4H/X7Oe+JVs5U+8b06JU1tbz8ze2MTQ+nAdnZtgdzjnunuEqCJgGTAJqgTUishlo639uR9+stsaTtbmPiCzE0YTFwIEDXQ5WdW8ZCZGs3lXKmfomeod0vsnlaOUZXvz8IEu/OsKZhiZE4Krh8Xx36iC+kR7n9SGPRRW1bC86xYjESFL7hV3w+glXHD5Zy8/f2Ma4AX15+LqL/6OUnZnEL97KY2dxFaOT+1708dzlaOUZfvvOLj7YWUJabB8W3TWZbwyLO7c9NCiQx24dy/D+EfzH+7spLD/NC/OySOzb27aYjTH84u08TtbU89L8SV1qMvQUdyeNIuAzY8wJABF5D5iAo5/DeTjGAKDlPLZURBKNMces5qcyp2OltLPPeYwxzwPPA2RlZfnmzxzldSMSI2i2phMZlxLl8n77Sqt57rMDrLSaWmZlJjEnK4V1+Sd4/asj3PnnTQzqF8b/u2QQs7MGEBUW4qm3QMXpelbtOMaK3KNsKqg4V96vTwgTBkUz0bqNSe7bqT8sdQ1N/Oj1zQjw1HcmEBp08X+UrhudyCMrdrJyW7FPJI2zjU28+Pkhnvz7fgB+fs1wvnd5WpvvVUT4/hWDGRLfh/uX5DLrqXU8/92JjB/ovYkAna3cVsyq7cf4+TXDfeKzdObupPEh8C8iEgbUA98AHrcSQrWITAE2AvOAJ619VgLzgUet+xVO5feKyFLgEuBUSzOWUq7ISLBGUJVUuZQ0cgrKee6zA3y8u4zewYHMm5rKgsvTSI5y/OK8ZHA/7r06nQ92lvDq+gJ+/95u/vujvcwal8S8qamMGeCe/9xn6ptYvbuUFVuP8tm+4zQ2G4bGh/Ozbw1j6pBY9pVWs7mwgs2FFazeVQpAcKAwOrkvEwdGk5UazYRB0cRHtD+W//erdpN3tIoX5mWREhPmlrij+4Rw5fA4VuYW89DMDNsuPgNHX8AjK3dy6MRpZo5K4Jc3jGBAdMfv8+qM/rz5o0tZsGgTc57fwGO3jiU707ut4sWVZ/jl23lMHBTND64Y7NXXdkWHSUNEluAYBRUrIkXAIzg6xp8E4oBVIpJrjLnGGFMhIn8ENuFoSnrPGLPKOtQPcYzE6g28b93AkSyWi8gC4DAw2yp/D7gOyMfR1HXnxb1V1dMMjAkjLCTwgnNQNTcbPtlbxrOfHiCnsILosGB++s1hzJs6iOg+Xz+DCAkKODdWfvexKl7dUMhbW47y181FZKZEMW/qIK4bk9jp5oTGpma+yD/BitxiPtxZQm19EwmRvbhrWhrZmUmMTIw81xw1cVA0cyc7mmFP1px1JJDDFWwprGDxhkJe/OIQACkxvckaFOM4IxkYzfCECAIDhJXbinl1QyE/uGIwM0a6d6W3WZnJfLy7jI2Hypk6pJ9bj+2Ko5Vn+N27u3g/r4TUfmG8cuckrhwe3/GOTob1j2DFPdO4+y+b+fHSXPaX1vDAjGFeSYLNzYaf/XUbTc2GP942jiAvT0boCvHVTquuysrKMjk5OXaHoXzETc+sIyQwgGU/mHpeeUNTM+9sK+a5zw6wr7SG5KjefP/yNG6blNLh1BmtVdU18LfNRby6oZCDx08T0yeE27JSuOOSgRf8FW+MYeuRSlZsPcq7249x8nQ9kb2CuG5MItmZyUxOiyGwk3+ozjY2sbO4ii2FFeQUVJBTWMGJmrOAY0TZ+IFRbCmsYERiJEsWTnH7DKm19Y1k/e5jsjOT+M+bx7r12BfS0hT11N/zMRjuvWoo379i8EU1u9U3NvPrFXks3XSEa0b198pCRy99cYjfvruLR28ew+2Tvds/KyKbjTFZHdbTpKG6s4ff3MF7O46R++sZiAi19Y0s23SEFz8/xNHKMwzvH8HdVw7mhrFJF/0H1BjDuvyTLF5fwMe7SzHA9Ix4vjs1lcuHxp77pZpfVsPK3KOs2FZM4claQoIC+OaIeLIzk7lyeJxb+hecYyqqOENOYbnVpFVJbX0jSxdO8VhH70+WbuWTvcf56hfT3fpe2rN233F+s3InB0+c5ppR/fnVDSNdaopyhTGGP68r4HerdjE8IZIX52eda650t32l1dzw5BdckR7LC/OyLnqgQ2dp0lAKeHV9Ab9asZNV90/jo52lLF5fQEVtA5NTY7j7ysFcNTzeI/85iyvP8PrGwyzddJgTNfWk9gvjmlEJrDtwgryjVQQIXDokluzMJK4ZnUBkL/sv2nKXT/aUcecrm3hhXpbbm7+clZyq49/f3cl7OxxNUY/MGsVVnWyKctWne8u47/WthAYH8H/fzXL7Snn1jc3c+PQ6Sqvq+OAnVxAXEerW47tCk4ZSONa3nv3cegIEmg18c0R/fnjlYCYOivHK659tbOKDvBIWry9kc2EFYwf0JTszmW+PTSQ+0r5J5zypoamZS/5jDZcO6cdT3/HM7D81Zxu54YnPKamq496rhvK9ywd7fFhqflk131uUQ3FlHY/eMoabJ7g+P1d7jDEUn6rjhbUHeeXLAp7/7kS+NSrBDdF2nqtJwzsrkStlk1FJkYxJ7suw/hHc/Y3BLs3c6k6hQYFkZyaTnZlMbX1jp/tL/FFwYADXj0nkr5uPUHO2kXAP9AP8ekUeh8trWbpwKpPTvPMDYGh8BG/fcxk/em0LDyzfxr7SGn5+zXCX+51O1Tawp6SKvaXV7CmpZl9JNXtLq6mucyxLfPukFNsSRmd0/2+w6tHCQoJ4575pdocB0CMSRovszCRe3VDIRztL3PKL3NmK3KO8ueUoP56e7rWE0SIqLIRFd03mNyt38txnB8gvq+ZPt48/LzHWNTSRX1bjSAxWgthbUkVp1dlzdSJ7BZGREMmNmckMT4ggIyGCCTZdE9JZPedbrJTymgkDo0mO6s2K3GK3Jo3DJ2v5xVt5ZA2K5r6rh7rtuJ0RHBjA728aw/CECP7tnV3c8syXXDsmgb3WmUPBidO0TGEVEhRAenw4lw2NZXj/CCtBRNI/MtTrHd3uoklDKeV2AQHCrMwknl97kBM1Z4kNv/iO3YamZn68bCsi8KfbM22/hmHe1FQGx4Zzz+tb+N81+0nt14dh/cO5YWwSGQmOBDEoJsz2ON1Nk4ZSyiOyM5N49tMDvLfjmFvWs/7fj/ez9XAlT31nvNuG1F6saemxbPzX6RhDl+Y380fdKwUqpXxGRkIkGQkRbpkuff2Bkzz9aT63ZQ3ghrH2rlzXWq/gwB6TMECThlLKg2ZlJrG5sOKippWvOF3PT5flktavD7+Z1blVBZX7adJQSnnMt62zgpXbuna2YYzhwb9t5+Tpszwxd3yPGoHmqzRpKKU8JiUmjKxB0by9tWuLM7228TAf7SrlwZkZPjdFeE+lSUMp5VHZmUnst65b6Iz9pdX89t1dXDEsjrsuS/NQdKqzNGkopTzqujGJBAZIpzrE6xqauG/JViJ6BfE/s8fZujaHOp8mDaWUR/ULD+WK9Fje2VZMc7NrTVSPvr+HPSXVPDZ7nC2T96n2adJQSnlcdmYyRyvPkFNY0WHdNbtLeeXLAu66LM1js9aqrtOkoZTyuBkj+9MrOIAVuUcvWK+sqo6fv7GdkYmRPHjtcC9FpzpDk4ZSyuP6hAYxY2QCq3Yco76xuc06zc2GB5Zv40x9E0/MHe+VBZxU53WYNETkZREpE5E8p7LZIrJTRJpFJMupPFVEzohIrnV7zmnbRBHZISL5IvKEWLN1iUiMiKwWkf3WfbRVLla9fBHZLiKemZhfKeUV2eOSqKxt4Iv8421uf/7zg3yRf4JHvj2SofHhXo5OucqVM41XgJmtyvKAm4G1bdQ/YIzJtG53O5U/CywE0q1byzEfAtYYY9KBNdZzgGud6i609ldK+akrhsURFRbc5iiqbUcq+e8P93LdmATmTEqxITrlqg6ThjFmLVDeqmy3MWavqy8iIolApDFmvXFc4bMYuNHanA0ssh4valW+2DhsAKKs4yil/FBIUADXjUnko52l1NY3niuvOdvIj5duJT4ilP+8aazfThneU3iiTyNNRLaKyGcicrlVlgwUOdUpssoA+htjjgFY9/FO+xxpZ5/ziMhCEckRkZzjx9s+9VVK2S97XBJnGppYvav0XNkjK3ZyuLyWP90+nr5h3Wet9O7K3UnjGDDQGDMeeAB4XUQigbZ+OnQ0YNvlfYwxzxtjsowxWXFxcZ0KWCnlPZNSY0js2+tcE9WK3KP8bUsR917t/VX4VNe4NWkYY84aY05ajzcDB4BhOM4SnJfvGgC0NGyWtjQ7WfdlVnkRkNLOPkopPxQQIMwal8TafcfZdqSSX1qr8N1v0yp8qvPcmjREJE5EAq3Hg3F0Yh+0mp2qRWSKNWpqHrDC2m0lMN96PL9V+TxrFNUU4FRLM5ZSyn/NykyisdnwnRc2gI+swqdc1+E8wyKyBLgSiBWRIuARHB3jTwJxwCoRyTXGXANcAfy7iDQCTcDdxpiWTvQf4hiJ1Rt437oBPAosF5EFwGFgtlX+HnAdkA/UAnde1DtVSvmEkYmRpMeHs7+sxqdW4VOuka5MV+zLsrKyTE5Ojt1hKKUu4JO9ZRwoq+F7lw+2OxRlEZHNxpisjurpiiZKKa+7ani8zivlp7QhUSmllMs0aSillHKZJg2llFIu06ShlFLKZZo0lFJKuUyThlJKKZdp0lBKKeUyTRpKKaVc1u2uCBeR40Ch3XF0IBY4YXcQLtA43ctf4gT/iVXjdJ9BxpgOpwnvdknDH4hIjiuX69tN43Qvf4kT/CdWjdP7tHlKKaWUyzRpKKWUcpkmDXs8b3cALtI43ctf4gT/iVXj9DLt01BKKeUyPdNQSinlMk0aHiAiKSLyiYjsFpGdIvLjNupcKSKnRCTXuv3ajlitWApEZIcVx9dWsLKW3H1CRPJFZLuITLAhxuFOn1WuiFSJyE9a1bHtMxWRl0WkTETynMpiRGS1iOy37qPb2Xe+VWe/iMxvq46H43xMRPZY/7ZviUhUO/te8HvihTh/IyJHnf59r2tn35kistf6vj5kQ5zLnGIsEJHcdvb12ufpVsYYvbn5BiQCE6zHEcA+YGSrOlcC79odqxVLARB7ge3X4VieV4ApwEab4w0ESnCMK/eJzxTHUscTgDynsv8CHrIePwT8oY39YoCD1n209Tjay3F+CwiyHv+hrThd+Z54Ic7fAD9z4btxABgMhADbWv/f83Scrbb/D/Bruz9Pd970TMMDjDHHjDFbrMfVwG4g2d6oLko2sNg4bACiRCTRxnimAweMMT5zEacxZi1Q3qo4G1hkPV4E3NjGrtcAq40x5caYCmA1MNObcRpjPjLGNFpPNwADPPX6rmrn83TFZCDfGHPQGFMPLMXx7+ARF4pTRAS4DVjiqde3gyYNDxORVGA8sLGNzVNFZJuIvC8io7wa2PkM8JGIbBaRhW1sTwaOOD0vwt4keDvt/0f0lc8UoL8x5hg4fkgAba1v6muf7V04zirb0tH3xBvutZrRXm6nuc+XPs/LgVJjzP52tvvC59lpmjQ8SETCgb8BPzHGVLXavAVH88o44EngbW/H5+QyY8wE4FrgHhG5otV2aWMfW4bdiUgIMAv4axubfekzdZUvfba/ABqB19qp0tH3xNOeBYYAmcAxHE0/rfnM5wnM5cJnGXZ/nl2iScNDRCQYR8J4zRjzZuvtxpgqY0yN9fg9IFhEYr0cZkssxdZ9GfAWjlN8Z0VAitPzAUCxd6L7mmuBLcaY0tYbfOkztZS2NONZ92Vt1PGJz9bqgL8BuMNYDe6tufA98ShjTKkxpskY0wy80M7r+8rnGQTcDCxrr47dn2dXadLwAKst8yVgtzHmj+3USbDqISKTcfxbnPRelOfi6CMiES2PcXSK5rWqthKYZ42imgKcaml2sUG7v9585TN1shJoGQ01H1jRRp0PgW+JSLTV3PItq8xrRGQm8CAwyxhT204dV74nHtWqH+2mdl5/E5AuImnWWentOP4dvO2bwB5jTFFbG33h8+wyu3viu+MNmIbjlHg7kGvdrgPuBu626twL7MQxumMDcKlNsQ62YthmxfMLq9w5VgGexjEqZQeQZVOsYTiSQF+nMp/4THEksmNAA45fuwuAfsAaYL91H2PVzQJedNr3LiDfut1pQ5z5OPoBWr6rz1l1k4D3LvQ98XKcr1rfv+04EkFi6zit59fhGLF4wI44rfJXWr6XTnVt+zzdedMrwpVSSrlMm6eUUkq5TJOGUkopl2nSUEop5TJNGkoppVymSUMppZTLNGkopZRymSYNpZRSLtOkoZRSymX/HwHBRjXFhsS/AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "nrange = 1+np.asarray(range(N))\n", "plt.plot(nrange,mse)\n", "#plt.fill_between(nrange, mseminus, mseplus, alpha=0.2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Again, we might draw a similar conclusion as for PCA.\n", "Altough the MSE is minimized for 14 components, it is fairly low for other values as well.\n", "\n", "**Task**: Finally, we want to take a look at the declared variance in the response in terms of the number of compontens used in the PLS regression.\n", "You can copy your code from above, and only have to change the `scoring` option.\n", "What do you observe?" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Observation**:" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution**:" ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.3780667415235278\n", "0.3607969231863876\n", "0.3547150074082742\n", "0.3548712309562639\n", "0.31615670731167844\n", "0.31286545476024114\n", "0.314770810251885\n", "0.3292632601937663\n", "0.32483001852412297\n", "0.32200714771909256\n", "0.3219131888346342\n", "0.32325327123174763\n", "0.3230365880832327\n", "0.3343035025073656\n", "0.3304193531453876\n", "0.3279709712678428\n", "0.3308248403083828\n", "0.3304385050644538\n", "0.3245713522246282\n" ] } ], "source": [ "r2 = []\n", "for i in range(N):\n", " pls = PLSRegression(n_components=i+1)\n", " cvs = cross_val_score(pls, X, y, cv=10, scoring='r2')\n", " r2.append(cvs.mean())\n", " print(cvs.mean())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Observation**: Most of the variance in the response variable is declared by taking only one component of partial least squares." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.6.6" } }, "nbformat": 4, "nbformat_minor": 2 }