{"cells": [{"cell_type": "markdown", "metadata": {}, "source": ["# PCA (Principal Component Analysis)\n", "\n", "This notebook shows how to plot a PCA with sciki-learn and statsmodels, with or without normalization."]}, {"cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": ["%matplotlib inline"]}, {"cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": ["import matplotlib.pyplot as plt\n", "plt.style.use('ggplot')"]}, {"cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [{"data": {"text/html": ["<div id=\"my_id_menu_nb\">run previous cell, wait for 2 seconds</div>\n", "<script>\n", "function repeat_indent_string(n){\n", "    var a = \"\" ;\n", "    for ( ; n > 0 ; --n)\n", "        a += \"    \";\n", "    return a;\n", "}\n", "// look up into all sections and builds an automated menu //\n", "var update_menu_string = function(begin, lfirst, llast, sformat, send, keep_item, begin_format, end_format) {\n", "    var anchors = document.getElementsByClassName(\"section\");\n", "    if (anchors.length == 0) {\n", "        anchors = document.getElementsByClassName(\"text_cell_render rendered_html\");\n", "    }\n", "    var i,t;\n", "    var text_menu = begin;\n", "    var text_memo = \"<pre>\\nlength:\" + anchors.length + \"\\n\";\n", "    var ind = \"\";\n", "    var memo_level = 1;\n", "    var href;\n", "    var tags = [];\n", "    var main_item = 0;\n", "    var format_open = 0;\n", "    for (i = 0; i <= llast; i++)\n", "        tags.push(\"h\" + i);\n", "\n", "    for (i = 0; i < anchors.length; i++) {\n", "        text_memo += \"**\" + anchors[i].id + \"--\\n\";\n", "\n", "        var child = null;\n", "        for(t = 0; t < tags.length; t++) {\n", "            var r = anchors[i].getElementsByTagName(tags[t]);\n", "            if (r.length > 0) {\n", "child = r[0];\n", "break;\n", "            }\n", "        }\n", "        if (child == null) {\n", "            text_memo += \"null\\n\";\n", "            continue;\n", "        }\n", "        if (anchors[i].hasAttribute(\"id\")) {\n", "            // when converted in RST\n", "            href = anchors[i].id;\n", "            text_memo += \"#1-\" + href;\n", "            // passer \u00e0 child suivant (le chercher)\n", "        }\n", "        else if (child.hasAttribute(\"id\")) {\n", "            // in a notebook\n", "            href = child.id;\n", "            text_memo += \"#2-\" + href;\n", "        }\n", "        else {\n", "            text_memo += \"#3-\" + \"*\" + \"\\n\";\n", "            continue;\n", "        }\n", "        var title = child.textContent;\n", "        var level = parseInt(child.tagName.substring(1,2));\n", "\n", "        text_memo += \"--\" + level + \"?\" + lfirst + \"--\" + title + \"\\n\";\n", "\n", "        if ((level < lfirst) || (level > llast)) {\n", "            continue ;\n", "        }\n", "        if (title.endsWith('\u00b6')) {\n", "            title = title.substring(0,title.length-1).replace(\"<\", \"&lt;\")\n", "         .replace(\">\", \"&gt;\").replace(\"&\", \"&amp;\");\n", "        }\n", "        if (title.length == 0) {\n", "            continue;\n", "        }\n", "\n", "        while (level < memo_level) {\n", "            text_menu += end_format + \"</ul>\\n\";\n", "            format_open -= 1;\n", "            memo_level -= 1;\n", "        }\n", "        if (level == lfirst) {\n", "            main_item += 1;\n", "        }\n", "        if (keep_item != -1 && main_item != keep_item + 1) {\n", "            // alert(main_item + \" - \" + level + \" - \" + keep_item);\n", "            continue;\n", "        }\n", "        while (level > memo_level) {\n", "            text_menu += \"<ul>\\n\";\n", "            memo_level += 1;\n", "        }\n", "        text_menu += repeat_indent_string(level-2);\n", "        text_menu += begin_format + sformat.replace(\"__HREF__\", href).replace(\"__TITLE__\", title);\n", "        format_open += 1;\n", "    }\n", "    while (1 < memo_level) {\n", "        text_menu += end_format + \"</ul>\\n\";\n", "        memo_level -= 1;\n", "        format_open -= 1;\n", "    }\n", "    text_menu += send;\n", "    //text_menu += \"\\n\" + text_memo;\n", "\n", "    while (format_open > 0) {\n", "        text_menu += end_format;\n", "        format_open -= 1;\n", "    }\n", "    return text_menu;\n", "};\n", "var update_menu = function() {\n", "    var sbegin = \"\";\n", "    var sformat = '<a href=\"#__HREF__\">__TITLE__</a>';\n", "    var send = \"\";\n", "    var begin_format = '<li>';\n", "    var end_format = '</li>';\n", "    var keep_item = -1;\n", "    var text_menu = update_menu_string(sbegin, 2, 4, sformat, send, keep_item,\n", "       begin_format, end_format);\n", "    var menu = document.getElementById(\"my_id_menu_nb\");\n", "    menu.innerHTML=text_menu;\n", "};\n", "window.setTimeout(update_menu,2000);\n", "            </script>"], "text/plain": ["<IPython.core.display.HTML object>"]}, "execution_count": 4, "metadata": {}, "output_type": "execute_result"}], "source": ["from jyquickhelper import add_notebook_menu\n", "add_notebook_menu()"]}, {"cell_type": "markdown", "metadata": {}, "source": ["More about [PCA](https://en.wikipedia.org/wiki/Principal_component_analysis): [Implementing a Principal Component Analysis (PCA) in Python step by step](http://sebastianraschka.com/Articles/2014_pca_step_by_step.html)."]}, {"cell_type": "markdown", "metadata": {}, "source": ["## Download data"]}, {"cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [{"data": {"text/plain": ["'auto-mpg.data'"]}, "execution_count": 5, "metadata": {}, "output_type": "execute_result"}], "source": ["import pyensae.datasource\n", "pyensae.datasource.download_data(\"auto-mpg.data\", url=\"https://archive.ics.uci.edu/ml/machine-learning-databases/auto-mpg/\")"]}, {"cell_type": "code", "execution_count": 5, "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>mpg</th>\n", "      <th>cylinders</th>\n", "      <th>displacement</th>\n", "      <th>horsepower</th>\n", "      <th>weight</th>\n", "      <th>acceleration</th>\n", "      <th>year</th>\n", "      <th>origin</th>\n", "      <th>name</th>\n", "    </tr>\n", "  </thead>\n", "  <tbody>\n", "    <tr>\n", "      <th>0</th>\n", "      <td>18.0</td>\n", "      <td>8</td>\n", "      <td>307.0</td>\n", "      <td>130.0</td>\n", "      <td>3504.0</td>\n", "      <td>12.0</td>\n", "      <td>70</td>\n", "      <td>1</td>\n", "      <td>chevrolet chevelle malibu</td>\n", "    </tr>\n", "    <tr>\n", "      <th>1</th>\n", "      <td>15.0</td>\n", "      <td>8</td>\n", "      <td>350.0</td>\n", "      <td>165.0</td>\n", "      <td>3693.0</td>\n", "      <td>11.5</td>\n", "      <td>70</td>\n", "      <td>1</td>\n", "      <td>buick skylark 320</td>\n", "    </tr>\n", "    <tr>\n", "      <th>2</th>\n", "      <td>18.0</td>\n", "      <td>8</td>\n", "      <td>318.0</td>\n", "      <td>150.0</td>\n", "      <td>3436.0</td>\n", "      <td>11.0</td>\n", "      <td>70</td>\n", "      <td>1</td>\n", "      <td>plymouth satellite</td>\n", "    </tr>\n", "    <tr>\n", "      <th>3</th>\n", "      <td>16.0</td>\n", "      <td>8</td>\n", "      <td>304.0</td>\n", "      <td>150.0</td>\n", "      <td>3433.0</td>\n", "      <td>12.0</td>\n", "      <td>70</td>\n", "      <td>1</td>\n", "      <td>amc rebel sst</td>\n", "    </tr>\n", "    <tr>\n", "      <th>4</th>\n", "      <td>17.0</td>\n", "      <td>8</td>\n", "      <td>302.0</td>\n", "      <td>140.0</td>\n", "      <td>3449.0</td>\n", "      <td>10.5</td>\n", "      <td>70</td>\n", "      <td>1</td>\n", "      <td>ford torino</td>\n", "    </tr>\n", "  </tbody>\n", "</table>\n", "</div>"], "text/plain": ["    mpg  cylinders  displacement horsepower  weight  acceleration  year  \\\n", "0  18.0          8         307.0      130.0  3504.0          12.0    70   \n", "1  15.0          8         350.0      165.0  3693.0          11.5    70   \n", "2  18.0          8         318.0      150.0  3436.0          11.0    70   \n", "3  16.0          8         304.0      150.0  3433.0          12.0    70   \n", "4  17.0          8         302.0      140.0  3449.0          10.5    70   \n", "\n", "   origin                       name  \n", "0       1  chevrolet chevelle malibu  \n", "1       1          buick skylark 320  \n", "2       1         plymouth satellite  \n", "3       1              amc rebel sst  \n", "4       1                ford torino  "]}, "execution_count": 6, "metadata": {}, "output_type": "execute_result"}], "source": ["import pandas\n", "df = pandas.read_fwf(\"auto-mpg.data\", encoding=\"utf-8\",\n", "                    names=\"mpg cylinders displacement horsepower weight acceleration year origin name\".split())\n", "df[\"name\"] = df[\"name\"].apply(lambda s : s.strip(' \"'))\n", "df.head()"]}, {"cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [{"data": {"text/plain": ["mpg             float64\n", "cylinders         int64\n", "displacement    float64\n", "horsepower       object\n", "weight          float64\n", "acceleration    float64\n", "year              int64\n", "origin            int64\n", "name             object\n", "dtype: object"]}, "execution_count": 7, "metadata": {}, "output_type": "execute_result"}], "source": ["df.dtypes"]}, {"cell_type": "markdown", "metadata": {}, "source": ["We remove missing values:"]}, {"cell_type": "code", "execution_count": 7, "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>mpg</th>\n", "      <th>cylinders</th>\n", "      <th>displacement</th>\n", "      <th>horsepower</th>\n", "      <th>weight</th>\n", "      <th>acceleration</th>\n", "      <th>year</th>\n", "      <th>origin</th>\n", "      <th>name</th>\n", "    </tr>\n", "  </thead>\n", "  <tbody>\n", "    <tr>\n", "      <th>32</th>\n", "      <td>25.0</td>\n", "      <td>4</td>\n", "      <td>98.0</td>\n", "      <td>?</td>\n", "      <td>2046.0</td>\n", "      <td>19.0</td>\n", "      <td>71</td>\n", "      <td>1</td>\n", "      <td>ford pinto</td>\n", "    </tr>\n", "    <tr>\n", "      <th>126</th>\n", "      <td>21.0</td>\n", "      <td>6</td>\n", "      <td>200.0</td>\n", "      <td>?</td>\n", "      <td>2875.0</td>\n", "      <td>17.0</td>\n", "      <td>74</td>\n", "      <td>1</td>\n", "      <td>ford maverick</td>\n", "    </tr>\n", "    <tr>\n", "      <th>330</th>\n", "      <td>40.9</td>\n", "      <td>4</td>\n", "      <td>85.0</td>\n", "      <td>?</td>\n", "      <td>1835.0</td>\n", "      <td>17.3</td>\n", "      <td>80</td>\n", "      <td>2</td>\n", "      <td>renault lecar deluxe</td>\n", "    </tr>\n", "    <tr>\n", "      <th>336</th>\n", "      <td>23.6</td>\n", "      <td>4</td>\n", "      <td>140.0</td>\n", "      <td>?</td>\n", "      <td>2905.0</td>\n", "      <td>14.3</td>\n", "      <td>80</td>\n", "      <td>1</td>\n", "      <td>ford mustang cobra</td>\n", "    </tr>\n", "    <tr>\n", "      <th>354</th>\n", "      <td>34.5</td>\n", "      <td>4</td>\n", "      <td>100.0</td>\n", "      <td>?</td>\n", "      <td>2320.0</td>\n", "      <td>15.8</td>\n", "      <td>81</td>\n", "      <td>2</td>\n", "      <td>renault 18i</td>\n", "    </tr>\n", "    <tr>\n", "      <th>374</th>\n", "      <td>23.0</td>\n", "      <td>4</td>\n", "      <td>151.0</td>\n", "      <td>?</td>\n", "      <td>3035.0</td>\n", "      <td>20.5</td>\n", "      <td>82</td>\n", "      <td>1</td>\n", "      <td>amc concord dl</td>\n", "    </tr>\n", "  </tbody>\n", "</table>\n", "</div>"], "text/plain": ["      mpg  cylinders  displacement horsepower  weight  acceleration  year  \\\n", "32   25.0          4          98.0          ?  2046.0          19.0    71   \n", "126  21.0          6         200.0          ?  2875.0          17.0    74   \n", "330  40.9          4          85.0          ?  1835.0          17.3    80   \n", "336  23.6          4         140.0          ?  2905.0          14.3    80   \n", "354  34.5          4         100.0          ?  2320.0          15.8    81   \n", "374  23.0          4         151.0          ?  3035.0          20.5    82   \n", "\n", "     origin                  name  \n", "32        1            ford pinto  \n", "126       1         ford maverick  \n", "330       2  renault lecar deluxe  \n", "336       1    ford mustang cobra  \n", "354       2           renault 18i  \n", "374       1        amc concord dl  "]}, "execution_count": 8, "metadata": {}, "output_type": "execute_result"}], "source": ["df[df.horsepower == \"?\"]"]}, {"cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": ["final = df[df.horsepower != '?'].copy()\n", "final[\"horsepower\"] = final[\"horsepower\"].astype(float)"]}, {"cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": ["final.to_csv(\"auto-mpg.data.csv\", sep=\"\\t\", index=False, encoding=\"utf-8\")"]}, {"cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [{"data": {"text/plain": ["(392, 9)"]}, "execution_count": 11, "metadata": {}, "output_type": "execute_result"}], "source": ["final.shape"]}, {"cell_type": "markdown", "metadata": {}, "source": ["## PCA with scikit-learn"]}, {"cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [{"data": {"text/plain": ["PCA(copy=True, iterated_power='auto', n_components=2, random_state=None,\n", "  svd_solver='auto', tol=0.0, whiten=False)"]}, "execution_count": 12, "metadata": {}, "output_type": "execute_result"}], "source": ["from sklearn.decomposition import PCA\n", "X = final[df.columns[1:-1]]\n", "Y = final[\"mpg\"]\n", "pca = PCA(n_components=2)\n", "pca.fit(X)"]}, {"cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[536.44492922,  50.83312832],\n", "       [730.34140206,  79.13543921],\n", "       [470.9815846 ,  75.4476722 ],\n", "       [466.40143367,  62.53420646],\n", "       [481.66788465,  55.78036021]])"]}, "execution_count": 13, "metadata": {}, "output_type": "execute_result"}], "source": ["out = pca.transform(X)\n", "out[:5]"]}, {"cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [{"data": {"text/plain": ["(array([0.99756151, 0.0020628 ]), 55.14787750463889)"]}, "execution_count": 14, "metadata": {}, "output_type": "execute_result"}], "source": ["pca.explained_variance_ratio_, pca.noise_variance_"]}, {"cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [{"data": {"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAD8CAYAAAB6paOMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztnX14FOW5/78zuwGJgbDJBsJLoCaAEU8oYACJSNSsaNF6aKoWtPoTX/j5i61XjXourS/1FPGK1Yjigdr6Qo9WT4GrQnt+9Sc2RJNiQKIQ5JAGTRALGgjJLnnhpXmZ5/fH7M7OzM4ms7szu7OZ+3NdXrKTmdl7nn3mvp/nfu77fjjGGANBEARhW/hEC0AQBEEkFjIEBEEQNocMAUEQhM0hQ0AQBGFzyBAQBEHYHDIEBEEQNocMAUEQhM0hQ0AQBGFzyBAQBEHYHDIEBEEQNseZaAH08u2330r/drvdaG9vT6A0+iFZzYFkNQeS1RwSJevEiRN1nUczAoIgCJtDhoAgCMLmkCEgCIKwOWQICIIgbA4ZAoIgCJtDhoAgCMLmkCEgbAFraYLw3hawlqZEi0IQliNp8ggIIlpYSxOEyseB/n4wpxP8g0+Dy8tPtFgEYRloRkAMe9ihA0B/P8AEYKBf/EwQhAQZAmLYw11YADidAM8DDqf4mSAICXINEcMeLi8f/INPgx06AO7CAnILEYQKMgSELeDy8skAEEQYyDVEEARhc8gQEARB2BwyBARBEDYn5jWC9vZ2rF+/HqdOnQLHcfB4PFi6dCk2b96MHTt2YMyYMQCAFStWYO7cuQCArVu3orq6GjzPY+XKlZg9e3asYhAEQRBRErMhcDgcuO2225Cbm4uzZ8/ikUcewaxZswAA1113HW644QbF+ceOHUNdXR1eeOEF+Hw+rF69Gi+99BJ4niYnBEEQiSBm7etyuZCbmwsAGDVqFCZNmgSv1xv2/Pr6ehQVFSElJQXjxo1DdnY2mpubYxWDIAiCiBJDw0fb2trw1VdfYdq0aWhqasL27dtRW1uL3Nxc3H777UhLS4PX68X06dOlazIyMjQNR1VVFaqqqgAAFRUVcLvdQaGdTsVnK0OymgPJag4kqzlYXVbDDMG5c+dQWVmJO+64A6mpqViyZAluvPFGAMCmTZvw5ptvoqysDIwxXffzeDzweDzSZ/l+n7RXqTmQrOZAspoDyTo0cd2zuL+/H5WVlbj88suxYMECAMDYsWPB8zx4nkdJSQlaWloAAJmZmejo6JCu9Xq9yMjIMEIMgiAIIgpiNgSMMbzyyiuYNGkSrr/+eum4z+eT/r1nzx7k5OQAAAoLC1FXV4e+vj60tbWhtbUV06ZNi1UMgiAIIkpidg0dOnQItbW1mDJlCh5++GEAYqjoxx9/jCNHjoDjOGRlZWHVqlUAgJycHCxcuBDl5eXgeR533XUXRQwRBEEkkJgNQX5+PjZv3hxyPJAzoEVpaSlKS0tj/WqCIAjCAGgoThAEYXPIEBAEQdgcMgQEQRA2hwwBQRCEzSFDQBAEYXPIEBAEQdgcMgQEQRA2hwwBQRCEzSFDQBAEYXPIEBAEQdgcMgQEQRA2hwwBQRCEzSFDQBAEYXPIEBAEQdgcMgQEQRA2hwwBQRCEzYl5Y5r29nasX78ep06dAsdx8Hg8WLp0KXp6erB27VqcPHkSWVlZeOCBB5CWlgbGGDZu3Ih9+/Zh5MiRKCsrQ25urhHPQhAEQURBzDMCh8OB2267DWvXrsWaNWuwfft2HDt2DNu2bUNBQQHWrVuHgoICbNu2DQCwb98+HD9+HOvWrcOqVavw2muvxfwQBEEQRPTEbAhcLpc0oh81ahQmTZoEr9eL+vp6FBcXAwCKi4tRX18PAPj000+xePFicByHGTNm4PTp04qN7gmCIIj4ErNrSE5bWxu++uorTJs2DZ2dnXC5XABEY9HV1QUA8Hq9cLvd0jWZmZnwer3SuQGqqqpQVVUFAKioqFBc43Q6FZ+tDMlqDiSrOZCs5mB1WQ0zBOfOnUNlZSXuuOMOpKamhj2PMRZyjOO4kGMejwcej0f63N7eLv3b7XYrPlsZktUcSFZzIFnNIVGyTpw4Udd5hkQN9ff3o7KyEpdffjkWLFgAAEhPT5dcPj6fD2PGjAEgzgDkDdLR0REyGyAIgiDiR8yGgDGGV155BZMmTcL1118vHS8sLERNTQ0AoKamBvPmzZOO19bWgjGGL774AqmpqWQICIIgEkjMrqFDhw6htrYWU6ZMwcMPPwwAWLFiBZYtW4a1a9eiuroabrcb5eXlAIA5c+Zg7969uP/++zFixAiUlZXFKgJBEAQRAzEbgvz8fGzevFnzb08++WTIMY7jcPfdd8f6tQRBEIRBUGYxQRCEzSFDQBAEYXPIEBAEQdgcMgQEQRA2hwwBQRCEzSFDQBAEYXPIEBAEQdgcMgQEQRA2hwwBQRCEzSFDQBAEYXPIEBAEQdgcMgQEQRA2hwwBQRCEzSFDQBAEYXPIEBAEQdgcMgQEQRA2x5DN6zds2IC9e/ciPT0dlZWVAIDNmzdjx44d0l7FK1aswNy5cwEAW7duRXV1NXiex8qVKzF79mwjxCAIgiCiwBBDcMUVV+Daa6/F+vXrFcevu+463HDDDYpjx44dQ11dHV544QX4fD6sXr0aL730EnieJicEQRCJwBDtO3PmTKSlpek6t76+HkVFRUhJScG4ceOQnZ2N5uZmI8QgCIIgosCQGUE4tm/fjtraWuTm5uL2229HWloavF4vpk+fLp2TkZEBr9drphgEQRDEIJhmCJYsWYIbb7wRALBp0ya8+eabKCsrA2NM1/VVVVWoqqoCAFRUVMDtdkt/czqdis9WhmQ1B5LVHEhWc7C6rKYZgrFjx0r/LikpwbPPPgsAyMzMREdHh/Q3r9eLjIyMkOs9Hg88Ho/0ub29Xfq32+1WfLYyJKs5kKzmQLKaQ6JknThxoq7zTFuh9fl80r/37NmDnJwcAEBhYSHq6urQ19eHtrY2tLa2Ytq0aWaJQRAEQQyBITOCF198EY2Njeju7sa9996Lm2++GQcPHsSRI0fAcRyysrKwatUqAEBOTg4WLlyI8vJy8DyPu+66iyKGCIIgEoghhuBnP/tZyLGrrroq7PmlpaUoLS014qsJgiCIGKGhOEEQhM0hQ0AQBGFzyBAQBEHYHDIEBEEQNocMAUEkGNbSBOG9LWAtTYkWhbApppaYIAhicFhLE4TKx4H+fjCnE/yDT4PLy0+0WITNoBkBQSQQdugA0N8PMAEY6Bc/E0ScIUNAEAmEu7AAcDoBngccTvEzQcQZcg0RRALh8vLBP/g02KED4C4sILcQkRDIEBBEguHy8skAEAmFXEMEQRA2hwwBQQwDhlsI6nB7HqtDriGCSHKGWwjqcHueZIBmBASR5Ay3ENTh9jzJABkCgkhyhlsI6mDPQy4jcyDXEEEkOXpCUFlLU9KEqIZ7HrXLqPeXLwPuCQmWdnhgiCHYsGED9u7di/T0dFRWVgIAenp6sHbtWpw8eRJZWVl44IEHkJaWBsYYNm7ciH379mHkyJEoKytDbm6uEWIQhG0ZLAQ1GXzuakOl9Txql1HfwX1AMRkCIzDENXTFFVfg5z//ueLYtm3bUFBQgHXr1qGgoADbtm0DAOzbtw/Hjx/HunXrsGrVKrz22mtGiEAQRBis7nMPGCq27W3x/2HcPmqXUcrFc+Is6fDFEEMwc+ZMpKWlKY7V19ejuLgYAFBcXIz6+noAwKefforFixeD4zjMmDEDp0+fVmx0TxCEsVh9DUGvoQq4jLh/vRX8g09jRL61niOZMW2NoLOzEy6XCwDgcrnQ1dUFAPB6vXC73dJ5mZmZ8Hq90rkEQRiL1ctYcBcWgDmdwED/kIaKsrDNIe6LxYyxkGMcx4Ucq6qqQlVVFQCgoqJCYTycTqfis5UhWc2BZI0Q9yJgwaIhT0uIrO5F6P3ly+g7uA8pF8/RPdK3RLvqxOqymmYI0tPT4fP54HK54PP5MGbMGADiDKC9vV06r6OjQ3M24PF44PF4pM/ya9xut+KzlSFZzYFkNYeEyeqeABRPwFkA0Pn91K5DM3HiRF3nmZZHUFhYiJqaGgBATU0N5s2bJx2vra0FYwxffPEFUlNTyS1EEASRQAyZEbz44otobGxEd3c37r33Xtx8881YtmwZ1q5di+rqarjdbpSXlwMA5syZg7179+L+++/HiBEjUFZWZoQIBEEQRJQYYgh+9rOfaR5/8sknQ45xHIe7777biK8lCNuQTAlhRPJBmcWEbUkW5ZoMCWFGkSy/yXCDDAFhS9TKlVt+D9DTZUkFpBVnbzUZjcBOBs9qkCEgbIlCufb3gb3zCiAwSyqgSOLskxm7GDwrQoaAsCUK5QoOEASAsZgUkFluDXlCGNLGSJm3w01J2sXgWREyBIQtCVGuf3g1JgVkpFtDy6AE/j+cXSdWz4AezpAhIGIimRf35OUK2KSpMT1HiFtjV7V0P7iHzuiV7jOIQYnUdWL13yacwbOirMMdMgRE1Aynxb1YFZDCrcE7wHZWAYIQcd38wZR9JK4Tq/82VpfPbtAOZUTUWL28cTyRV8ZEUYl/zUFWN1/vfQapFKquvjnobMDiv43V5bMbNCMgooYW95QEZhWspQnCrmqpXVIuniPW0NF5j8H85HpnLlb/bawun93gmFY5UAvy7bffSv+mYlPmEI2sifJDW71d5e2StWBRQmSN5rcxq121ZIm171i9D8ixetE5mhEQMUGLe9pYoV2sIAMQfj3AKvIRtEZAJCGspQmn//hm2C0NCWtB6wHWh2YERFIRGF32+H3LVow2EWq3g332MbhLLgO/+JpEi5NwaD3A+pAhIJIKxegS1itDINRuB3trPQCANTZAAIDSWxMqU6KhRDHrQ64hIqmw/Ebsn3086Ge7wuXlg19605AJcMJ7W8jllwBoRkAkhGgjRgKjy9Rjh3Fmcq7lRpfcJZeBNTYED4wYid6mA7oTyiLB6pnDkUAJZonFdENw33334bzzzgPP83A4HKioqEBPTw/Wrl2LkydPIisrCw888ADS0tLMFoWwCFovPQDdSo3Ly8f5CxbhrD8cz0oKkV98DQQA7G8fAEcPA/vr4WtsAF++2lDZzKptFEk5DCOJZ+VRPf3FSn0qHsRlRvCLX/xC2rweALZt24aCggIsW7YM27Ztw7Zt2/DjH/84HqIQFkD90gt11cCu6qiUmhVHkvziayD0dIF93RIsc22wYjNKcarbL5JyGEYSrwVlPf3Fin3KbBKyRlBfX4/i4mIAQHFxMerr6xMhBpEgQvz8HKIOL7RqaKLiGZ0poqIz0Adu1FqJuv0iKYdhJJGUz4gFPf3Fqn3KTOIyI1izZg0A4Oqrr4bH40FnZydcLhcAwOVyoaurKx5iEBZBHUUCAKyuOqrRoFVDE+XPOHb+InR2dho6yjQqEkfdfpGUwzCaeCSY6ekvVu1TZmJ6iQmv14uMjAx0dnbi6aefxsqVK/GrX/0Kv/vd76RzVq5ciY0bNyquq6qqQlVVFQCgoqICvb290t+cTif6+/vNFNswhpusvU0H0HdwH1IunoMR+ca9IJHeVy6rWTIZhdPpROemN9Dzzm/FYnQ8j7RbVuH8H96eaNEAKNsv9V/mWKq/ymUDoPido3239PQXo/tUovTAiBEjdJ1n+owgIyMDAJCeno558+ahubkZ6enp8Pl8cLlc8Pl8ivWDAB6PBx6PR/osr9NBNUbMYShZ5b5TGO07dU8AiieIo1Edi8AKWTWutRJutxtnJucCDicAcZR5ZnKutNidcGTtN6K/3zL9VdHfeB7gOGBgQOp7Uddw0tNfDOxTrKUpYVFulqg1dO7cOTDGMGrUKJw7dw6ff/45brzxRhQWFqKmpgbLli1DTU0N5s2bZ6YYhEHEO7JDqHwc6OsD43lwt947aJau1aM8rJZUZfX2AlT9TWDiVqJA0G+/IDERTpGQDJnwgMmGoLOzE88//zwAYGBgAIsWLcLs2bORl5eHtWvXorq6Gm63G+Xl5WaKERPJ8MLEi3C+UzPaiB06APT1AWCAMAD2zitgk6Zq3j9ZojzC+cDj3ccGC9/tnb8oIVFDWij6G+efEQgDSeW3t3omfABTDcH48ePx3HPPhRwfPXo0nnzySTO/WiKWlyxZFEy80BrVmtVG3IUFYDwvvviAuNtXmO0f4zlTMZpE9DHNbTXrxPBd3182g/vR3UBPV8IHP5pBBUk2KEuWhedhnVkc60umeGH6+yD8+R3wN9ySNJ3QDNSjWrOUMJeXD+7We8HeeUVcYHU4w27/aNTLlojZXyztF5AXaWMiUtzq9mIMQRn6+sDefgVgzBKDH3V/k+9lcLrmMJgFs8vlWD0TPsDwNgQxKinphenvE/2TjfshfNmY8JfDSpg54uEXXyNtKs86TgJ/+0AZ717sNwQG+N8TNfvT237qiqZBef19E5x4n6IS8EVXDSq71khb2lGN44ABAQCz7OwqWfzuAdSZ8FZkWBsC9UuGtDEQ3tuif+Tkf2GEP78DNO4HwMQR065qS3e8eGL2Iqje7R9jjUFPlHtJT/tpVjTt6fLLG4j+ZqJRqN0OYVf1kMpR3V4BGdLGT0D36y/qMuyJWj9LFr97MjG8DYHsJUPaGLA/vBrRiE+o3S7WjEkZATh4MXQNDGxnFdhC5agrYS+FBRaz45IIpFKYI/ILDAkVVbhXEuTLHar9tCqa8jfcopytBv8a3ezXL0Oq243TY926avEkav0sWfzuycSwNgRAsIML722JaMQnH4WJN+JkfxxQXJ+olyKR35sI4zOkwoywmBgARV4Et/weSyySqlFXNOUuuSxkkIN/tIhrKEyIWTmGa2e5eyo4I4n/Aj2Xlw9u+T1IOVCPvoJ5lvqtkpVhbwgC6B1FBBQF27db9Qem+Dc7eyb4MYxbwewFrUS4M2JegDfJiERTTIwrukrRfujpAr/0JsNkMgqpoqlq17OQhfuFV5lmoNXuKUyfKSZ5McR9VM5amsD+8Cp6B/qBgw1hw4oJ/djHEOjwxbKWJgjPPyZOt3mH6ga83yfp54OtYLMXiC+j3MjwDrCOk+KL84dXI17QikRRGjlF1mu0Yo1yMWsGo0eukHMYEuYOihR+8TXAENtemumiC9lg58tGwJkCXLYE3NQ8RWE2s2eLtEZgPLYxBIAO18KuatEIAGL8+vSLgb5eYGwmuIJLwH7/66AxEAQIddVw+O/JP/i0GI+9swr42wdgOznxXDE2T1dnDRmxhnFVyI2FEQu1kURhxGJ8zJzBRFNMjC+6CigybxSdaIycfYVsuAMAwgA4DsG1N1kZCDNdlbRGYDy2MgRDEVJ+b0IOHLeVSR8HDnwGNARdRvJlAy4vX1R0guA3Fpx/6sx0d1Z13oIYQ6+M59YaVfNLb5JKHEeVOBfBCCuWKCEzX2A9cklRYHXV0m8Xj4XuRGD07EtyT+38APj6MACxXytyEAZkM2YTXZXJEpufTJAhkMEXXQXh4ypxNsA7xBGj/O/XlkL4n8/Eju8Qt3JgLU1SRwxRdMvvwflsQHdnVVwPzm9UlFEgQl11sPSCrFZ6LC99pAo6WuUZTlkbNXIN8ZmHu++uarD+frC6ocMsrcpQbWbG7CvgngpZcA/MpOUjKd5h6kg9GWLzkwkyBDK4vHzwDz8junhUswOp869YJe48VbcDrPYDSZkAgPD+u0DmOCB7MvhrS8Hl5SOlvRVsz07p/kN+vzzc9b9+66+86JA2NsHHf4W4QgeA48XjsSbOmTDCCqeotJS1GesG4e5rpIJMZMjwUG1m9uxLnYOgyLXhOKCoRJrBDlfX23CCDIEGgborgcQcQBVmWHQV2MCAslbLzir/SB7A8WMQjh8DZi+Ar/r/ikloOpVc4CVjLU2Qtorw/1+oq/bnMvgpuCQ4qo7xpTdyhBWJcjdr3SDcfSNVkOqM3mie0Wj0tFkkLrxY9yzm8vLB33ALhC8bFesvVKsreSBDoCLsNnWDRJswhqARCHD8GPD+seDnCJVcUOkzgAnitZw0FwAAcOniLm/RvvRxieoYokaTWSPXwbLK9baVVkZvwBiE6yfxGP3qbTM9Ljyj9izW6oOK3J3+Pim4wk4E3jcrVXXVggyBirCllgeJNgEgri2ojYHixnxI2eZwxcJYSxNQtwNqFxAAsJ1/FQ0ExwFT8oK3D8wM/ApJy1ceTXZ1NITUaPp7+BpNsSw+DyqD2s2mem49+QJaGb2BEE5NQxOn0a+RbaY2aPIaTtHIpXADBirIDvjXuup2gA1RBynemDkwkhtZ3182gy9fbalnl+N46qmnnkq0EHro7u6W/p2amoozZ84Mcnb0cBlucPmzAPd48N9fIXbuwDGeByZfAC7nAvH49IvFv2W4wV30XbAT3wLek2FuzIG7zAP42iH86lHg7w3A5/XAoQNgn9QAo9PBDu4T8xAOHQD+7ve3ggO3eAn4K74H+NpFt1UghPXgPnAXfRdchlvsdM/9HGjcD7brQ+XxyseBRv/3BWYZwgCQlQ1u+sWGt2ugvZj3JNB+wu/aYoB7vOL7FOf721IOa2kC2/2huEai+pseWQP3ZQf3if5rJgwqRwh9vWKbBZg0FezoV2Db3wXSxoC/plTqJ2g9GvY7zOiv4dpMjboNQ9qUd4B98hECUUDn/+hO/DN1tGEy4lQH8HVz8KC/XQb7bfUSa7sG3439YJ98BC5/VtSyaN5/94fBPsEi6HcGMnq0vt+SZgQahJtSq9cOQhZB/2UuWPPflYlnAYQBCH/8nZiaL8j8/EwsFiaV/uV5oKBQmbW5UIxeYocOKGcdslIXQmDzdwAY6Ifw/h/huO8x/4hPXY/G/9VnI3uJIkp20/AbR5RzYKB/OVr3U0jIZMMnAD4R5WtsAJt9qRQUAMS+TmM0Wnkp4szIv+vbLeKub2bUcArALbxKfG9k7WKVtQOzM/MV/c6ZYok+EQ4yBDrRtUB3YQGYwwH0axgCAGj+u6ZCFndeCpT+FcRcBY4Dpk4Dd/kSVXhqSjDpDZzo9mhpAv7Rorzn/noxtFWang8gBPU14Z69pUk0NHU7IkoWisWFYeRLGqkcioS9xddA6OkCO6LRVg27IRzcK7WFOkch0YS04WcfBwcFAwNgbwd3fYunGyvSul9mYXZimvzZx85fhC5aIwiloaEBGzduhCAIKCkpwbJlyxIlii50dxotRT/U32bNBw7uFV0R8nOPfAn2zdeKWircZSVg3x4FWkSjwv7rt2J0UYiiZ2CHDoBfehO4W+4VZxyC8hzuksuGeGrZqDKQuwAEI6V0KNaocw4MfkmHkiMQHYQpecCO/1aOVgOyyH+fAGpFZqEchZB1jCl5fpejn0AQgskyaq4dWGD2ZNb6lPo7uLx8jHC7DZ1pGU1CDIEgCHj99dfx+OOPIzMzE48++igKCwsxefLkRIijCz2dRsosHvRGnFi3KKCUnSngry0FK7hEWe00QF8fhPf/CIxxAYEFaZ4XN/MGE0dW4rcH7++/b+AFU2zwcvYM8I+WkHDIcEijSileiRP9yrLdwsxQeEMlnxkZhaGoNNvYAICDPGGPX3qTJAs7ewZo+hw4+hXUlT4TUQRwMDQXzOWDkQS5K/Qq4HhEuIV1A9ss/yEhhqC5uRnZ2dkYP348AKCoqAj19fWWNgSAapMUjXIOIdEyWjAGXDRLHJlJtYggrh1oX+D3TcsQBEibefMO8R4BpRSmPpEhI3OOB7fII4os2y3MzFIC4ZLPfP/9B3CXlYBbGHsUSkhBNQ5i+8qUvJYsakVhxEjXaAUUkFtyxwBiv7nou1Ftu2pWFrjW9yQsT8MiaxjxJCGGwOv1IjMzU/qcmZmJL7/8UnFOVVUVqqqqAAAVFRVwu4Or+U6nU/E5nvQ2HYDvhSfEBTdnClz/vk5cYAMA9yL0/vJlnPvo/+Fs9V/87hQNGhsgjbCFAaQeOwxu/AR0c9zgriUZvHsc+LEZGFVyPZxT89B3cB9SLp4TlCWM7IOdN9DciFH7PwWXNgasp0s8b4H4TPLrepsOwLf7Q9HgOVMwdv4iceprMqdrDouF8ZgA9AtgNdvBdn2o/A2i4MziJeiWFVQbefnV4EelggFIOdUOVvNeaJu5FwELVMlX7mBbcWljwI4dRkp6OpzZ2br666B9K0Z65y+C7y+bpd/Mddv/0bz3YO+WmfIF7h9ou3/u/gi90npGP1KPHcb5qvY2Sw8o+lmY746UROosPSTEEDANZcepVtg8Hg88Ho/0uV3mX3O73YrP8UTYs1NU8P4kmVN7doKXuyjcE4Ab7wQ/pwhCXTX4g3shdLSp7qJ8/p72k8CmN3QbAQAQTh6HcPI4ug9/AVzmEfepdU8I64dU+Pp5Htyt94Zmyr7whOgH9++Bi5SUYGmGyZ1ge3aC6+wUp/blqyWXwynZcTNhk3NFXzcLzLgY0NcL31u/jmp0K3HJIuDaFmD7uwBj+OfH1f4qmv04p9EWg+KeADa5M7ifMM+j/54HceaSoRXJkH0rFtwTpN+Mu7BAXLjU6CuDvVsK+Xr/Cd9vnge//B5DfvfQPZj9cBzgcOLM5NyQrHez9IDUz9Af9rsjJVE6a+LEibrOS4ghyMzMREdHh/S5o6MDLpcrEaJETETT/13VEPrVswK//zkAY5ICioqB/iH3qWUtTWItmMBipzAA9s4rikXo0DBTMaxVWp8IEzGkZwpthDtBusfye8RaT7t2BPfsbQyfsKYXblQqWOC3EQZUv0dk2z8q2nJgAN2vVoIf6x46yioOUSwxuXPkCWIAcORLCM/9HPzDz8TunpPWV1TvwdRphhobvQEOZi8iW42EJJSNHTsWW7ZsQWFhIUaOHInf/e53+MEPfoD09PSw18QroWwotBLOtAgmk/iLcI2bAEz/F+DEt1DPCAxBGADzngQ3boIiKUYaaZ34ViUgU57PO8TENvVi9/FvxISgQHirLFlKkTADBtbdBfbxX4G+XnBTpym/P4qknUDSEWs9Bvbq8+J3NTaA/9FdGPu9Upz79ihw8oRCLvCO6BKV5IlVvMOfx+FvC44TF/W/v0LfPXmHuLeFXKmpkve00Nu3zGSwd0szQQwMjHcA3xzBiUjkAAAbFUlEQVSJKTlMan91ZNt354O/4nsRy6om0n6oN2FPL4nSWXoTyhJiCHieR3Z2Nl5++WW8//77uPzyy3HppZcOeo1VDAGgs5PIFYszBfx9j4GbdpFS2QCI3ij4jcvpnuChkyfAPqlRdHJJWaujigCg/QTY7o/AfB3gci7A2Gt/gH+OHguclyrWStL8Wh78v94SmpXK8UDrP4CTx8Vs3LEZ4KZOCzEWejNL5ZnSOPCpbJQu3mNMyXU4m5YuGi9/Viw3t8hvMBpERTzGJRmkIVtTpoS5Sy4DDnwmfp/DASxaAv5Hd+lXzL52sO5O4ESrOAFMGQH++8u1M3u15DBQAeklIFfKyPMGzyw+fzTY7o+CAwbeARw7IrZ5XTWQrr/N5UjZ6IH7Mf9786O7w7ZFRIYgTD+MF2QIwjBhwgR873vfw9KlS3HRRRcNeb6VDIEeAh37/Km56L/2h8pSFWplo+a8VOCSy4C2b0NH6N+ZDsyaD27xNcCemtDrmaAcfQZGpwHj8935wHmjgE6feK0gAF+3gH1Sg5TcC9F7ugdczgXBqKbQJxOzUtNGi8/kfx5425WRT/194BdeGVLCgP/+CrHMxhCjM+Evm4EjXwQeShyh+/3F/PdX4PzJU3D2vFTFCFos89Dgj6JiwP/sBTdztm6lqihJ8ffPxe/lOHBzLgW/oFjzGq0SDkLl40DrMb8RuRque8rRO/kC00saRItcrnO1HwwqV6CcChwOMeFxah7w1ZcA/G1+cC+4i/S3ufre/HfniffPygZ/w+CzooAe0FWuQqMfxrPtrW4IKLPYRLi80NLOinC+MOsC3ILF4H9cBtayNFji2h8eGvCXCu9t0c5ZYEwsZge/kgq4KBgD2IAYiup0ikYhUHfIvx7Q/WplsKDd/GJwE3PEuHn5GgYTNNck2J//S/kM/mS1sJmlgQS1/j7tLG1VpVXMmgfughlD728gz6Jmgu7EN8V36/DVB9pWnU+hyCUQAC4zSyrbEFKV1SLVOEN2xhtiLUTe5qylCexvfw26dAQh5kqsekJLA7kkrLNT1zqVHf3+kUCGIEGo4/MBiArMGawtJOUtLFRVOn1vi6jseYdmxVP2dYssSkidDSsufHKBcsoBIwNOpkAZ8MlHwG33wbH0JghZ2arM5ODiKQBxIVru2519abBcc0BhMtFLAsBvqGR7LfgNl6J9Fl4F9vEOaQ8I3j+rGrRN8/KDWdR+wxlN4ttQSkOzbQf6g+Ul5HWiZEYkkmqc8UxoGqomzmCycHn54G6917+tqiBuaG9iJVahdru/Lwrw/WUzuIVX6k7iM7OURrJDhiBBqJUNAM2d0QLnSolssg1ykHsh8OXB0AvqdoBxkGUdq+B5KRErYGSQNgbs9xsUs5RA2WUpM1k1O5Fe+MDoPrCoem2peH1LE4TnH5NqIwkfV4F/+BnRhRTImeA4zWQ6Li8f/ENrIlaGiizqjpOaiW96lOxgSiOYbR042Z8tXrcDbKBf/DxrXojx4vLyxYS8mvfFA7KigYr7xzmhSd4X1TVx9Mgib3MjdswLB2tp8hsc/6Cjrzdkb5BEF3ZL1oxkMgQJJGSKPUh1UyB0Co/mv8vvBmmULQwEX5CQvWTFqpMAghu1+Gvzj+w+hTPb3g6eOyVPWYDtx2UQpuRJO3ahpytYfkIjW1UKowwQUHyB4nkGbKwy2HWspUncT9fgypeKETTvELdl5ABWuz3ohvv8U+DaH4Zeq1GNU00iSlUE2kxdEyecLGqFF+KiUynnaBSk+pqQEi48L+6ENjXYJxOpfJM5I5kMgUXQXd1UcidxwXhujhMXgQ/u09w8B2ljwL5uEQeufreTVocd/b/uw5lz54JrAn/9E4SqP0v5A8Eyxv1gXzaK5SzkG/aokrpCqqX6916Ol7827PpEjEpW676spSm4aRCAcAXd9Dy7VYqyhZNlKIWnNduNVEFqfoe8P/E8Rt/zIE4Dij4pz42JN1arNRUJZAgsgp6XX/6CSUXEAkr42h8C1/4wtP6NRkccTBmGTaySyhgHr0NP19AKvagE6DoFjBkrZj/L5IrlJeltOgBhz86Iq58apWQVs45A3Sn5+sQgBd2GenYrLWxGa0zlzxiN8dVSqvLif9yFBUhdsAg9b/46KuVrhgvHSgY8UsgQWAS9L7/CnSTzyw6m+EPuobXN4ntb0Dt/kfJv4ADO7+pwOMFdchmYaqMZLaWm3r8ABk+TWUuTWPOmry/iKbiRSlZr1Ipb7xU3sknPiOg+4X7HwIJ8oo3BUMZ00AXlKBRkuGuMMOxmuXCsZMAjhWNahX8syLffBjNjE1lrKFKsKqvWPsZISQFfvhoAZIq8X1pX4BdfE3a/ZUn5d50SE8Dk0Uw8D+7yJUBGlmZJ6aFeGvV5wntbwP70tugv5nlw/3qrrj2IjUZ4bwvYtrfF0aj/GcVIJ78rzOEE//AzyFqwKGwfUAcASLWdwhyXX2eGwtHbX+XfD2BQWaOVd6hrArJGeu+Q3y0O/YdqDREJJ9yIM2Sa748h55feBO7QAbCBgeDipz+yR5ruq7dAfOc3mqGs4kVcSBin+h6D1kkazFecyI1NVKNRxqC5OB5SpVRG2MXYQfzNVliUjNT1E40rUO81kd47mV04ZkGGYJgz5MJemBjywab/8J4M3QJRayvMADm5wNctCkUBQHEPoa4anMaoLpyv2PXv68SqpwmcgmsuitbJZgT+xXE1CsMczgUyiLKy2qKkGYrVzDDMZHbhmAUZgmHOUEpD/lLIY8jVLwv75mspkQcOB8BzgCDGz3OXXAZ26H+CMwLeAQRKDjhTRJfJN68q1iTwj5Zg4pUUgx9a3VRKwhLEekYBJTMiv8C4Es0xoB6N8g+tCSbQTc0TM2DT06Xd1LQMs5ZSGkxZWW1EO5is0bqEzJ7xxBqsEC1WzTMgQzDM0RuNxGnEkMsjYxSJPIHtMv0JYdykqeAefkbKrA2EqCrCK/0L24o1CYcDuGyJPwZ/kB3P5MlnFkcRTRTYTe0vm8GXrw7r8uGX3qSpFMIpKyuOaMMGDagX0zF06QmrzXhiIZpyGImADMEwRD3qiFVphCTygAt+7hc3sucWXgUuMytkVCtdobUm4a/Fw11YEDbJih06EHQ7hcnEtSLh6vcYHb5qZUIU+q5qKWly0LpAFpvxRIt6MBBJOYx4Q4ZgmCEv68CcKeAfWhOz0lAszgZG5UJwAxu286+6a/poveQBYyXV6hni/GQg7NqLBUfzZqG9mD60IoxHG8XDRRMyGLBYOQw5CStDHSnJVoY6QLxlZe9tBr7yl2/2+/O5WfN0XRtOVql8dlY2kDVBXPhV77IWqPPOBM0NckLuFSgdDYglhP95Dnj/j8BXXyj2ElCU7p5bJJaa5h04f/KUiNtVV7lig5DL7frx/0bv5AuUf0vAngN6MLK/hmy0k3OB7lLQetooWlnjVg5cvSfJzXeBK7oqrhsPURlqmxKyPYFBWSIK3/eu6tAaRvIv1Ll1JPvm6+B6ARdwN4nbO7Lfb4AAsaCZVshq7y9flhZg9ZCIkMtway9mYNVFyJDFdAvMhuK1BjFYIIZuWeP0u5pmCDZv3owdO3ZgzBixxPCKFSswd+5cAMDWrVtRXV0NnuexcuVKzJ492ywxbAdfdBWEj6vEhV3eIdYcMhB552ZffQk07Pb/gRNnDG3HAX+Zaq2QUEWiFB+ol+TfIF6+nTNjYG8H91VWv7x9B/cBxREYgiRZgBRqt0sF1AKlvIfCCnkFepEbhkQZr0gzo2P6rhgGA/H8XU2dEVx33XW44YYbFMeOHTuGuro6vPDCC/D5fFi9ejVeeukl8Dxvpii2gcvLB//wM6a+YIrZwcG9wSqcE6eIO5UJA2FDQhUKWeDE6CO/qwAl3w/ZBCfcImvKxXNwJoKXNxnWGoTa7WBvrQcAsMYGaUY0FMli5ORIW5EODIA5HOAffiai8NJAJE4ks8IAmvkfFjSk8fxd4+4aqq+vR1FREVJSUjBu3DhkZ2ejubkZM2bMiLcow5ahFoeNGv1IL1Rgn4L99cFyEtAOCQ1RyMvvUZSqkDbBURVui/XlTYZFWvbZx8rPf/sA0GEIksHIqRECUWKANHvUs1tbuLDcSIk0MzoRxPN3NdUQbN++HbW1tcjNzcXtt9+OtLQ0eL1eTJ8+XTonIyMDXq/XTDEIGUZPN6VRviCILxKDVFNIKyR0KIWs3uREKwQVAPpq3ov45bV6yCV3yWVgjQ3BA/84jIG3NiiqtmpelwRGTo16K1K9KSLqSJxwGekRyWJRQxrP3zWmonOrV6/GqVOnQo4vX74c06dPl9YHNm3aBJ/Ph7KyMrz22muYMWMGFi9eDAD49a9/jTlz5uDSSy9V3KOqqgpVVVUAgIqKCvT2BrcFdDqd6A+3+5bFsJqsp//4Jnre+a1UsC3tllU4/4e3A4he1t6mA/D94n5xAdmZAte/r8OI/AL0Nh1A38F9SLl4jrhvr4776D1/oLkR7Y/dF/KdViSSdvVVPILeT2qDBzgOSBkRt+eLV3/tbToA35M/kQrVuX75H7r7iNTXeIe4rjQwEHMfiLSvRkqi9MCIESN0nRfTjOCJJ57QdV5JSQmeffZZAEBmZiY6Ojqkv3m9XmRkhJbs9Xg88Hg80md55T6rVvTUwmqyssm5oj8e4ujnzORcnPXLF7Ws7gngy1dLI5cu9wRxYcw9ASiegLPAkAtlQ1XbDPnKaTO1v9OCRNKu7Mrrgb27g9t/Mgb09+HUnp1xKakRt/7qngD+wTWR/36yvjbyTDfOffBnaXYQUxtF0Fejur1dq4/6fD64XC4AwJ49e5CTkwMAKCwsxLp163D99dfD5/OhtbUV06ZNM0sMQoVZ081Y3S7RLIxZ3dUTDSHrLv79oY1yV1gpzDTa3y9w3aj2Vpyrfs9yLp1kxDRD8Pvf/x5HjhwBx3HIysrCqlWrAAA5OTlYuHAhysvLwfM87rrrLooYijNWVKCDVtu0kPKKB1JU1sKrDH3uZAoz1cOI/IKkWxuxKqYZgp/+9Kdh/1ZaWorS0lKzvppIQsLNVIab8ooEow12MoaZDoUVBzXJCGUWE5ZBs4LlMFReicKq0TFE4iFDQFgaSXn19wHgxL0MiKiIZX3Ibu45u0GGgLA0XF6+fytMcVMc9odXwSZNBdzht38kwhONK8XO7jm7QKu0hPXp6RLLXjOm3OqSkGAtTeKm7C1Nxt9bwz1HDC9oRkBYHoV7iLOmeyiRrhOzR+y0tjD8IUNAmIKRilFyD739CjAguod6Z86KquCYGSTadWL2gnoylrAgIoMMAWE4ZihG9nVLcL+C/j6c3vQ62LU3RuXvNlqhJTqyKdIRezRtQGGawxsyBIThGK0YWUsTULcDUpkyxtC7vx442BCRkTFr5J5o10kkI/ZEz14Ia0KGgDAcoxWjYgN76SCL2MiYNXK3gutE74g90bMXwpqQISAMx2jFqDAsvCO4R3KkRiZtjLgrmsAZPnJPFtdJomcvVseu+RJkCAhTMFIxam1Kk3rsMM5Mzo3ILcT+8Kq4NSbPg1t+j61e9ABWmL1YFTu7zcgQEEmB2rCcv2CRVD5bD5JLBEz8r6crIaM/K4w4k2X2Em/s7DYjQ0DYArVLBGljFKM/9ZaZZmDnEadVGMwQ29ltRoaAsAVql4h6y0OxhAUDcziAopIht4eMBjuPOK3AUIbYzm4zMgSEbVC7RKTRHzj/nssM6BeA2u0QdlVThu4wQ48htqvbjAwBYUvkoz+kjREXkgPbQyLy0NRIv9NuI04rQIY4PDEZgl27dmHLli345ptv8MwzzyAvL0/629atW1FdXQ2e57Fy5UrMnj0bANDQ0ICNGzdCEASUlJRg2bJlsT0BQUSJfPTHJk01bXvIcN9JxBcyxOGJyRDk5OTgoYcewm9/+1vF8WPHjqGurg4vvPACfD4fVq9ejZdeegkA8Prrr+Pxxx9HZmYmHn30URQWFmLy5MmxiEEQMWPW9pCEtSBDrE1MhiCcAq+vr0dRURFSUlIwbtw4ZGdno7m5GQCQnZ2N8ePHAwCKiopQX19PhoCwDKQoCDtiyn4EXq8XmZmZ0ueMjAx4vd6Q45mZmfB6vWaIQBAEQehkyBnB6tWrcerUqZDjy5cvx7x58zSvYYzpPs5xnOa5VVVVqKqqAgBUVFTA7XYHhXY6FZ+tDMlqDiSrOZCs5mB1WYc0BE888UTEN83MzERHR4f02ev1IiMjAwAUxzs6OuByuTTv4fF44PF4pM/tsixSt9ut+GxlSFZzIFnNgWQ1h0TJOnHiRF3nmeIaKiwsRF1dHfr6+tDW1obW1lZMmzYNeXl5aG1tRVtbG/r7+1FXV4fCwkIzRCAIgiB0EtNi8Z49e/DGG2+gq6sLFRUV+M53voPHHnsMOTk5WLhwIcrLy8HzPO666y7wvGhz7rzzTqxZswaCIODKK69ETk6OIQ9CEEZihZpABBEvYjIE8+fPx/z58zX/VlpaitLS0pDjc+fOxdy5c2P5WoIwFaoJRNgNU1xDBJHMaJUiIIjhDBkCglDBXVgAOJ0Az1MpAsIWUK0hglBBpQgIu0GGgCA0oAxjwk6Qa4gghhmspQnCe1vAWpoSLQqRJNCMgCCGERTxREQDzQgIYhhBEU9ENJAhIIhhBEU8EdFAriGCGEZQxBMRDWQICGKYQRFPRKSQa4ggCMLmkCEgCIKwOWQICIIgbA4ZAoIgCJtDhoAgCMLmkCEgCIKwORwLt9M8QRAEYQuSckbwyCOPJFoE3ZCs5kCymgPJag5WlzUpDQFBEARhHGQICIIgbI7jqaeeeirRQkRDbm5uokXQDclqDiSrOZCs5mBlWWmxmCAIwuaQa4ggCMLmWLL66K5du7BlyxZ88803eOaZZ5CXlwcAaGtrwwMPPICJEycCAKZPn45Vq1YBAA4fPoz169ejt7cXc+bMwcqVK8FxHHp6erB27VqcPHkSWVlZeOCBB5CWlma6rACwdetWVFdXg+d5rFy5ErNnzwYANDQ0YOPGjRAEASUlJVi2bJn0fC+++CJ6enpwwQUX4Kc//SmcTnN+os2bN2PHjh0YM2YMAGDFihWYO3duVHLHG6vIIee+++7DeeedB57n4XA4UFFREbbvMcawceNG7Nu3DyNHjkRZWZmpboMNGzZg7969SE9PR2VlJQBEJdtHH32Ed999FwBQWlqKK664Ii6yWrGvtre3Y/369Th16hQ4joPH48HSpUst265DwizI0aNH2TfffMN+8YtfsObmZun4iRMnWHl5ueY1jzzyCDt06BATBIGtWbOG7d27lzHG2FtvvcW2bt3KGGNs69at7K233oqLrEePHmUPPfQQ6+3tZSdOnGA/+clP2MDAABsYGGA/+clP2PHjx1lfXx976KGH2NGjRxljjFVWVrKdO3cyxhj7zW9+w7Zv326orHI2bdrE/vSnP2k+T6RyxxOryKGmrKyMdXZ2Ko6F63ufffYZW7NmDRMEgR06dIg9+uijpsp28OBB1tLSonh3IpWtu7ub3Xfffay7u1vx73jIasW+6vV6WUtLC2OMsTNnzrD777+fHT161LLtOhSWdA1NnjxZGvXrwefz4ezZs5gxYwY4jsPixYtRX18PAKivr0dxcTEAoLi4WDputqz19fUoKipCSkoKxo0bh+zsbDQ3N6O5uRnZ2dkYP348nE4nioqKUF9fD8YYDh48iEsvvRQAcMUVVxguqx4ilTveWEUOPYTre59++ikWL14MjuMwY8YMnD59Gj6fzzQ5Zs6cGTILjlS2hoYGzJo1C2lpaUhLS8OsWbPQ0NAQF1nDkci+6nK5pBH9qFGjMGnSJHi9Xsu261BY0jU0GG1tbfi3f/s3jBo1CsuXL8dFF10Er9eLzMxM6ZzMzEx4vV4AQGdnJ1wuFwDxx+vq6oqLnF6vF9OnT5c+Z2RkSDKpZf3yyy/R3d2N1NRUOByOkPPNYvv27aitrUVubi5uv/12pKWlRSx3vNH6rRMhhxZr1qwBAFx99dXweDxh+57X64Xb7ZauC/TXwLnxIFLZ1O0ej/4px8p9ta2tDV999RWmTZuWdO0aIGGGYPXq1Th16lTI8eXLl2PevHma17hcLmzYsAGjR4/G4cOH8dxzz6GyshLM5MCnaGQNJ5PWcY7jYhMwDIPJvWTJEtx4440AgE2bNuHNN99EWVmZJeQeDKvIoWb16tXIyMhAZ2cnnn766UFntFZ9BiAy2eIls5X76rlz51BZWYk77rgDqampYc+zYrvKSZgheOKJJyK+JiUlBSkpKQDEmNzx48ejtbUVmZmZ6OjokM7r6OhARkYGACA9PR0+nw8ulws+n09acDJbVrVMXq9Xkkktq8vlwujRo3HmzBkMDAzA4XAozo8WvXKXlJTg2WefjUrueKP1WydCDjXy/jZv3jw0NzeH7XuZmZlob2+Xrk3EM0QqW0ZGBhobG6XjXq8XM2fOjIusY8eOlf5tpb7a39+PyspKXH755ViwYAGA5GpXOZZcIwhHV1cXBEEAAJw4cQKtra0YP348XC4XRo0ahS+++AKMMdTW1qKwsBAAUFhYiJqaGgBATU1N2BG80RQWFqKurg59fX1oa2tDa2srpk2bhry8PLS2tqKtrQ39/f2oq6tDYWEhOI7DxRdfjN27dwMQIwkCz2AGcp/0nj17kJOTE5Xc8cYqcsg5d+4czp49K/37888/x5QpU8L2vcLCQtTW1oIxhi+++AKpqalxNwSRyjZ79mzs378fPT096Onpwf79+6UIHbOxYl9ljOGVV17BpEmTcP3110vHk6ld5VgyoWzPnj1444030NXVhfPPPx/f+c538Nhjj2H37t3YvHkzHA4HeJ7HTTfdJP3ALS0t2LBhA3p7ezF79mzceeed4DgO3d3dWLt2Ldrb2+F2u1FeXm5o+Gg4WQHg3XffxYcffgie53HHHXdgzpw5AIC9e/fiP//zPyEIAq688kqUlpYCEI2bOnw0MAMympdffhlHjhwBx3HIysrCqlWrJGUUqdzxxipyBDhx4gSef/55AMDAwAAWLVqE0tLSsH2PMYbXX38d+/fvx4gRI1BWVqYIOzaaF198EY2Njeju7kZ6ejpuvvlmzJs3L2LZqqursXXrVgBimOOVV14ZF1kPHjxoub7a1NSEJ598ElOmTJFcOStWrMD06dMt2a5DYUlDQBAEQcSPpHINEQRBEMZDhoAgCMLmkCEgCIKwOWQICIIgbA4ZAoIgCJtDhoAgCMLmkCEgCIKwOWQICIIgbM7/BwFvVWkTE02aAAAAAElFTkSuQmCC\n", "text/plain": ["<Figure size 432x288 with 1 Axes>"]}, "metadata": {}, "output_type": "display_data"}], "source": ["import matplotlib.pyplot as plt\n", "plt.plot(out[:,0], out[:,1], \".\");"]}, {"cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[ 1.79262233e-03,  1.14341275e-01,  3.89670355e-02,\n", "         9.92673415e-01, -1.35283460e-03, -1.33684138e-03,\n", "        -5.51538021e-04],\n", "       [ 1.33244815e-02,  9.45778439e-01,  2.98248416e-01,\n", "        -1.20752748e-01, -3.48258394e-02, -2.38516836e-02,\n", "        -3.24298106e-03]])"]}, "execution_count": 16, "metadata": {}, "output_type": "execute_result"}], "source": ["pca.components_ "]}, {"cell_type": "markdown", "metadata": {}, "source": ["## PCA with scikit-learn and normalization"]}, {"cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [{"data": {"text/plain": ["Pipeline(memory=None,\n", "     steps=[('norm', Normalizer(copy=True, norm='l2')), ('pca', PCA(copy=True, iterated_power='auto', n_components=2, random_state=None,\n", "  svd_solver='auto', tol=0.0, whiten=False))])"]}, "execution_count": 17, "metadata": {}, "output_type": "execute_result"}], "source": ["from sklearn.decomposition import PCA\n", "from sklearn.preprocessing import Normalizer\n", "from sklearn.pipeline import Pipeline\n", "\n", "normpca = Pipeline([('norm', Normalizer()), ('pca', PCA(n_components=2))])\n", "normpca.fit(X)"]}, {"cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[0.02731781, 0.00012872],\n", "       [0.03511968, 0.00666259],\n", "       [0.03247168, 0.00632048],\n", "       [0.0287677 , 0.0060517 ],\n", "       [0.02758449, 0.00325874]])"]}, "execution_count": 18, "metadata": {}, "output_type": "execute_result"}], "source": ["out = normpca.transform(X)\n", "out[:5]"]}, {"cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [{"data": {"text/plain": ["(array([0.86819249, 0.08034075]), 4.332607718595102e-06)"]}, "execution_count": 19, "metadata": {}, "output_type": "execute_result"}], "source": ["normpca.named_steps['pca'].explained_variance_ratio_, normpca.named_steps['pca'].noise_variance_"]}, {"cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [{"data": {"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYYAAAD8CAYAAABzTgP2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztnX10VNW997/nzCTSMECTCRB5EWFUUjEUTFSIVCCkarV6U3pZvt/V6rIPleoCtU9rxS6r5cpaVVgt2GpdltXHRkWvgu19vKIxBYSABDWIYnxIoAiCAkmQhICZmbOfP87MmbNnzkwmc+Y9389aLDIze+bsfSbZ371/b1sRQggQQgghAdRMd4AQQkh2QWEghBAiQWEghBAiQWEghBAiQWEghBAiQWEghBAiQWEghBAiQWEghBAiQWEghBAiQWEghBAi4cx0BxLl8OHDme5C3JSWluL48eOZ7kbSyLfxAPk3pnwbD5B/Y8rEeMaMGRNXO+4YCCGESFAYCCGESFAYCCGESFAYCCGESFAYCCGESFAYCCGESFAYCCEkBxDtrdBefxmivTXl18rZPAZCCBksiPZWaE8sBXw+CKcT6n2/heIpT9n1uGMghJAsR3y6G/D5AKEBfp/+OIVQGAghJMtRJlcATiegqoDDqT9OITQlEUJIlqN4yqHe91uIT3dDmVyRUjMSQGEghJCcQPGUp1wQgtCURAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRILCQAghRCIpR3u2tLRgzZo10DQN8+bNQ11dnfS61+vF6tWrsW/fPgwbNgyLFy/GqFGj0NbWhqefftpot2DBAlx66aXJ6BIhhJAEsS0Mmqbh2WefxdKlS+F2u/HAAw+gqqoK48aNM9o0NjZi6NChWLVqFbZu3Yr6+nosWbIE48ePx/Lly+FwONDV1YWf//znqKyshMPhsNstQgghCWLblNTW1oaysjKMHj0aTqcT1dXVaG5ultrs3LkTc+bMAQDMmDEDH330EYQQOOusswwR8Hq9UBTFbncIIYTYxPaOobOzE26323jsdruxd+/eqG0cDgeKiorQ3d2N4cOHY+/evfjTn/6EY8eO4e677+ZugRBCMoxtYRBCRDwXvvKP1eb888/HihUrcOjQITz55JOYNm0aCgsLI9o3NDSgoaEBALB8+XKUlpba7XracDqdOdXf/si38QD5N6Z8Gw+Qf2PK5vHYFga3242Ojg7jcUdHB4qLiy3buN1u+P1+9Pb2wuVySW3GjRuHIUOG4ODBg/B4PBHXqa2tRW1trfH4+PHjdrueNkpLS3Oqv/2Rb+MB8m9M+TYeIP/GlInxjBkzJq52tn0MHo8HR44cwdGjR+Hz+dDU1ISqqiqpTWVlJTZu3AgA2L59O6ZMmQJFUXD06FH4/X4AwLFjx3D48GGMHDnSbpcIIYTYwPaOweFw4Pbbb8eyZcugaRrmzp2L8ePHY+3atfB4PKiqqkJNTQ1Wr16Nu+++Gy6XC4sXLwYAtLa2Yv369XA4HFBVFXfccQeGDx9ue1CEEEISRxFWDoAc4PDhw5nuQtxwC5z95NuY8m08QP6NKa9NSYQQQvILCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJCgMhhBAJZzI+pKWlBWvWrIGmaZg3bx7q6uqk171eL1avXo19+/Zh2LBhWLx4MUaNGoUPP/wQ9fX18Pl8cDqduO2223DRRRclo0uEEEISxPaOQdM0PPvss/jVr36FlStXYuvWrTh06JDUprGxEUOHDsWqVatw7bXXor6+HgAwbNgw/OIXv8ATTzyBRYsWYdWqVXa7QwghxCa2haGtrQ1lZWUYPXo0nE4nqqur0dzcLLXZuXMn5syZAwCYMWMGPvroIwghMHHiRJSUlAAAxo8fD6/XC6/Xa7dLhBBCbGBbGDo7O+F2u43HbrcbnZ2dUds4HA4UFRWhu7tbavPuu+9i4sSJKCgosNslQgghNrDtYxBCRDynKMqA2hw8eBD19fV48MEHo16noaEBDQ0NAIDly5ejtLQ00S6nHafTmVP97Y98Gw+Qf2PKt/EA+TembB6PbWFwu93o6OgwHnd0dKC4uNiyjdvtht/vR29vL1wul9H+8ccfx6JFi1BWVhb1OrW1taitrTUeHz9+3G7X00ZpaWlO9bc/8m08QP6NKd/GA+TfmDIxnjFjxsTVzrYpyePx4MiRIzh69Ch8Ph+amppQVVUltamsrMTGjRsBANu3b8eUKVOgKApOnTqF5cuX46abbkJ5ebndrhBCCEkCtncMDocDt99+O5YtWwZN0zB37lyMHz8ea9euhcfjQVVVFWpqarB69WrcfffdcLlcWLx4MQDgjTfewBdffIFXXnkFr7zyCgBg6dKlGDFihN1uEUIISRBFWDkAcoDDhw9nugtxwy1w9pNvY8q38QD5N6a8NiURQgjJLygMhBBCJCgMhBBCJCgMhBBCJCgMhBBCJCgMhBBCJCgMhBBCJCgMhBBCJCgMhBBCJCgMhBBCJCgMhBBCJCgMhBBCJCgMhBBCJCgMgwzR3grt9Zch2lsz3RVCSJZi+zwGkjuI9lZoTywFfD4IpxPqfb+F4uEBSYQQGe4YBhHi092AzwcIDfD79MeEEBIGhWEQoUyuAJxOQFUBh1N/TAghYdCUlCJEeyvEp7v1ybd0Vqa7AwBQPOVQ7/ut0S+akQghVlAYUkC4Lb/vkVVA6dmZ7hYAXRwoCLmHeaHB74+kGgpDCgi35Xs//gCYnR3CQHIPBg2QdEMfQwoIt+UXTJme6S6RHIZBAyTdcMeQAsJt+YXlFcDx45nuFslRlMkVEE4n4PcxaICkBQpDiqAtnyQLBg2QdENhICQH4EKDpBP6GHIAlrEghKQT7hiyHEakEELSDXcMWQ4jUggh6YbCkOWwjAUJkiqTIk2VJByakrIcRqQQIDkmRavs6fDPVW68E+g5yd+1QU5ShKGlpQVr1qyBpmmYN28e6urqpNe9Xi9Wr16Nffv2YdiwYVi8eDFGjRqF7u5urFixAm1tbZgzZw7uuOOOZHQn72BECrEyKQ7kdyKasEif6/NC1D8FCGG0CV6bQjG4sC0Mmqbh2WefxdKlS+F2u/HAAw+gqqoK48aNM9o0NjZi6NChWLVqFbZu3Yr6+nosWbIEBQUFuOGGG/DZZ5/h4MGDdrtCSN5iN8ktQli2NerPuYbrpkq/D1AUwK8BEIDfB62pEdjWyMCHQYhtYWhra0NZWRlGjx4NAKiurkZzc7MkDDt37sSCBQsAADNmzMBf/vIXCCEwZMgQlJeX44svvrDbDULyGrsmRUlYVAfElgZA0wCT+Qiu4RAvPhMSHwUQNnYpJHexLQydnZ1wu93GY7fbjb1790Zt43A4UFRUhO7ubgwfPtzu5fOaaBU14620mU0VObOpL7mKHZOiWVhExzHgnTeNCR89J6Feoy/cxNgJoXLxAERTI0txDEJsC4MQIuI5RVEG3KY/Ghoa0NDQAABYvnw5SktLB/T+TOJ0Ogfc377W3eha8ZBu93UWoPg3f0BheUXU5+N9fybGk8q+JItEvqNsxnI8pbOAy2bp38f2fwI+L+AswDcvnYXCYNtAmyB9j6yC9+MPUDBlesa/s0HxHWUJtoXB7Xajo6PDeNzR0YHi4mLLNm63G36/H729vXC5XAO6Tm1tLWpra43Hx3OoKF1paemA+6vt2AJ4vYZT8MSOLVBLz476fLzvz8R4UtmXZJHId5TNxBxP6dlQ733U2BmcLD07epHH0rOB2WfjNJDxQpCD6jtKEWPGjImrne08Bo/HgyNHjuDo0aPw+XxoampCVVWV1KayshIbN24EAGzfvh1TpkwZ8I5hsBEtf8F4XlEAKLrzcADvzwTx9oXx9OlD8ZRDvWYBzXrEEkVY2XkGyPvvv4+//vWv0DQNc+fOxfz587F27Vp4PB5UVVWhr68Pq1evxv79++FyubB48WLDWb1o0SL09vbC5/Nh6NChWLp0qeS4jsbhw4ftdjttJLoyiGaX1zZvgHj+qYDzsCBqWGHw/XANT2pseiLj6c/HYA6nRAYiYLgazX7ybUzZvGNIijBkgsEgDFaI9lZof38e2LMLgNBX4bOuNMIKwydVOxNutMk8Fb/Q2usvQ6yv181Nqgrl324xHKLpgJNO9pNvY8pmYWDmcxYR/6raC0Do5iSrsMJAjLoyuSLhxKiBZNomI+JosB9Gw6gtkk1QGLIA0d6qT+bvvAVoGoTDAfXn/xkxQYQm+YAofOvbUK+/GeLzA4CqAJoixagbJQ4SmHDjFZRkVX8dzKU/WEGXZBsUhjDSvXIzJgVvX+jJQNapI+z64atq9fqb9c948Rk9Y1VVgYsuBnY1yzHqCUy48a7gwwVEa2qEkmgS1iAt/RFLhLmTIJmAwmAimSs38x80SmdFbxecFMKwCtqyWlVrr78ceL/Q/w0vjtghJDLhxruClwREUYGmtyH8fq58B4BxD31e/YsPRJpxJ0EyBYXBhN1CZcbnhP1B9z2ySo8HR+QKUJoUgnEADieUmTWWnx0+yUfsIqprgOqapKwy4xEUs4Cg8xjE5jdt37/BhuIph3LjnXoBO78G8eIzRgZyMn4fCRkoFAYTyXKAhv9Bez/+AJh9dtQVoHrfbwM+hjf1EFT0HyhmFhirlX06J5CggIj2VpZQSJSek4GFgQiJwCB3yJPMQWEwYdcBKuUNmP6gC6ZMx2lE35EonnKIbY0BUQDg90Nsa4wZBRQuMOkM7YzGYHYg2yWqCMysAU6eAIZ/M7MdJIMKCkMYiTpAw/MFzAeeFJZXAMePx1wBhmeTxMouyWYTw2B1ICeFmTVQFBhmRCM0WQgACrRtjXH5GUR7K05t2gcxbhK/C5IQFIYEsDwJK2yyNlesDBJrRa1W10Db2gBofkB16L6CaNdzDQ+VxEiiiSGeCJhcmnRijScd0T5GGLLQv994d4DKzBo5NFlvFdciIPhZPUGfEx3WJAEoDAMkmp8g2m4gfCIN/tM2b4D29+ehVF4O9Yqr9Odv/l8Q722FUnm58cdsLn8hnAW6k/LFZ0IZwjfemZQ//HgiYGJNOtkWVhlrPFavxYocS/j6jz8YSEYEtK0NlrkpQJQdYHhQghLfIkD6LGTXbpLkDoNKGJIxecXyE0SEkgYm9R5NA1QHlFm1+mrw8wMQzz2pf96eFmgAlLET9Anf54PYuwdi7AT99fqn9F0EoJetfm9raCUphO60TALxmKeiTTrZGFYZMzfA4jVzqenkXd8bekLzR52kI0J+O48BQCjaawC1ruiwJslg0AhD0jJ0Y/zhme3ror1VntT9PohNb0Bs3gAMGyH3bcubwLQZkZMVEHJIA/oOofJyiL17ImLeEyWaw9xqQom6K8pCn0fM7ykNk6d+jYKQOKiOqNcJLiq0pkZg61sQmzZAbAnsMAJlTeJdzAQ/q+jQPvTmgLmPZCeDRxisVomIrEhq+d6wnUY8kTfi093ypG68IPQoEzMjSqJPugWByUVRody8EOoVV0EDImLe++2/ha07lsPccmUbZdLJxlVqrO8p1mvJMokpnnKo9y+Ly8cQbI+mRsAfWkhozzwOdB0HhIAIVNGNVxyGXjYLp/Oo4BxJL4NGGMInL7iGx7WDiOpT6O8P1DUccKj6H7qixAwzUiqqok5WlhOYVcx7LMdmFFt3PA7ziL5aTDrZGqYa63uyei3ZJrGBRmgpSlgGS8fR0M8+b1bsxMjgwPZBPblCcPJS/u0W3dnYc9LadBNGtJ1GLER7q6l+kQOOc86L0TEF4kC7XtoCiDg8RfFEHqgykEN4otm6B/o5/WHVz1xjoN91sg8WUmbWAM6C6A1smg0JiZdBs2MAIldw8Zg/EjGThOofCUBo8H/WFqVDejVUvb6QD0JVge/WQflGkbHylmouIWT6slqhW5lBYtm602FSySUG8l2nwuEumZ+2vBUyKwH6r1IcZkNCksGgEgYz8Zo/EjKTmPMMoASieIxPBKZdCpSNAz5r189C3vux/pLfD7zxCoSiyKGpPh+EqgCaHokkCgqkbGfR3qo7Lpve1rOmw0xesWzd6TCppJJkCthAvutUOdyNTPiZNfp3dqAdONCmmw6zxLFP8p9BKwxA/DbgeNqZo3uMPIOgOEiGYwF8uBPY/Z6+MrX8sMAkYISmaoBp8Wi2N4fKdgcO7wEiJpCB2rqzMcrIilSt2uP6nUixw90QiOAYs8ixT/KfQS0MyUKK7lEV3bcQnKSFRWSS5o98LuJDIYemhjuvg6WZzWYrAOZsaMsM7ThW2NkYZWRFJgUsfHcB6MeTDnTn0t/3ka2OfZLfUBhsYpzBHFyxa0roMAVFDQhFlJ1BLFwuqFdcBTF2gslMFEpsEy/8Wbc3hyVHBZPoAERm91o8l8uTUaYFLGJVHyxp0U/Yb5B4dzwD3fHZYTD6lkgkFAYbWJ6+JoQe5aNperjqxPOBvXvkNyrBYDCh/yy0yB1BTw9EeysUTzkcnnKI6hpoa58B9u/VX/d5IbY1Qr31roiVa/BsBMsIm+DuIyz80SrXIdsnhnQKWMy6S+adi8+r55gI0f9Z2Vlmsssl3xJJLRQGG2hNjYGdQvgLoSQlmEMZFVX3R2sB/8PosUDZWP0ozvAzGIQG7Y1X4Fj0oP5WTzkw3hMSBoS0JGLl6g1kRauq/rFB09LnB0JvEkI+KSzOuj7ZRjoEzHx/hLMA6v3LIsKHQ7s2kymxn8k+0zuecLJNqEjmGPTCkOjWWbS36uad8AlddQR2DIHENnPI4bkeYH8gdFUI4ItDwLEjFplNAXY1Q9u8IWSWmOCRmikTPHKfPt0d2r0IoW9IvnMllJk1oYQ2KHqfFcWosxQr14FmBehnZQTvj88Lbe0zUGZdKZmLzHWN9ByW/if7bDPZZZtQkcwxqIXBztZZfLpbnvSDTJhkTBridC/wxiuh1yZPBQ6062GnQfx+oNgNdHVYXESTzBKYWRPKojZN7AbhCVBCACUj5ZyGggI5+/v1lwN1kuRch1iZ4X2tu6Ht2JIVk1k6iEha378XIrBzE6oK5ZafGhVyARjHcsZV9C6LTHbZJlQkcwxuYbCzdQ7mKkRMGm0Qhw7ozt5Pd0OYVujKN4ogpl4CtLwrv8dKFADA4ZDMEooCPVkt2oouXChUVZ78e04ajlFjZRusk3TTT4DP2g0fQ7R7I9pb0bXiIcDrzTs7dLTdo3FWhlUQgaZB/O1PUuJZNk32AyWX+06Sx6AWhkS3zkbJC6Hpk3fVd4Adm0L1iwKOXbiGB+z8mr4ij1XSwFxPSVGAc8/Tdxhv/yPUv5k1xiEuVis6fUdQaBTdw3frApN/6BQwBJLjLOsk3XpXKFnuZJc+Ng2IrKTqzTs7dF/r7tAOSQ1FdxlJgj//z1B0mC/MryS0vLkPhACDXRgS3DpLp2sJAXSfkHcOgZW6+UAdzLtOf2yOYDI6osjCIITuZD50QF7hB8thhBW6M690zeMRn+6GiHYKmIUoivZWaL/7VWhl7HBIPgogrMRGHtmhvR9/YBJKTS+R3tQoZZA7POXQJngg/vYnOT9lgPeBIaEk2xnUwgAktnWOmFSlMxL08th6uQuTI/iz9kAimgnVAVyp10YSze8Ah/4lvx5YySuTKwZ0GplZOKRTwALXjFY+XHv9Zdlv4tckH0XwfhX/5g84kWc+hoIp0/WCguZ7ZbUjCjfXnXs+1AGcomd5ehzo5CfZxaAXhkSwmlTNDkcA0J5/OvQGhyMkHn6fLgjV84x8AW3zBuDQAYsLqcbKfyCnkUnF9268U/cdvPOWHiZr8qSGi6IyuQLC4QjtGJzWK+HC8gqopWcn4U4mFzsr8cLyitBhOU1v61FlFjuB8EXBQEQBsPi+tjVCNDUyd4BkFRSGBLCsYmrOJfj786GVt6LgGzXXoi94wE7wTOeAkIjPD+hnOlvEqyqzakMCEM0XYi7YZzYJmQ7gwcyakP8jhj3cbEtXFEgmpGwnGclZ5mTCqH4cT7le3DDsbO64rxEmLEKAuQMk60iKMLS0tGDNmjXQNA3z5s1DXV2d9LrX68Xq1auxb98+DBs2DIsXL8aoUaMAAOvWrUNjYyNUVcWPf/xjTJs2LRldShnxHTIfKI+hKICzAM6JF+D0c38MVT/9fx9DCKGv4NWwXIcgBYVGaYtovhDJCa6q+u4AkEt0GNFM8TnZg5NjrpHM5KxY5kXjnpvO5h7IdSxrLG1rZO4AySpsC4OmaXj22WexdOlSuN1uPPDAA6iqqsK4ceOMNo2NjRg6dChWrVqFrVu3or6+HkuWLMGhQ4fQ1NSEFStWoKurC48++ih+//vfQ1Wz9/wgaQLy9kF78RnDnCA5pRUF+Na3oVReju6//B7o64OxK/CZHJfR6uldXC2bhKyczubraRrE7p2hKKSgMJmimYI7gXykvwizaMebhtOfOSoUlRVZViSe90f0m7kDJAuxLQxtbW0oKyvD6NGjAQDV1dVobm6WhGHnzp1YsECf1GbMmIG//OUvEEKgubkZ1dXVKCgowKhRo1BWVoa2tjZccMEFdruVMpTJFfqBOv7A5P6vvdB+96vQwe1m+/P1N5uykaMf7WnJjk3QLpgSWp1amEekvggBfNgcipQKCJN6/c2GwGBbI4TPJ0Xb5DLxnsUdreQHACOs+JTwQ1McMe83AN10Z1FWxLhOP47lhI+KJSSN2BaGzs5OuN1u47Hb7cbevXujtnE4HCgqKkJ3dzc6Oztx/vnnG+1KSkrQ2dlpeZ2GhgY0NDQAAJYvX47S0lK7XU+M0ln4qvb7OLNhfeg5zY+iQ/sw9If/gb5HVsH78QcomDIdheUV6D1xHN0mh686sgzasS/6v44ACnY3o88fMo8UHdqHoZfNkvrSVVmNvh2bA+8JFPATAnAWoPi2n6KwXF85n9q0Dz2mzzrrgyaID5qgABgy53tGu3hwOp1pu/+9G9bj6+0bcdaMOSi6KmSi7GvdrSfaBeoXFf/mDyi8bBZgvj8BTm3ahx5zRJjmR+EHTTjz9v81xKInWFsq6KC3ut8ATgm/3jYgvkOFH0MD98LqHp/+5/9IffQekttYXSMZpPM7Shf5NqZsHo9tYRAWh9wrYfaKaG2sno9GbW0tamtrjcfHTYfRpxsxvRowTSpQHTilOHDq//xJXxnOvgan9U5C+zJYC0kfq2Y+4N2SQKZ0QQG8FZcAH7cA0HcgveMm4bRp3KK9Fdp7TaG3OgIZzIEQ15OlZwOB9mLcJL0MBvTy3Kcb/tuIPjrd8N8DKppXWlqalvuvbd4A8dyTAIC+lh3oOXUK6hVX6a/t2KL7UQIVTU/s2ALlq6+sdwzjJoVCUQFAdeDMl0fkRDXD5xMQ1rD7bT6ICabsc6mN+R47nDhz5kxkHydXSG3Cv9Nkka7vKJ3k25gyMZ4xY8bE1c62MLjdbnR0hEo6dHR0oLi42LKN2+2G3+9Hb28vXC5XxHs7OztRUlJit0spR/HIx2UqEzxRTRBSNnL4MZ/m+kTGhytQrrgqVPguRt0d8elufTILvA/V84yJ06rPwWgaFJ4ll+UIFM3LNnOGeG9r5OPA+MLNdkZtJ683UL9ooXwvqucBJ08Aw7+pl/zY1hhm3NMDBazOUgiP8op23oKVY1k0yY5l+hRILmBbGDweD44cOYKjR4+ipKQETU1NuOeee6Q2lZWV2LhxIy644AJs374dU6ZMgaIoqKqqwh/+8Ad8//vfR1dXF44cOYLzzjvPbpfSgtkurL3+slSPX2tqNCJ7zAlh4nQvsGEdQk5oi5LdCCt8F8P+HOHTqK4JfUqY/d0cTQOHQ/8XjIYKJL1limgOW6Xycog9LdJj4+ewCVb35QQPS/JDPP8UxNgJACBN6lIU2da3A6fuqXB6LoB/Ro2lsFqWDwkLBDD3S6qzZCEC9CmQbMe2MDgcDtx+++1YtmwZNE3D3LlzMX78eKxduxYejwdVVVWoqanB6tWrcffdd8PlcmHx4sUAgPHjx2PmzJm49957oaoq7rjjjqyOSIpGhBO46W0IU+RLYXmFbuZ4Yqm8Y4j4oMCqtZ9JWtu8wYijjxbGGu7glCY3TS/HLb7qAk506mUvMjRRxQr/VcNyP8In7fAJVgTLnQN6lJZxOFFkGKu069vSAF/7p8C/2iPCT0V7q37oUdjZFubXozm8o0WTEZLtJCWP4eKLL8bFF18sPXfDDTcYPxcWFuLee++1fO/8+fMxf/78ZHQjYyieciizaiE2vaE/YWGaiTjzwHhzdBOGFWa7u9jTAty2KEYYq2kyDDe9nOMBguavz58ZcDx+sugv/0C94irDfBQLxVMO5ZaFerKgpkkCGy2M1Qgx1rSAYIdllZtNSA4HcPmVUqirtnmDcT2hOuRjVX/3Kz1nxeHImUOPCAnCzOckocysibAnS68buwpT4sK0y6BMvGBgBfxi2N2la4VNhpamlyzIuJX6qjogOo4ZR5rGi7Rq//ljESv4WGGsouNYwKwWqD5rDj8N32UpCO1CAP2sDNNpfcHCe5gyPVRWxO+TTIuE5AKOhx9++OFMdyIRuru7M90FCaWkFEr5VKB0NNTrbpImoKKiIpweUgQMLwY+el9/sqAQ6n/8DOpls4Gu4xDb/6nb+kv6CV/z9un5CsHrXv1DKBNkv0y0viglpVDOn6JfQ3VAvLsRgB6Bo153U//XNo2nt7c3rrb9YfRVVYHP9gEH2iDe3QilfKrUH9HeanmPjFX9nl36+6proF42W2ojjTv8ff8KhFYr0M2Ae1pC1zbfI9Wh9++TDyHe3ah7if7VFjkgoQHf+AZwIhR2rZx7HpSplyTlfsVLMr+jbCHfxpSJ8QwbNiyudtycpO+gAAAWSUlEQVQxJJH+nIrqFVdFRBkNtMZPLLt7uL3bMnM3SonujJiRTH1ByciQSSe8UGCsMiQJ7nzEtkZTCfTggUsi0g8RuEei4xjwzpvGdRQF+ml4VmXUJ08FDv4rFMU0syayDSFZDIUhzUQ4TAc4sYn2Vj0qJpDRbH6+P4Hpr0R3Ognvi3LjnXqugYUpLtY9ilUKQ8o9CPhvAF0UxDtvhTqjqvo/v18/4KjzmHTWtnrNAr2/5ppGwTIjb7wih/4KAbz9DymfhP4FkmtQGDLMQE6Ri7lyNq+AowhMvyW60ziJWYaARtnBxLpHsQoMGgUNA6fXCafTyGqWqKhC8Q23o+t/1ukRZZs3hN4TOPEu2nWUiRdAtOyAVPKkn5BWQrIdCkMcpHLiHEjCU8yJfUtDqKHqMM56liaxaKe22SxXnQjRHORRy4HHuEdW75MKDOrPRMkbAbD7PeCG26G4R0L4/fJ7wkxLEdeeXKGblIICpCj9Cjwh2Q6FIUCsePRkTpzRznIAQhEvUX0DwQiawDnMwclfdByTM6AvutgyE9tqgtWT80LVQrWmRigp3j1IhwjFaW4ZaFKYITzmE9mi4ffpFXBnzJVPcYtjkjffU7PJiuYjkstQGJAa5+ZArmNpbzdNMFI8vapC+c6VwDkeOZPZlICF4cXRbfLhE2x4tdCtDRCalrLdQ8RYZtX2/6YEMCbsbY16It9XnXoUURSR8O3dA/yrLXTi3VddwPBiI5s8fPdlHk/UkiU825nkKBQGJO7ctHUdbx+0tc9AveFO+XmfV4+PF8I6a1lAj+DpORmRyYySkZGHvwScqVFzA3pOwijcB0jx96nIbZDvtWbE/iciQvFMvMFjM3UHsyPSv2DG74M40A5sazQiikSMOljxHdrEIztJ7kFhQGLOzXgITlx9l84CSs+OPMth/15ojz+oR7AEr68ogddN9u3g+zRhnAMNhGX0hh3DaT6/WGx+U68LZDpnWhp7QUFIRBQl6nnHycDSxJOACMUVhRUuqOd4QnkLgC4Wnm8B7Z+Eyl0ogDAvEt7bOuDztvt7jZBsh8KAxJybQfr1TXi96Pr781BuXghl7ISQLyCIT49gMSqfnuMB3v5HpEgZpbsFxLZGKDNrIo+INJk7FE85lE93B5ypmn5q3OYN0LbJq/PwsQNIqfkjeL2gaCUqQuG7rPDsYiufjPKdKyE+P6CLkqICV9bp91oTgEM/GlUZO0HOYK+8HGLvHutFQ6wFRRRHP01LJBegMARIpOJlLN+AVO3TH6j2efl3I23cTt2JbPgL9u6J8DFor78cKqXh9+vml61vQ71/WSjG3upUsODkFOxHWJRNtLGnetJSPPq50qK6JrICbJwTZ0Thwq0N8AMhn0CYTya8jDlcwyHeeTMU4qtpetLg2AkRi4Ropc9jLSisBJemJZIrUBhsEOEbeP4pQBMhkQir9mlkywbPZvj2pVCvnt9vWefQBG/KsvV59Z2Dpzym2UKprtEdqbvf01/PolDKoCCJ9lZof/ujHnIb5viOWpLbE1a40O/Td0RNbwPjzg2JYcAnEx4BZj7qE4AuLp/sgrZ3T0TiX6xFQ7yvSaXZaVoiWQ6FwQaSuQCKdCwkek6aqn0KozSCMrMmYqITnx/QTUWwDo80zC8vPiPZyIObj9j5CV5doL5bB+UbRWkJQxWne4HP2vst2RF8Tje5mUQvMHEC/ayyz/GEjuPUP00fr1H/yPp+Rq10K6x3VMkgmUEMhKQaCoMNwmPYxYvPRCRsibETUHRoH3rHTZJWraK9VV9FBt8nAuU7p0yPfq0b79TLOWt+QHVIB/NgZg0UBYbJRMpP8PuBN9dD+d+PpVQUwid4sacFGqDb7QPnHph3BACg/f15eZI2TeYxnbvtrRAv/DnSZ2Pm3POg3nBnpM9ocgWE+fQ81QE4AiUxrIQkCb4BO0EMhKQbCoNNzOYCK1u04inH0MtmRZ7VHLSBK4HjPoXQ/7XsgPbxB5Y2aMVTDvXn/xmzCF+wYFtEmW+hpdR8YUzi4c+/86a+IwrfEWxr1J28QZNPQBDCI6ciju8MONgtV/1njweOHg6N+bP9ln1VPKZDer7qgjKiGK5vTUXPl0dSmuCYiB+LkExAYUgSA1lVSithKIECboEQ1SgO4iDhk0u0VbXiKYdy88JAToQW18lwdrD0gwDAN0uAA+2mhroACIGAkARE4VvfjigMCMDYCZkT+oSqApMmR3biy88BTzmw9xP9c2OIoSE8TyyF8PnQve2fUO99NHbIa5LNTIxSItkKhSEJDHRVGWFvDmbbbmkYsIM4lu3aqsx3qjCbSsw+BmXsBGgff2AcxBPcEQCQqpX2Vy1WAaTEOOz9OLITmga0teqRXnGEwYYnHAad+dK4UuQbYAIcyWYoDElgoKvKaPZmYeGY7g87ORjJJtq1ovUvVh5GxKT96Ud6ToIvhk8BACCA6nl6f5R++ju5AiL4mUJAbGmACEsUTJVvgAlwJJuhMCSBRFaVVpNoopN4ttiupYN3ICfKic8P6I7mczxSdJTkJ/F6IVQVyi0LI7PEvzikm9zM5TscTuDyWqBoKPDm+pDJLFDGQvh8EE2NUl5IeL9QPQ/YvAGxTE+puL+MUiLZDIUhCcS7qsxnm7JkGnE4dEd6IAIJ864D3nhFb7inRfekFBSGchXMyYCaH6L+T8B3rgQqKuVDcLSATwbQtwOX18Jx21369addFiofEi2/JKxf6n2/hVpdIx/Ak6YJmlFKJJuhMCSJ/laV+W5Tlidj+TwDvN8U+Qazo9yoBRVKBhSb3gCcBbpfIvg8ENg1BM6pNoXrRjjlg/WYgFCUki8oLCEHv3rNAqj3/TYipDgdZMtOj5Bw1Ex3YLBgZVPOJ5TJFbrTV1X1/4OlwB1OYKJFBJFpda54yqHcslB/jxm/DzhnkukiCpTvfBfKv90SU1gVT3ko8zy8BEnwGM+w6w/94X9wkiYkAHcMaSLfbcpWtYGkxD9AN/9cOhvKmPGRjuhABJV0hrIQQPlU4PMDUavIRqXnpJ5xLnVSF5ZgeXIKASHWUBjSRHAVK97bqodxDmBSyhXfRLhpJJSBHUh8UxQoY8ZHPQtZ8ZRHnqH8Wbt+3w609xtlJH1WeHlvRdEd0/EKCyGDGApDihHtrTi1aR80xRFK0Nq7B2LshLgmqFzzTYSL2EB3ShFnKH+yC+L/fayHk2qa5aE+0Y5L5ZGbhCQGhSGFBCf1nmCRPREqshdv3HouxbtHE7Hg+QvxrPiN9n9/Hvhkl36/jFIbkfculnDSuUtIYtD5nEKMSV3T9IndwvHZH5JTN8t9EzEd7NsaITa/Ce2JpRDtrTE/R/GUQ73+5kBUUqQz23wPIg7s+fvz/X4+ISQ23DGkEMvSFwM0aeRSvHs0s1Eiu55ozuyIcxnCfQl7QmcqJHqvcsWnQ0iqsCUMPT09WLlyJY4dO4aRI0diyZIlcLlcEe02btyIV199FQAwf/58zJkzBwDwwgsvYPPmzejp6cFzzz1npytZSXBysxsjnwsmEWMytRC/qILRzwRs5cy2amOYnvbsQn9FCOMZRy75dAhJBbaEYf369aioqEBdXR3Wr1+P9evX49Zbb5Xa9PT04L/+67+wfPlyAMAvf/lLVFVVweVyobKyEldffTXuueceO93IaqzKbucb2uYNgexivSRF+GRqtetJdjlr9fqboUU5m3kg5JJPh5BUYcvH0NzcjNmzZwMAZs+ejebm5og2LS0tmDp1KlwuF1wuF6ZOnYqWlhYAwAUXXIDi4mI7XSAZoK91N7TXX9ZX/O2temlvvz/gKPZaJu8pnnKo1ywIOY2TnPAXFJ/+kt/6/Zwc8ukQkips7Ri++uorY2IvLi7GyZMnI9p0dnbC7XYbj0tKStDZ2WnnsiSDiPZWdK14SC9453RCqa6RT1FT1fiKCKYg4S8ZJrdc8ukQkir6FYZHH30UJ06ciHj+xhtvTPiiykAylQI0NDSgoaEBALB8+XKUlpYmfP1E6GvdDe/HH6BgynQUlg9sEnM6nWnvb6o4tWkfenxeY6U/ZMgQnC4s1A/oUVUMu/M+FF02q//7VToLfY+sSvieJhvpOyqdBVw2K6P9sUs+/c4FybcxZfN4+hWGhx56KOprI0aMQFdXF4qLi9HV1YXhw4dHtCkpKcGePXuMx52dnbjwwgsH3NHa2lrU1tYaj4+n0WYvHcUZxR4ey5FaWlqa1v6mEjFukh5G6vMCDie+nl4NdXq1MfZeTzlOvbtFul+Yd51xcI96xVXyvZp9DU4DQAruz0Cii/LpOwLybzxA/o0pE+MZM2ZMXO1smZKqqqqwadMm1NXVYdOmTbjkkksi2kybNg0vvPACenp6AAC7du3CzTffbOeyaac/h+RgimRRPOUo/s0fcGLHlohM4yDhuQXBkttiTwv8x74A3v5Hyu/VYPpOCEk2tpzPdXV1+PDDD3HPPffgww8/RF1dHQCgvb0dTz31FADA5XLhhz/8IR544AE88MAD+Pd//3cjpPVvf/sbFi5ciL6+PixcuBAvvfSSzeGkhv4cknYdqaK91XDmWj3ONgrLKyRHcjjS/Qrn/aa0VJnN92q2hKQSWzuGYcOG4de//nXE8x6PBx6Px3hcU1ODmpqaiHa33nprRHhrNtLv8Zk2HKkRZxvfeGeoplKOrnTN90uc7g0d0gMAF1frO4YUV5nN92q2hKQSZj7HSayIFzuRLBEr2/e25kUcvfl+aSPLjKqy6hVXSaetpWpsjC4iJHEoDEkimnAEq6uKKJnPESvbysshkpColU2oV1wFXHGV8Thdmdy5kDFOSDZCYUghUnVVh7VZyDIreOyErF3p9id0hJDch8KQQiQzEaKbhaxqAmXjpBuP0BFCch+W3U4h+VZeQSojzkgfQvIW7hhSSLKqq2YLjPQhZHBAYUgx8VRXzZX6//kmdIQQaygMGSbXMnRTWUY8mkDminASki9QGDIM6//rRBPIXBNOQvIBOp8zTL45qBMlWgkLlrYgJP1wx5BhUpWhm2vml2iObTq8CUk/FIYsINl5C7lofokmkCxtQUj6oTBkEcla5eeq3yKaQGZrwh8h+QqFIUtI5iqf5hdCiB0oDFlCMlf5NL8QQuxAYcgSkr3Kp/mFEJIoFIYsgat8Qki2QGHIIrjKJ4RkA0xwI4QQIkFhIIQQIkFhIIQQIkFhIIQQIkFhIIQQIkFhIIQQIqEIIUSmO0EIISR74I4hDfzyl7/MdBeSSr6NB8i/MeXbeID8G1M2j4fCQAghRILCQAghRMLx8MMPP5zpTgwGJk2alOkuJJV8Gw+Qf2PKt/EA+TembB0Pnc+EEEIkaEoihBAiweqqSaKnpwcrV67EsWPHMHLkSCxZsgQulyui3caNG/Hqq68CAObPn485c+bg66+/xooVK/Dll19CVVVUVlbilltuSfcQAAAtLS1Ys2YNNE3DvHnzUFdXJ73u9XqxevVq7Nu3D8OGDcPixYsxatQoAMC6devQ2NgIVVXx4x//GNOmTcvEECJIdEwffvgh6uvr4fP54HQ6cdttt+Giiy7K0ChC2PmOAOD48eNYsmQJFixYgOuvvz7d3Y/AzngOHDiAP//5zzh9+jQURcFjjz2GwsLCTAxDItEx+Xw+PPXUU9i/fz80TcMVV1yBH/zgB+kfgCBJ4bnnnhPr1q0TQgixbt068dxzz0W06e7uFosWLRLd3d3Sz2fOnBG7d+8WQgjh9XrFQw89JN5///209l8IIfx+v/jZz34mvvjiC+H1esX9998vDh48KLV54403xNNPPy2EEGLLli1ixYoVQgghDh48KO6//37R19cnvvzyS/Gzn/1M+P3+tI8hHDtj2rdvn+jo6BBCCHHgwAHxk5/8JL2dt8DOeIL87ne/E0888YR47bXX0tbvaNgZj8/nE/fdd5/Yv3+/EEKIkydP5vzv3DvvvCNWrlwphBDizJkz4q677hJffvllegcghKApKUk0Nzdj9uzZAIDZs2ejubk5ok1LSwumTp0Kl8sFl8uFqVOnoqWlBWeddZaxEnU6nZg4cSI6OjrS2n8AaGtrQ1lZGUaPHg2n04nq6uqIcezcuRNz5swBAMyYMQMfffQRhBBobm5GdXU1CgoKMGrUKJSVlaGtrS3tYwjHzpgmTpyIkpISAMD48ePh9Xrh9XrTPQQJO+MBgB07dmD06NEYN25curtuiZ3x7Nq1C+eccw7OPfdcAMCwYcOgqpmf0ux+R2fOnIHf70dfXx+cTieKiorSPQT6GJLFV199heLiYgBAcXExTp48GdGms7MTbrfbeFxSUoLOzk6pzalTp/Dee++hosLe0Z6JEN4/t9sd0T9zG4fDgaKiInR3d8c1tkxgZ0xm3n33XUycOBEFBQWp73QM7IznzJkzeO2117BgwYK09jkWdsZz5MgRKIqCZcuW4Re/+AVee+21tPY9GnbGNGPGDAwZMgQ/+clPcNddd+G6666zNEmnGvoYBsCjjz6KEydORDx/4403JvyZiqIYP/v9fvz+97/H9773PYwePTrhz0wUYRGgZu5frDZWz2cDdsYU5ODBg6ivr8eDDz6Y/A4OEDvjeemll3DttddiyJAhKevfQLEzHr/fj9bWVjz22GM466yz8Mgjj2DSpEkZWVSZsTOmtrY2qKqKp59+GqdOncKvf/1rVFRUpH0+oDAMgIceeijqayNGjEBXVxeKi4vR1dWF4cOHR7QpKSnBnj17jMednZ248MILjcdPP/00ysrKcO211ya343HidrslE1ZHR4exCwpv43a74ff70dvbC5fLFfHezs5OwwyTSeyMKdj+8ccfx6JFi1BWVpbWvlthZzxtbW149913UV9fj1OnTkFRFBQWFuLqq69O9zAi+hpkoL9zF154ofG3Nn36dOzfvz/jwmBnTFu2bMG0adPgdDoxYsQITJ48Ge3t7WkXBpqSkkRVVRU2bdoEANi0aRMuueSSiDbTpk3Drl270NPTg56eHuzatcuI3HnxxRfR29uLH/3oR+nstoTH48GRI0dw9OhR+Hw+NDU1oaqqSmpTWVmJjRs3AgC2b9+OKVOmQFEUVFVVoampCV6vF0ePHsWRI0dw3nnnZWAUMnbGdOrUKSxfvhw33XQTysuz4yxuO+N55JFH8OSTT+LJJ5/ENddcgx/84AcZFQXA3ni+/e1v47PPPsPXX38Nv9+PTz75JCt8J3bGVFpaavgbzpw5g71792Ls2LFpHwMT3JJEd3c3Vq5ciePHj6O0tBT33nsvXC4X2tvb8dZbb2HhwoUAgMbGRqxbtw6AHq46d+5cdHR04Kc//SnGjh0Lp1PfxF199dWYN29e2sfx/vvv469//Ss0TcPcuXMxf/58rF27Fh6PB1VVVejr68Pq1auxf/9+uFwuLF682FjNvPrqq/jnP/8JVVXxox/9CNOnT097/61IdEyvvPIK1q9fL+0Uli5dihEjRmRwNPa+oyAvvfQShgwZkhXhqnbGs3nzZqxfvx6KomD69Om49dZbMzwanUTHdObMGfzxj3/EoUOHIITA3LlzM/IdURgIIYRI0JRECCFEgsJACCFEgsJACCFEgsJACCFEgsJACCFEgsJACCFEgsJACCFEgsJACCFE4v8DsYc4JdkD3UAAAAAASUVORK5CYII=\n", "text/plain": ["<Figure size 432x288 with 1 Axes>"]}, "metadata": {}, "output_type": "display_data"}], "source": ["import matplotlib.pyplot as plt\n", "plt.plot(out[:,0], out[:,1], \".\");"]}, {"cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[ 0.00415209,  0.92648229,  0.11272098, -0.05732771, -0.09162071,\n", "        -0.34198745, -0.01646403],\n", "       [ 0.01671457,  0.0789351 ,  0.85881718, -0.06957932,  0.02998247,\n", "         0.49941847,  0.02763848]])"]}, "execution_count": 21, "metadata": {}, "output_type": "execute_result"}], "source": ["normpca.named_steps['pca'].components_"]}, {"cell_type": "markdown", "metadata": {}, "source": ["## PCA with statsmodels"]}, {"cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [], "source": ["from statsmodels.sandbox.tools import pca\n", "xred, fact, eva, eve = pca(X, keepdim=2, normalize=False)"]}, {"cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[536.44492922, -50.83312832],\n", "       [730.34140206, -79.13543921],\n", "       [470.9815846 , -75.4476722 ],\n", "       [466.40143367, -62.53420646],\n", "       [481.66788465, -55.78036021]])"]}, "execution_count": 23, "metadata": {}, "output_type": "execute_result"}], "source": ["fact[:5]"]}, {"cell_type": "code", "execution_count": 23, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([732151.6743476 ,   1513.97202164])"]}, "execution_count": 24, "metadata": {}, "output_type": "execute_result"}], "source": ["eva"]}, {"cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[ 1.79262233e-03, -1.33244815e-02],\n", "       [ 1.14341275e-01, -9.45778439e-01],\n", "       [ 3.89670355e-02, -2.98248416e-01],\n", "       [ 9.92673415e-01,  1.20752748e-01],\n", "       [-1.35283460e-03,  3.48258394e-02],\n", "       [-1.33684138e-03,  2.38516836e-02],\n", "       [-5.51538021e-04,  3.24298106e-03]])"]}, "execution_count": 25, "metadata": {}, "output_type": "execute_result"}], "source": ["eve"]}, {"cell_type": "code", "execution_count": 25, "metadata": {}, "outputs": [{"data": {"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAD8CAYAAACVZ8iyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztnX98FPW19z8zu4sCIXGziwk/IpYlCHrRhCYIEYLAqq3t7cOltU+0ffoUFNoXtt6C1LYi9EdE0xdGrL7EthZ53aI+V3hK8D6veouNqUQIXiIQyoUGSKIWbEogGwgRaHYz3+eP2Z2dmZ3Nzu7O7Mwm5/168SIzOztzZvY73/P9nnO+53CMMQaCIAhiWMNbLQBBEARhPaQMCIIgCFIGBEEQBCkDgiAIAqQMCIIgCJAyIAiCIEDKgCAIggApA4IgCAKkDAiCIAiQMiAIgiAAOK0WIBn+9re/SX97vV6cP3/eQmn0Q7IaT7bICZCsZkGyJmb8+PG6j6WZAUEQBEHKgCAIgiBlQBAEQYCUAUEQBAFSBgRBEARIGRAEQRAgZUAQQw7W3grhrR1g7a1Wi0JkEVm1zoAgiMFh7a0Qap8AQiEwpxP8o0+C802zWiwiC6CZAUEMIdiJo0AoBDABGAiJ2wShA1IGBDGE4G6aATidAM8DDqe4TRA6IDMRQQwhON808I8+CXbiKLibZpCJiNANKQOCGGJwvmmkBIikITMRQWhAETnEcINmBgShgiJyiOEIzQwIQgVF5BDDEVIGBKGCInKI4YghZqLNmzfj0KFDyMvLQ21tLQCgr68PmzZtwrlz5zB27FisWrUKOTk5YIxh69atOHz4MK655hqsXLkSkydPNkIMgjAEisghhiOGzAzuvPNOPP7444p9u3btwowZM/D8889jxowZ2LVrFwDg8OHD+Pvf/47nn38eK1aswG9+8xsjRCAIQ+F808Dfex8pAmLYYIgyuPnmm5GTk6PY19zcjPnz5wMA5s+fj+bmZgDABx98gMrKSnAch6lTp+LTTz9FT0+PEWIQBEEQKWJaNNHFixfhdrsBAG63G729vQCAQCAAr9crHefxeBAIBKRj5dTX16O+vh4AUFNTo/ie0+lUbNsZktV4skVOgGQ1C5LVWDIeWsoYi9nHcZzmsX6/H36/X9qWF5SmYtjmkC2yZoucAMlqFiRrYsaPH6/7WNOiifLy8iTzT09PD3JzcwGIMwH5Q+nu7tacFRAEQRCZwzRlUFZWhj179gAA9uzZg/Lycml/Y2MjGGM4efIkRo0aRcqAsA208pgYrhhiJnruuedw/PhxXLp0Cd/+9rfx1a9+FYsXL8amTZvQ0NAAr9eL1atXAwBKS0tx6NAhPPLIIxgxYgRWrlxphAgEkRasvRVsfwPY3npAEGjl8RCFtbdSyHAcDFEG3/ve9zT3r1+/PmYfx3F46KGHjLgsQaRFpGNATi7Yv78MBPujH4ZXHhvdYVBnZB2UZmRwKDcRMSyRdwzgODH1RASOM2XlMXVG1qKVZoSefxRSBsSwQpoNBM5FOwZwYuoJxgDeAVQsAl+x0PhZAXVGlsLdNAPM6QQGQpRmRANSBsSwQTEb4HnA4QAEiB1D1XKgr9dU8w11RtaiJ83IcDbjkTIghg2KkTkDcMfd4DxjM/biU84j6xms8M9wN+ORMiCyAiNGbNxNM8B4HhAYwPGmmIISymCDKmTDefQ7GMPdjEfKgDAFIzsc9YiNq1qOT9kA2MTJyZ+b40TfQJxV70Od4T76HYzhbsYjZUAYjtEdjmLEFuwHe/Ul9HEAHMmdm504CgwMiBvCgKEjv2wZbQ/30e9gDHczHikDwnCM7nC4m2aAORxAKBz+GbH5I7lzmzXys8toW49CGu6j30TYwYxnFaQMCMMxusPhfNOAikVA426EtUBKawHMGvnZYbStVyEN99EvER9SBoThmNHhcJN8ovOXCYDDiZGLvoB/lFYkfW4zRn52GG0no5BSeQbZYgYjUoeUAWEKRna6rL1VTBfBBIDnwd2/ArlLvmab9MV2GG2bqZDsYgYjzIWUAWF7oqNeJv7r67VapBistjWbqZDsYAYjzIeUAWF/cnIBngMEc3IGDRXMUkh2MIMR5kPKgLA1koloIGwiqlpOo9I0SMX2bwczGGE+pAwIWyOZKMDEfzY0EWUL6dj+rTaDEeZjWqUzgjAC7qYZgNMZTixHJop00LL9E0QE02cGDz/8MK699lrwPA+Hw4Gamhr09fVh06ZNOHfuHMaOHYtVq1YhJyfHbFEIm6NlwhjOJgqjwzmHqu2fwl6NISNmoh//+MfIzc2Vtnft2oUZM2Zg8eLF2LVrF3bt2oWvf/3rmRCFgD1fnsFMGMPNRGFWCc6hqFgp7NU4LDETNTc3Y/78+QCA+fPno7m52QoxhiWRl4ftek383yaF31MxYQzF4vXS77PnD+II3mCTDuebBv7e+2zRYRrx+5HpyzgyMjPYsGEDAOCuu+6C3+/HxYsX4Xa7AQButxu9veQUzBRmxoyrZxzJzECSNWH0tx4dkiNCtr9BWYvZpBKcVsPaWyE8sxYIBcGcLvBrNuhPOBhuV/2z5hpm+rLjbDnTmK4MqqurkZ+fj4sXL+LJJ5/E+PHjdX+3vr4e9fX1AICamhp4vV7pM6fTqdi2M3aStX/WXPT8fjsQCgJOF66bNRcjDHiu/a1H0fPsOunlHrPsX3HplV9I2+6fPo8R02ZIxwaPHYbrllJpH7xz0f+zFxA8dhhcTi7YmQ648vKin6u48p5y5DzqTAdG3z43+QeSAfQ+0/7Wo+jZ+0fFvhGz5mH04q/FfQ7q78c8V5NkTZeL/7cJV0NBcSMUxIjDTcjT8fvJ29mF32/HdT99Hgi3m1TvW9125W3VKNJ5r9L9TfViujLIz88HAOTl5aG8vBxtbW3Iy8tDT08P3G43enp6FP4EOX6/H36/X9qWpx/wer22SUeQCFvJ6h0HfnW1NArq9Y4D0niukREV6z4HBINSmulLf6iLboeCuHBgL3jvOGXpSWdsuUmB+zPYr2sBQQBcrrgj/tzptwEOJwBxRHh54mRcSSC3VaM/vc9UOLBXXE8hIzh+UsxvpIX6uaY6U8pUWx24clWxffXKVQT1PqNwu2KRdnXvfcD8cbgCJHxOic4pb6tGkspzNeI3TWbwbaoyuHr1KhhjGDlyJK5evYo///nP+MpXvoKysjLs2bMHixcvxp49e1BeXm6mGMOORJ2eUQ5ZRWN1OMKFYyCmjDjdoawxHJ6+K8xUoSDYa78EGBOn+ov+Gdi9U/w+IH4uM2PJ72vE7XOTcoZmg6NRStU9EBJ3OF26zR7ZljKCr1gIYV89IAwAvAN8xUJd31OYhZJ4PmoU74hNo6wy/ZuaqgwuXryIZ555BgAwMDCAuXPnoqSkBD6fD5s2bUJDQwO8Xi9Wr15tphjDCtbeCmHj48DAAJjDAf77T5nWgJSNlQHX5QM93eEPGXDHopgaw4oXj+PCI2Emmq3kigAAOD6qRFSdef/PXkhKqWVDZ8n5poH//lMQmhrAcQA3R39ZTrt2aPGI3Gs6q6GvmzVXnDUlidbAwI5RVpn+TU1VBgUFBdi4cWPM/jFjxmD9+vVmXnrYIjQ1REeWAyEITQ1wGNy4JdPQlcvRnEGMRRVB2OmpVWNY/jIjJzecaiIEgBNNQ9KBHLivfTs6K1B15sFjh4H5+juCbOksOd+0lH6vbAkbVc9aU5Ez8r0RXm9KZiGtgYFZEVZyZzeSVFyZ/k0pHcUQI2KpkW8bQaRRSx14KKgcxcsvOP028F96IHEE0YRJ4NSKIRQUcxA98G3wlfdEr919TmF2ct1SKtqIdWJFZ5lOR5AKdl+PYRdTXaYGBvL77fn9dvCrq8X9SbTBTP6mpAyGGNychWD73ok6aOfos8UOhsI3wMtMO4oL80CkLvEgikBzin7vfeJnEybFvCjKa/Pg5t0Nbs5CMbIiyVFhJl8srY5Ay/dh587bCOT3moypzsxnlKmBgdo/JjQ1APsbLFeG8SBlMMTgfNPAr9mQcry/FsqXGGKeIEGlDG6bBe4zxQmvkaztXnE8A5A/1lYvUDxiHOXh+7TL6DgTqO+Vq1ou5plKMCLPxDPKxMAgxtnNAczGfitSBkMQeUM34sWKKUgPDii+BWj7i9iwnS5wMz4L9nG7OPIJo6WQkJOr2SHEkzN6bQbwDtva+tXEi3oxypGdDbML9b2ir1fXiDwbnP16UDu7L168CBbx6dnQb0XKYIiT7osldeL/NBNoOQCAiefqaBV9BrwD8H8J7PVfiecHwPb+UTQbDQxII0LRH6BcW4CcXCl9wKByRnwTWj4KmyHlFWKQ7lMe9aJlr062Y7fj7EIzyaDGveoZkWeLs18Pcmc3d/68rZ38pAyGOHpeLPmLDO9cxX7FOgKHQ+yspcifcI2Bv7YDAwPREw4MQLQnAQgGwd57O2aEyN00Q5cJgZ04Gr0WEyQl0d96FMKBvbZ6qeQpFgCAOZzgv/+Uwr+htlcDiNuxx1MSdhs5x53VpWibz4RN36qZlZ2d/KQMhjiJXiyt+P1I5Iui0xEAbt7dQP5YZUiowwnus3eAnfjvaEir8griAjSeF3VHZDSs04QQbyTd8+w6UdFYPDKOdZAGox8KA+I+VZoFeYcgvLVDs2MfNIurzUbOgymndMNHzcCOMys7QMpgGKD1Yknmn8C5mPh9NvFi1L4vt9eHF0Gx9lawOQsVC6PYhElitMTpduDDNiiijRgDN+8u0fkrjxTSYULQUmZiBxq0fGSs7SB1RRWCDh9HvI49UQdrtbkh2RW8qYzE9X4naTObzWZWdoGUgY1IxXacSoegDtdUpI3IyVV+BqW9PqYDDIeuRhZLRT8Pr0OIZN1UraZNpkNTKwmx8wl3uimMjI0yEWjNbrj7V4DtfRvIywf/uS8nto/HeQ6JOlgrzQ3JruBNZSSu9zupnNvMmVU2OPbjQcrAJiTT+JWLtKLHA/oWtMSEa95xt5Q2gp3pkJmGmMx5K0RzxctHVfsbYlaUylcYqx3FaoWQqgnB/dPncSEFn4GRJgJ1p8KuXAbefF1KsofPfVn3/WjKoJp9WUVMavIkV/CmMhKP9x21fyuVc5s1s8p28xMpA5ugp1HHXfwVCkL4w++AY4c1lUOkU5Y6bFUnJk8bwV2QLeTieMDBi8nE5CaMyHd5h2Y1LoVNvHG3mIxOEMBkWUjTHUGNmDYjpcyShod2yiOjXv+l+KyAmCR7mt+Nc+/xZl9WoNnBJWsWSmEkHs9XFJl1Mp7H5eWPpjzKN2Nmle3mJ1IGNkHXC6Zw6HLh3BPh0fuRA9FRfDgnkbjaMZI2ggPjeWDSZHBz79YcGQmNu3HptZfCEUMAOIC7f4VCkbD2VmnEyhiA996O2/hZe6tmBwnEj6AxmphcOAaYCGKirCoWgQucU+VW4uN2konu3U6dSrxZQLJmoWRH4oP7ihgwMIBLL9eC//7TYujywX3gPnuHpZ2v3Rz7yULKwEZwFQvF0H2NBG9AbGPDjcXAqWPihxFFELHRS6sdI45cJnbKH54C+/AU8L8eBn/vfVLpQeTkSiN4CUEQo3wix726WTET4KqWg+2PZCDlgcA5sPbWqINYCguN3ACfVFqCdHwi8UxpRpgIlKuLBaBxt7gwziHLxnprmfZ35ZXM4ty7nTqVeLKoR9aJ0k6kkghO01fE89EwZkFQpng4dRxswiTLFIIdHPvpQMrAQhSd1v/5tVR9DHFyu6vt8ez1X2mclIGrWg5uwiRxtWOchHLsvbch/LVd6tzBcdERfIRwNEzU1CNfSxANBxWaGoCmd8Aad4Pt/aOUZE7h6A0nn1NEEoWC4nVzYosbpWp/lafwFjOqhtdCyDredE0EUgcZDEJaa8EE0ffCQXymR5ohHDssmuvCazdYe6v4mer5xpzfok6Ftbfi0z0dYBMnJxW5pC5hyd2/QlfaiWThfNPAPfDtcA0MAXCNsF2KBzuvI0gEKQOLkDq7YCQuXVbQZX9D3FFXpLEJb+2I7bwjx3/cDr7yHuklZh+eBFr+S3nQ6Q6wj07JdnDhqKKwYvBNAzf+Bggt/6VRZyBal5fzTQN34ijYQEiavrPXXpJGaPyjT0ZX5E6YJH7dN02cVbz2S2BAAPv3l2NGdIrRZVB8JpH96sVxchQpvAUWTqCXfB3hwWYlivvaWy/KGPa9SLMheYH28DoDxUyJ40TzUpyOI9OdSqQ99kX8SPJ1DQlkYRFzJCD+/9d205QZX3mPlNAw2RQP2RzpkwlIGZiMerQVUyZSnf0TADvSDOGG3eLofn8DWKNol2ccD9w2C8i9Dtwknxj6ORCrECJpqyXF0bgb7OgH4rEcB0y5WcwrJP+C0wWuajlGswF8yjnEDrrtL7GzCo4H5t0tmbKk9NLgovciCAqFJs5QQhD2N0Q7mb7e8LmZ5ohOmQ+Jgb33RzEbazjFRWRxXIyiVKXw1ptAT/2bJZqVRJ4tm7MwpoNRONi7z6G/9SjgHafpuLcLCuWL5DKKqpsIYwCfhjLTW6kvmRQP2R7pkwksUwYtLS3YunUrBEHAokWLsHjxYqtEMQ31aEuRo8fhUPSfCi50g217MfYjJgAt74t/ggMm3gic+Uh5EqdLEX3C2lvDq4UFsSOfMg0YnQNpbQEvOkCljulMB9jRg1G7thqOUygCyZGqklaKSI3jH9ATR4+KRUDjbkhmmFBY8ckWx8W84KoU3vznliTtbxD+4/Xo/Qf7IfzH63HTcsfYtcOzBuEPvwOONAONu9Hz/p+kFNZ2tSnr9VOoI3q4B76dcgnLuOeXmZz4NRsSPic9syg7OeXtiiXKQBAEbNmyBU888QQ8Hg9+9KMfoaysDBMnTrRCnKTRvTJSEf3TD1b/H7IGqerqeYdoZ+3/h14pgDMfRjcdDnBz74qJSZdkQDjq6NTx6PHh2gDyjr0vXtEa2XWlmG/5/XGcqGzAFB3CYA7IRB0jX7EQwv6GqIM64tcIF7cR/rMuOruSR7qoUnjrJWq6UynC4y0QTh1PbjR59INoVFawX7e/wsqcOfyjT2LUmQ5clvkMYuQ7cVT2zAfAXvsluMeeTqmEZcy521shvPGywuSkNpmmip2c8nbFEmXQ1taGwsJCFBQUAAAqKirQ3NycFcogmelmNPoh3Cl0ng53mJGTKVM26FcEWoIxzVz/0kug7uQHBMXx0Y5d5Ru4bRaQ6waa3olZbxDzgoVj7hUraQfp9BN1jOrvApD8D6GP20WZIjMSWShnqvb2qOLUIInRZEwUFa8dZhrzPYtNGZxvGkbfPhdXBisalJMLxSwwnH8p3bKR8RQx+7hdEaGWKqnMyoabj8ESZRAIBODxeKRtj8eDU6dODfIN+5BsWCSuLwQ6z8g+CHcSMYZWASkziINUMls0NQD7/hj1MTiVx0eVhmyk73RJKRVYRaxtHEgcDhuRIdWXSf5d1t4q+R8u7fujrMPlwM31p/XCCo27wQ6/r6zpzMn+T2ZBkyqKaszyR3FZrxIx0ZRhSOfW16vcNqjGRFxF/HEbhNonDFGMybRDqxWzFViiDJiGGYLTKNZbX1+P+noxFK+mpgZer1f6zOl0KrYzRf+suej5/XYpDPS6WXPFwtwyLu/ehUsv10Yjc/TA88rRZAIcEyfBdUspXJ+ZiuCHJ8EBuDYvL0YWAGLkze1z0f/5f8HVd/8TDMDIOz8vplaWHdP/sxcw8JcjYKNywPp64bqlNHpM+BzSc2g9KmYODT+HvM//i/a1DeTTPR2i/4UJgMCHM6Ey8Xf4XOrXv7x7Fy5te1HaHnF7Ja4pnQ3W1wsuJzf2WQxCf+tRBM90gHvwe9L3Rv1TKUbFm3HIv6ujbaWK/PdiThfcP31e834SvVeSjMEgwHMYs/xRjFJlZU1JPvm9Oxxw3jgFobZwzYyBEEad6cBo1XXM7AMUbS3O9ZPBqv4qGSxRBh6PB93d3dJ2d3c33G53zHF+vx9+v1/aPi+bvnq9XsV2xvCOA7+6Whph9XrHKWrxsvZWCL+ujYZ96inI4nACM8ok57AeBrzjIIwag398+ilYw1tAKIQrDW9pmmrksuMrywAAQSC2hrB3HLz/MkN6rlfCxyjywQDKaKhwWccLB/amlB4iGdjEyeKzQkh0lP/Ph6R7Vf8OyTDQ+LZiu/9sJ4JnOxXPMPIsBpVP7lCPjCa94zAiFNLXVhO0rXQQDuyN/l7BfvRse0nTKZ7wvVLJeNk3DZeNkFF1XgEAap+QTJCXJ06OMV+Z2Qco2lqc6yeDVf3V+PHjdR9riTLw+Xzo7OxEV1cX8vPz0dTUhEceecQKUZJC3ilGirjHHKO2F8fDMxboCYhKg+OAwgnxj42MgCOKxeEEjn4AduSAMkdRsB/s1c2ir1iWB0iO0LhbWrrPV96j656l6TIfduIODIjRSKoaBWYjt/vKq4elfd7P3gF2vCW64/SHoq1aZh5QLBD8uB3o7QGX51Y47I0w8xi9viCm3GiwX2xHqTjFDZRRszKa6rxWRl7ZOfLLLCxRBg6HA8uWLcOGDRsgCAIWLFiAoqIiK0TRjW4bYk5u1NY8GN3non8PhIC3d2kf5xqhSIaGv7aDnfzvqB9C4MKdsqBUGKoi7Gx/A9jfTkvpK9jxFjGyNIFCUEZEyc4vK3Zj1ssyWIcxwus1bNTMV94DAQA7uA8YcW04z5Ns0RgghVPKf1cGgO17Rwp/tFvEinqmgkX/rFxAOEgyvUzJNdi7lImFd4kWFw4HJRDBsnUGM2fOxMyZM626fNLozSrK/v3l5J3BHKc9m/BcD+7e+6QOm7W3QoikrYjgcIo59N97G1CtKJYyPcpKMSrkPbgPSKAMFB2cKrxTT2plI2ouGOHASyQHX3kPUHmPeN1jh5TZMrUirSKEVG3BJmmngdg2i7+2Kw+Ik0zPcDl0pMC24lllwkkcuff+WXOlCoJ2hVYg6ySprKKDce0o4Opl1ck5wBEe3YMDRo0GLn8KBM4rUjWI51d16hWLwFfeg4GP25XK4LZycD5ZpkctbvAlvm+t8E6dnXs6L5uRHUZS4cCq+2WffAzW8n44EEBjxudUp1e2Pu10hJg2+9k7wE4dj8kVFa8GthGkmgI7E2Qieity7z2/3y4tPLQrpAw0iDeK1JVVNCYfggyHE7jz88Affqe6oLIsJDtxFGzXa5Dn5YmaIZRlFSOLu2JWgYYLq8R8Ry7vyFG6nofWKttESCt5pZTDyb1sRnYYyb70kfsVGneDyaKMUHI7uBllks8AALg8d0rXyARadu9IXh95SvJ4NbCNQCtLa6IU2JnCbKWkzG5rjUkuGUgZqNAayQBQ2l4rFmoqDKHlv+ImjwPHgXvgW6J9emyhaNY53SF2lBomF0Venr31Ys1h3zTA/6Wo3ZePLmDjfNM0V4Fyvmniitz9DWAXe4CjBxFJrGbWiEyZhC9xnH4830DaqabVztMkX3p2cJ9yR/8/lCa72ifAQiGwpoZw/ePMjXb1mt+0FLminamUWPDYYWC+McpgsCytdrDHm+0kVigbp8tyH1IiSBmo0BrhAVDsU+RQlykM7N4Z/8QcJy3Ykdun463MVeflkeR4e1fUVJFEmT/Fwi2TR2SKFBgcB0y/LW5un8FMOOl0GGrnqTzkFgCEt3Yk7khVUUbcZ++IvUdZ/eNMjXbNLN3puqVUDKNNQzbJxJZEllarMFMpmRX9ZhakDFTEmzoq9mnkUBcPijmbOFpUpXGQPh2kISry8sgdmXHSHNgpQiMmO2ccRQCYZ17R7KwjRXp0dqTyKCN1KK5WO8nUaDfZVfCDKSf16HjEtBkpR2nF+E1UsyU7ZWnNFGZEv5kFKQMNtHwDiqIyH7dDyvopVxgumW1+ynTwX/4mACRM/iVH8QJrlPNjLpdog+U44K7F0f1pdqpGRj0kM/02y24bV6kn+Zwis7iY81sYh64rmCFJp7kRkVqKtOyhYEZnS0T6kDKQoaiS5XAoKo5JfgGpID2vyPoJQFHIJTIKYieOwjVrLjgdHWy8BV7ycn5SURhBAN75f2Alt8ucy6l1qmZEPcTrYLQK9ZjRYcQ7r5HKxyq7t55nZqZDW70iPfpOyHKyMwbk5NrCN0Dog5SBDEWVrLBvwCF3tskjIwTtLKFSIZd99VJnrreDjbvAS/4yxykKk06nmqmoh3ijVaM6DC1Fo3leG60FSJVEz8ysGZf6N8SchbI2K8vDJfORmYFWiU4iPUgZyIiJCu08jYFN60UzzYRJYO/9MfpZeOQjbUpFUcLTZEF7NfCg1x9sgVectNHIyVU4Q1N5MZKNekh5IZnJo9VEZhE7rgUwCzNmXELjbrDdOxWhohyHhG3WaAYr0UmkzrBSBgnL6cmrZPG8InUDK745Nmw0PPKRV3+SImh4R/TF0BlWpn6BgdgFXvJjkJMrVU5LJ6okmaiHdCJZzIzr1rVCPE6k2FC1aRtpoolZcwGIZT0ZYiK1MhatpqNEJ6GfYaMM9Na1jVTJYoffV67oPfUX5Qk1nZLKUMrIZ8mElWnFhcc7RlxdbMxIW2/UQzqjezOdrnoUjeasKkty1utJMGhm2HDMmgu3B+i9CLz3Ntj+2JBgM7HLCuahxvBRBjo7MamjzckFU+T6UcWN3hEtpjJYKKWZYWVWvBTpXtMsh6IeRaM+JqZN7G+w5SxBPiqXJxhULKr7a7u4wEsQTFFsMZldJ00RazxbsOI68juOONyEq1euZuSaw4EhrwxSXYXKV96DgXN/jy4k4x3hGH+xmAqvijSyIoTOiutaGVKZKIeOHkWjPkZSbLzD1M40HdSjcnZwH9iESZqZVAEYkstJ/fuq11xwEyZBOHbY0tH51T/9JxAMQtjfYKvfK1sZ0spgsFWoehqO48v/W0wdsfdt4OMOAGIKCK5quaaJyaoww0xfN6ZDzcSqZhNy6MgVG+s+B7z3tiUj3YS+LI2V0ArTpOLg5Ep0qmW4XDAOwpbnNE1n6jUXVq4hkJI22igXVLYztJVBnFWour8fSUktL9LNmKkhc9mG2WmAFQuaBsmhk6pCiig21t4as+I7E+h5floroVl7a7hmdTDqq3I4xSy2g9Sjll9Xa63AJXmxpASdrJVrCBQJGMlvYAhDWhkvJCENAAAYU0lEQVSka9+O5tiJnDC1UddQxsxwUaFxN9jr4QV2Dqeispo8h44RCskq85fe56celaujypKZ8Q6+ViBcWxrM1m2d802D+6fP48KBvbbz8ajJljURpimD7du345133kFurhiLf//990vFbOrq6tDQ0ACe57F06VKUlJSYIkO6L7hCmfAO3aOuoUzMwi4TFzeJK63D4bwDIXCV90hpvuU5dIxSSJaY3NJ4fqnKq35eirUCqtrSVrd1rfrbEblGTJthet3tdMmmNRGmzgy+8IUv4Etf+pJi35kzZ9DU1IRnn30WPT09qK6uxi9+8QvwsnTMRpLOC26ls9SOxBuBm/GMNJPyxVkxnM2hhpYEAaif15yF4hobm2XXjJueJZIp2OBCPGaQTWsiMm4mam5uRkVFBVwuF66//noUFhaira0NU6dOzbQourDSLmo34o3AzXhG3E0zoon/uGhVLs1js0hp6ykEbzbxnpeZYdCpkCg9C263vzLIpoGKqcpg9+7daGxsxOTJk/GNb3wDOTk5CAQCKC4ulo7Jz89HIBAwUwzCIDLZsJPt4LNBaWei5q5esuF56UnPYnci7TiZzMVWwTGmVeVbH9XV1bhw4ULM/qqqKhQXF0v+gjfeeAM9PT1YuXIlfvOb32Dq1KmorKwEALz00ksoLS3F7NmzY85TX1+P+nqxUlJNTQ36+6NRPU6nE6FE9YZtwlCStb/1KILHDsN1S6lot7cIuZx2kSkeEVk//d1v0ff6r0XzF88j54EVGP3lb1gtnuL5jfqnUlu1VblsABS/c6rvlZ72YnSbsqoPGDFihO5j05oZrFu3TtdxixYtws9//nMAgMfjQXd3t/RZIBBAfn6+5vf8fj/8fr+0fV42ffV6vYptOzOkZPWOA+aPEyN5DLynZENDI3Kq15LYyUEXuaeIHZ5NnCxGRUGcWV2eOBlXLG4X6ueHn71gG58BAGV7AxRtL5X3Sk97MaNNWdUHjB8/Xvex5nhtAfT09Eh/HzhwAEVFRQCAsrIyNDU1IRgMoqurC52dnZgyZYpZYhA2hbW3Qnhrh9hhRuoJ73pN/L+9Vf954iSfsxr5PfX8+BGw9lbJZMD9j68Z0sHIn2HK51A9v+Cxw2nJlA5G3E/Ca+hoL3ZtU2Zjms/g1VdfxUcffQSO4zB27FisWLECAFBUVIQ5c+Zg9erV4HkeDz74oGmRRIQ9iUklXbEw5dBQuzro4tWIMMpWb5T/wegayKmSKX9KKgkN7dKmzMY0ZfDd73437mdLlizBkiVLzLo0YXNiRl4MSeWNkmPHSCLW3gr24SmAA8Q62PpSmCd1DQPXVhhVAzkdzFy8KCeVhIZ2aFOZYEivQCbsSUyW14qFQMVC3S+fekWnnSJjWHsrhGfWRmthczzGLPtXXDYj9DbLy3cqZMhwpNpQiE4zGlIGRMYZLM49EXZf0SklUIvuATMhl9VQG70OtfvJRkgZEJZgSCoFG67oVCRQAwDeYZodPptGr3qixbLpfoYipAyIrMLuzj3OF66Wt78BjAF8xULL7PB2wU6L7Yj4kDIgsopsWNFJI1wlmXIOE+lByoDIOjjfNIy+fa7lC7YIfdh9NkeIkDIg0iYTlc6yjUQlOjMtg5W/SzznsF3kI0RIGRBpQfbgWNIt0WlEJ2m330VtOrObfISJ6SiI4cFwXbqvRSSdgtDUkHKKh0SpOfSmbLD772J3+YYjNDMg0oLswSKK5GYOR9wSnQnPM4izNZnRtN1/F7vLNxwhZUCkBS0WElEWYgG4eXdrluhMxGCdZDJROXb/Xewu33CElAGRNtkaSjlYfd1k0SwlaUCuIPk5kh1N2/130ZKPnMrWQcqAGJYkqq+bbEdk5Eh30E6yarltitUbDTmVrYWUATEsSVRfN9VRvRmdl50L+BgJLU6zFoomIoYl3E0zxLTZPA/wDsDpEv+2oTNzuETeKH4TG/4OQx2aGRDDErVZB0jPZ2AmwyXyhpzK1kLKgBi2qM06du18hlMnaXen91AmLWWwf/9+7NixA5988gmeeuop+Hw+6bO6ujo0NDSA53ksXboUJSUlAICWlhZs3boVgiBg0aJFWLx4cXp3QBDDAOokCbNJy2dQVFSENWvWYPr06Yr9Z86cQVNTE5599lmsXbsWW7ZsgSAIEAQBW7ZsweOPP45NmzZh3759OHPmTFo3QBAEQaRPWjODiRMnau5vbm5GRUUFXC4Xrr/+ehQWFqKtrQ0AUFhYiIKCAgBARUUFmpub456HIAiCyAym+AwCgQCKi4ul7fz8fAQCAQCAx+OR9ns8Hpw6dSrueerr61FfXw8AqKmpgdfrlT5zOp2KbTtDshpPtsgJkKxmQbIaS0JlUF1djQsXLsTsr6qqQnl5ueZ3WCRmW8d+juPiXtvv98Pv90vb52VL+r1er2LbzpCsxpMtcgL2kFXvyl6rZE1l5bEdnqterJJ1/Pjxuo9NqAzWrVuXtAAejwfd3d3SdiAQQH5+PgAo9nd3d8Ptdid9foIg9GP3lb3JyGeHOhFDFVMWnZWVlaGpqQnBYBBdXV3o7OzElClT4PP50NnZia6uLoRCITQ1NaGsrMwMEQiCCGP3RWt65VOn9+5vtdd9ZDtp+QwOHDiAV155Bb29vaipqcGNN96ItWvXoqioCHPmzMHq1avB8zwefPBB8Lyod5YtW4YNGzZAEAQsWLAARUVFhtwIQRDa2H3Rml751EojeOwwMF9/0SBicDgWz8BvQ/72t79Jf5O90ByyRdZskRPIjKyJbO7Z5jPQklcyJ4WVhvtnL6A3iQpyVjIkfAYEQdgbPTZ3uy9ak8sX737UK7GTqRNBJIaUAUFkOUMt2+dg92N3pZbNUNZSgshyhlq2z6F2P9kCzQwIIssZaonshtr9ZAukDAhiCDDUzCdD7X6yATITEQRBEKQMCMJqWHsrhLd2gLW3Wi0KMYwhMxFBWIjdU0UQwweaGRCEhdg9VQQxfCBlQBAWQmGUhF0gMxFBWAiFURJ2gZQBQVgMhVESdoDMRARBEAQpA4IgCIKUAUEQBAFSBgRBEATSdCDv378fO3bswCeffIKnnnoKPp8PANDV1YVVq1ZJhRWKi4uxYsUKAEBHRwdefPFF9Pf3o7S0FEuXLgXHcWneBkEQBJEOaSmDoqIirFmzBr/+9a9jPissLMTGjRtj9r/88sv41re+heLiYjz99NNoaWlBaWlpOmIQBEEQaZKWmWjixIlJlVXr6enBlStXMHXqVHAch8rKSjQ3N6cjAkEQBGEApq0z6OrqwmOPPYaRI0eiqqoK06dPRyAQgMfjkY7xeDwIBAJmiUAQBEHoJKEyqK6uxoULF2L2V1VVoby8XPM7brcbmzdvxpgxY9DR0YGNGzeitrYWjLGkhKuvr0d9fT0AoKamBl6vNyq406nYtjMkq/Fki5wAyWoWJKuxJFQG69atS/qkLpcLLpcLADB58mQUFBSgs7MTHo8H3d3d0nHd3d3Iz8+Pex6/3w+/3y9tn5cVv/Z6vYptO0OyGk+2yAmQrGZBsiYmGTO+KaGlvb29EAQBAHD27Fl0dnaioKAAbrcbI0eOxMmTJ8EYQ2NjI8rKyswQgSAIgkiCtHwGBw4cwCuvvILe3l7U1NTgxhtvxNq1a3H8+HFs374dDocDPM9j+fLlyMnJAQA89NBD2Lx5M/r7+1FSUkKRRARBEDYgLWUwa9YszJo1K2b/7NmzMXv2bM3v+Hw+1NbWpnNZgiAIwmBoBTJBEARByoAgCIIgZUAQBEGAlAFBEAQBUgYEQRAESBkQBEEQIGVAEARBgJQBQRAEAVIGBEEQBEgZEARBECBlQBAEQYCUAUEQBAFSBgRBEARIGRAEQRAgZUAQBEGAlAFBEASBNIvbbNu2DQcPHoTT6URBQQFWrlyJ0aNHAwDq6urQ0NAAnuexdOlSlJSUAABaWlqwdetWCIKARYsWYfHixenfBUEQBJEWac0Mbr31VtTW1uKZZ57BuHHjUFdXBwA4c+YMmpqa8Oyzz2Lt2rXYsmULBEGAIAjYsmULHn/8cWzatAn79u3DmTNnDLkRgiAIInXSUga33XYbHA4HAGDq1KkIBAIAgObmZlRUVMDlcuH6669HYWEh2tra0NbWhsLCQhQUFMDpdKKiogLNzc3p3wVBEASRFob5DBoaGiRTUCAQgMfjkT7Lz89HIBCI2e/xeCQFQhAEQVhHQp9BdXU1Lly4ELO/qqoK5eXlAICdO3fC4XBg3rx5AADGmOa5tPZzHBf32vX19aivrwcA1NTUwOv1RgV3OhXbdoZkNZ5skRMgWc2CZDWWhMpg3bp1g37+7rvv4uDBg1i/fr3UsXs8HnR3d0vHBAIB5OfnA4Bif3d3N9xud9xz+/1++P1+afv8+fPS316vV7FtZ0hW48kWOQGS1SxI1sSMHz9e97FpmYlaWlrw5ptv4gc/+AGuueYaaX9ZWRmampoQDAbR1dWFzs5OTJkyBT6fD52dnejq6kIoFEJTUxPKysrSEYEgCIIwgLRCS7ds2YJQKITq6moAQHFxMVasWIGioiLMmTMHq1evBs/zePDBB8Hzot5ZtmwZNmzYAEEQsGDBAhQVFaV/FwRBEERapKUMXnjhhbifLVmyBEuWLInZP3PmTMycOTOdyxIEQRAGQyuQCYIgCFIGBEEQBCkDgiAIAqQMCIIgCJAyIAiCIEDKgCAIggApA4IgCAKkDAiCIAiQMiAIgiBAyoAgCIIAKQOCIAgCpAwIgiAIkDIgCIIgQMqAIAiCACkDgiAIAqQMCIIgCKRZ3Gbbtm04ePAgnE4nCgoKsHLlSowePRpdXV1YtWqVVH8zUgENADo6OvDiiy+iv78fpaWlWLp0qVQ7mSAIgrCGtJTBrbfeigceeAAOhwOvvvoq6urq8PWvfx0AUFhYiI0bN8Z85+WXX8a3vvUtFBcX4+mnn0ZLSwtKS0vTEYMgCIJIk7TMRLfddhscDgcAYOrUqQgEAoMe39PTgytXrmDq1KngOA6VlZVobm5ORwSCIAjCANKaGchpaGhARUWFtN3V1YXHHnsMI0eORFVVFaZPn45AIACPxyMd4/F4EioQgiAIwnwSKoPq6mpcuHAhZn9VVRXKy8sBADt37oTD4cC8efMAAG63G5s3b8aYMWPQ0dGBjRs3ora2FoyxpISrr69HfX09AKCmpgZerzcquNOp2LYzJKvxZIucAMlqFiSrsSRUBuvWrRv083fffRcHDx7E+vXrJUewy+WCy+UCAEyePBkFBQXo7OyEx+NBd3e39N3u7m7k5+fHPbff74ff75e2z58/L/3t9XoV23aGZDWebJETIFnNgmRNTCSIRw9p+QxaWlrw5ptv4gc/+AGuueYaaX9vby8EQQAAnD17Fp2dnSgoKIDb7cbIkSNx8uRJMMbQ2NiIsrKydEQgCIIgDCAtn8GWLVsQCoVQXV0NIBpCevz4cWzfvh0OhwM8z2P58uXIyckBADz00EPYvHkz+vv7UVJSQpFEBEEQNiAtZfDCCy9o7p89ezZmz56t+ZnP50NtbW06lyUIgiAMhlYgEwRBEKQMCIIgCFIGBEEQBEgZEMME1t4K4a0dYO2tVotCELbEsBXIBGFXWHsrhNongFAIzOkE/+iT4HzTrBaLIGwFzQyIIQ87cRQIhQAmAAMhcZsgCAWkDIghD3fTDMDpBHgecDjFbYIgFJCZiBjycL5p4B99EuzEUXA3zSATEUFoQMqAGBZwvmmkBAhiEMhMRBAEQZAyIAiCIEgZEARBECBlQBAEQYCUAUEQBAFSBgRBEAQAjiVbmJggCIIYcmTtzOCHP/yh1SLohmQ1nmyREyBZzYJkNZasVQYEQRCEcZAyIAiCIOD4yU9+8hOrhUiVyZMnWy2CbkhW48kWOQGS1SxIVuMgBzJBEARBZiKCIAjCpllL9+/fjx07duCTTz7BU089BZ/PBwDo6urCqlWrMH78eABAcXExVqxYAQDo6OjAiy++iP7+fpSWlmLp0qXgOA59fX3YtGkTzp07h7Fjx2LVqlXIyckxXVYAqKurQ0NDA3iex9KlS1FSUgIAaGlpwdatWyEIAhYtWoTFixdL9/fcc8+hr68Pn/nMZ/Dd734XTqc5P9H27dvxzjvvIDc3FwBw//33Y+bMmSnJnWnsIoechx9+GNdeey14nofD4UBNTU3ctscYw9atW3H48GFcc801WLlypakmhM2bN+PQoUPIy8tDbW0tAKQk27vvvoudO3cCAJYsWYI777wzI7Lasa2eP38eL774Ii5cuACO4+D3+3Hvvffa9rnqgtmQ06dPs08++YT9+Mc/Zm1tbdL+s2fPstWrV2t+54c//CE7ceIEEwSBbdiwgR06dIgxxti2bdtYXV0dY4yxuro6tm3btozIevr0abZmzRrW39/Pzp49y77zne+wgYEBNjAwwL7zne+wv//97ywYDLI1a9aw06dPM8YYq62tZXv37mWMMfarX/2K7d6921BZ5bzxxhvszTff1LyfZOXOJHaRQ83KlSvZxYsXFfvitb2DBw+yDRs2MEEQ2IkTJ9iPfvQjU2U7duwYa29vV7w7ycp26dIl9vDDD7NLly4p/s6ErHZsq4FAgLW3tzPGGLt8+TJ75JFH2OnTp237XPVgSzPRxIkTpdG/Hnp6enDlyhVMnToVHMehsrISzc3NAIDm5mbMnz8fADB//nxpv9myNjc3o6KiAi6XC9dffz0KCwvR1taGtrY2FBYWoqCgAE6nExUVFWhubgZjDMeOHcPs2bMBAHfeeafhsuohWbkzjV3k0EO8tvfBBx+gsrISHMdh6tSp+PTTT9HT02OaHDfffHPMbDhZ2VpaWnDrrbciJycHOTk5uPXWW9HS0pIRWeNhZVt1u93SyH7kyJGYMGECAoGAbZ+rHmxpJhqMrq4uPPbYYxg5ciSqqqowffp0BAIBeDwe6RiPx4NAIAAAuHjxItxuNwDxB+zt7c2InIFAAMXFxdJ2fn6+JJNa1lOnTuHSpUsYNWoUHA5HzPFmsXv3bjQ2NmLy5Mn4xje+gZycnKTlzjRav7UVcmixYcMGAMBdd90Fv98ft+0FAgF4vV7pe5H2Gjk2EyQrm/q5Z6J9yrFzW+3q6sKHH36IKVOmZN1zlWOZMqiursaFCxdi9ldVVaG8vFzzO263G5s3b8aYMWPQ0dGBjRs3ora2FszkgKhUZI0nk9Z+juPSEzAOg8l999134ytf+QoA4I033sBvf/tbrFy50hZyD4Zd5FBTXV2N/Px8XLx4EU8++eSgM1u73gOQnGyZktnObfXq1auora3FN7/5TYwaNSrucXZ8rmosUwbr1q1L+jsulwsulwuAGLNbUFCAzs5OeDwedHd3S8d1d3cjPz8fAJCXl4eenh643W709PRITiizZVXLFAgEJJnUsrrdbowZMwaXL1/GwMAAHA6H4vhU0Sv3okWL8POf/zwluTON1m9thRxq5O2tvLwcbW1tcduex+PB+fPnpe9acQ/Jypafn4/jx49L+wOBAG6++eaMyHrddddJf9uprYZCIdTW1mLevHm4/fbbAWTXc1VjS59BPHp7eyEIAgDg7Nmz6OzsREFBAdxuN0aOHImTJ0+CMYbGxkaUlZUBAMrKyrBnzx4AwJ49e+KO5I2mrKwMTU1NCAaD6OrqQmdnJ6ZMmQKfz4fOzk50dXUhFAqhqakJZWVl4DgOt9xyC95//30AYoRB5B7MQG6jPnDgAIqKilKSO9PYRQ45V69exZUrV6S///znP+OGG26I2/bKysrQ2NgIxhhOnjyJUaNGZVwZJCtbSUkJjhw5gr6+PvT19eHIkSNS5I7Z2LGtMsbwy1/+EhMmTMAXv/hFaX82PVc1tlx0duDAAbzyyivo7e3F6NGjceONN2Lt2rV4//33sX37djgcDvA8j/vuu0/6kdvb27F582b09/ejpKQEy5YtA8dxuHTpEjZt2oTz58/D6/Vi9erVhoaWxpMVAHbu3Ik//elP4Hke3/zmN1FaWgoAOHToEP7t3/4NgiBgwYIFWLJkCQBRwalDSyMzIaN54YUX8NFHH4HjOIwdOxYrVqyQOqRk5c40dpEjwtmzZ/HMM88AAAYGBjB37lwsWbIkbttjjGHLli04cuQIRowYgZUrVypCko3mueeew/Hjx3Hp0iXk5eXhq1/9KsrLy5OWraGhAXV1dQDEEMgFCxZkRNZjx47Zrq22trZi/fr1uOGGGySzzv3334/i4mJbPlc92FIZEARBEJklq8xEBEEQhDmQMiAIgiBIGRAEQRCkDAiCIAiQMiAIgiBAyoAgCIIAKQOCIAgCpAwIgiAIAP8fUuK/3Z1VUGIAAAAASUVORK5CYII=\n", "text/plain": ["<Figure size 432x288 with 1 Axes>"]}, "metadata": {}, "output_type": "display_data"}], "source": ["plt.plot(fact[:,0], fact[:,1], \".\");"]}, {"cell_type": "markdown", "metadata": {}, "source": ["## PCA with statsmodels and normalization"]}, {"cell_type": "code", "execution_count": 26, "metadata": {}, "outputs": [], "source": ["from statsmodels.sandbox.tools import pca\n", "from sklearn.preprocessing import normalize\n", "X_norm = normalize(X)\n", "xred, fact, eva, eve = pca(X_norm, keepdim=2, normalize=True)"]}, {"cell_type": "code", "execution_count": 27, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([3.65433661e-04, 3.38164814e-05])"]}, "execution_count": 28, "metadata": {}, "output_type": "execute_result"}], "source": ["eva"]}, {"cell_type": "code", "execution_count": 28, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[ -0.21720145,   2.87429329],\n", "       [-48.46551687,  13.57394009],\n", "       [ -5.89658384, 147.68504393],\n", "       [  2.99888854, -11.96508998],\n", "       [  4.79280102,   5.15588534],\n", "       [ 17.88981698,  85.8816515 ],\n", "       [  0.86125514,   4.75280519]])"]}, "execution_count": 29, "metadata": {}, "output_type": "execute_result"}], "source": ["eve"]}, {"cell_type": "code", "execution_count": 29, "metadata": {}, "outputs": [{"data": {"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJztnXt0VeXd57/PPiegwBsuSYtcBMspkJYBFIutMdoaYm2dvo5tp76A9V0u34GpdJWZVbGvVWvr4CVdwtsuu+rrjFpqtZbRVrxMWV5CliIEKwTRFEyK0aICvpggIhdNztnP/LHvt5Nz2cnZe+f7Wcsl57b385yd832e/bsKKaUEIYSQxKBUegCEEELChcJOCCEJg8JOCCEJg8JOCCEJg8JOCCEJg8JOCCEJg8JOCCEJg8JOCCEJg8JOCCEJg8JOCCEJI12pEx84cKBSp/altrYWPT09lR5G2SRlHkBy5sJ5RIs4z2Py5MkFvY87dkIISRgUdkIISRgUdkIISRgUdkIISRgUdkIISRgUdkIISRixEnbZ3Ql146OQ3Z2VHgohhESWisWxF4vs7oS69iYgm4VMp6FceytEpq7SwyKEkMgRmx277OoAsllAqkAuqz0mhBDiITbCLmbPBdJpQFGAVFp7TAghxENsTDEiUwfl2lshuzogZs+lGYYQQgIIRdiPHz+Oe+65B++88w6EELjmmmswa9asMA7tQGTqKOiEEDIAoQj7unXrcOaZZ+Laa69FNpvFJ598EsZhCSGElEDZNvYTJ07g9ddfR2NjIwAgnU5j9OjRZQ+MEEJIaZS9Yz906BCqq6tx9913Y9++fZgxYwauuuoqnHLKKWGMjxBCSJEIKaUs5wDd3d248cYbsXr1asycORPr1q3DqaeeisWLFzve19LSgpaWFgBAc3Mz+vr6yjlt6KTTaWSz2UoPo2ySMg8gOXPhPKJFnOcxYsSIgt5X9o69pqYGNTU1mDlzJgDgS1/6Eh5//HHP+5qamtDU1GQ+jlqh+zgX37eTlHkAyZkL5xEt4jyPIWu0MW7cONTU1JgdkTo6OjB16tRyD0sIIaREQomKufrqq3HXXXchm83i05/+NFasWBHGYQkhhJRAKMJ+xhlnoLm5OYxDEUIIKZPYlBQghBBSGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGBR2QghJGKEJu6qq+NGPfoTm5uawDkkIIaQEQhP2jRs3YsqUKWEdjhBCSImEIuy9vb3YuXMnFi1aFMbhCCGElEE6jIP89re/xXe/+12cPHky8D0tLS1oaWkBADQ3N6O2tjaMU4dGOp2O3JhKISnzAJIzF84jWiRlHvkoW9jb29sxduxYzJgxA7t37w58X1NTE5qamszHPT095Z46VGprayM3plJIyjyA5MyF84gWcZ7H5MmTC3pf2cLe1dWFHTt24JVXXkFfXx9OnjyJu+66CytXriz30IQQQkqgbGFfunQpli5dCgDYvXs3nnrqKYo6IYRUEMaxE0JIwgjFeWowZ84czJkzJ8xDEkIIKRLu2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2AkhJGFQ2Akhg47s7oS68VHI7s5KD2VYkK70AAghyUZ2d0JdexOQzUKm01CuvRUiU1fpYSUa7tgJIYOK7OoAsllAqkAuqz0mgwqFnRAyqIjZc4F0GlAUIJXWHpNBhaYYQsigIjJ1UK69FbKrA2L23EiZYWR3ZyTHVS4UdhIrkvpDTDoiUxe565Vk2z+FncSGJP8QydDjZ/tPyt8TbewkNtAJR8IkybZ/7thJbBCz50Km00Aum7gfIhl6omz7LxcKO4kNSf4hksoQRdt/GJQt7D09Pfj1r3+NI0eOQAiBpqYmXHLJJWGMjRAPSf0hEhImZQt7KpXClVdeiRkzZuDkyZO4/vrrMW/ePEydOjWM8RFCiANGRg1M2cI+fvx4jB8/HgBw6qmnYsqUKTh8+DCFnRASOoyMKoxQo2IOHTqEt956C5/97GfDPCwhhABgZFShhOY8/fjjj7F27VpcddVVGDVqlOf1lpYWtLS0AACam5tRW1sb1qlDIZ1OR25MpRCVefR1dqB/9yuomnMWRtSVFr3iN5cwjjvUROWalEsU5tF3TgM++PMjQLYfSFdh3DkNGFHkmKIwj8FGSClluQfJZrP4+c9/jvnz5+Mb3/hGQZ85cOBAuacNldraWvT09FR6GGVTqXnY7Z4AzNtluG6Xi7GPuudivw13HzfK8G8rXMq1sUdlHqUwefLkgt5X9o5dSol77rkHU6ZMKVjUSTIwfmAYUw25/l7T7inqG30z+sq1jyY5U5AUDiOjBqZsYe/q6sLmzZsxbdo0XHfddQCAJUuWYMGCBWUPjkQXx+5ZEUBOBSA1wZXQMvpyWUBJQfa+by0CZQgzE5QIKYyyhb2urg6PPPJIGGMhMcIh0qrQ0rIhgVQaSn0jUN8Iua0VcksL8OKzULe1QixeZgl+CcLMBCVCCoOZp6QkPLvnxcuAY0cdgiu7OgBVNXfoOHa0bGHmbXh5RC0GvJDxBL0nanOJEhR2UhKF7J79TCcU5soRtRjwQsYT9J6ozSVqUNhJyQwk0jSdRIuoOZ8LGU/Qe6I2l6hBYSeDCnfoQ8dApomoOZ8LGY+YPRdSUQBVAkIx3xO1uUQNCjshCaAQ00RU7qDsC5DfeNw5ERACkFL7v05U5hJVKOyEJIBCTRNh3kGV4rz0W4CUS74T+DrObQRyOe3FXBZqWytS+rl4NxgMhZ0kiuEaKRG2aWKg77FU52WgzdzIczj8vuN1IaCZYnKqtmtv2wRZ3xg4puF47f2gsJPEMJwjJcI0TRQUrVKi89JvAXImuynafzkJKCmIcxu1873wtHYANefvZDWO0d8PKQQw/xwoX/vWsLn+bijsJLa4d2jDPVKiHNOE/bss5Hss9Q7BbwFSNz5qS3aTAHRbul7GSpzbCNnWmvdcsqsD6O8HILXP7XoJ6l/boay6bcj+BqJ0x0BhJ7HEb1eJMdVaeQNVMFKiCNzfZSEZwuXcIbgXIMciIWzlKaQK2dUB5ZLvFJYzoSiAmrOezBa2uIchyFG7W6Swk1ji3lWqba3AtlZNFBQFYvGyiu+a4oL7uyw0Qzgs56V9kTALyrkWlUJyJsQV34N86G5zp4+09XlHpE1tgzX3kAQ5aneLFHYSSzymAAHIbBaA1P47drTSQ4wNUcgQtp9PTpnuLAG98dHCd9NmaKQCsWS5b5Zq3//6FVA7STtXSIIctbh6CjuJL3POAj48DNHwVYgp0we0wxJ/ohYTboh8Mbtp2d0J+ft7tNpE2jPm4u4W7/7drwBf1oTdLcgYU13cQmIbc5S+Qwo7qTjF2jhldyfUO2/QfowA5Dt/h1h1m/N2Xm+ZVukfWFyIYkx4Mbtps+CcgRKcpVo15yyc1N/mawYq0SwTpe+Qwk4qSik2TtnVYSWtAKaTzEh0KeR49iYh7qqUYSK7O3H8hTchp86Ifex1X2cH1Je3lFSJsRScTlUFOKzV9QfgOYeYPReyqkprmScUiKXfs15z7aZH1M0FbB2UDEF2ROf4xNjH4RoZUNhJRSnFxilmz4VMpcwdu8NJVkhhKVvMM6Clqst0VeiRDMZ5juk7RffxoxZJYR+XX5r/B//2Ey1OvMhKjKViCLLa1qolJm1+FnLrJs2Grqra30D9Iih6wpJYvAyyfSvE2edBueBia0xtrRBCC5v0K1vgWBzyxNhH6RoNBIWdVJRSnE4iUwflutsdP1hAc7JhTPWAoXqm+EOPnpByUCIZHIsMvMePWiQFkKdMbleHthsuoRJjOYhMHURXB2Qupx03a7RolkBWBTY/YzZxMc0oe/dATpkOAFDv/DGQy2ku9a2bIJYsx4eH9kNt+X++i0PeGPsQ5jRUu38KO6kopTqdRKYOCjQxkfv3QT78vzXzTCoFsfS/5zWvmItJtt8qLjUIDteBFq2oRVIAweKsjVU3deSrxBiCI9KDKz9BW4hzMCOgclnI9q2ecePw+y6TXT/kQ3fjYyMcEnAsDsYiFhhjX+Y1GsrdP4WdVJxSnE7Wj6TfilsGtB/1vm6krlyR93x2h1kYNnbfW3v9PKPefRMnfGzsUYukAIKFTGTqMP6Wu3Akj429EEdkkGkkCNndqce16/kJS5ZrEVBG20WpauM8+zzIvXsc41bbWn0OKL3PIfiOLegaldL5aSjv0CjsJJZYpgHvD9VW3TWQMCMY8u3ERKYOo7/YgJM2Z91gjSMM8i02I+rmQtHjv/N9PsiEAcCMZjJMIwOl/DvNZloIoxkOeW6jUzht8e/GHZ269Tnnrt2NUetdSQXuxt3XqOTOT0N4h0ZhJ/FkTLX/84pi2tzDoKCdWQRt5eUQxmLj64gMiGbKm1GaRww9ZhN7/LtuAsJ5FwGbnw4eqGmy99vJB3ykxM5PhZRGCAsKO4kd5u25+7eoKBBXXDMkO3E7UbSVDyX5zFAeE0ZANFPQsdzHAWxOcpsJzR6+Kv/wfzR7eroKYslyyzfgO3jV/H85FSoLfc9Q3aFR2EnscNyeCwFM/yzE9ExBNtuSzlNA84qwd2JxiZ0eyAzl3lG7o5kCTRzuUEZjJ25LTDPCVO0RMWZ1R0AT9xefBZouBbpe08Jb9+9z7s6Fon0GIvgu0EUh17vS/hMKO4kd7t2QMkgFv4rZiQ+VzT5qFGuGEpk6swNS3mP5RKuoRskI8wM+ETFu/r5X+y9wQEaJYBVy/b1QgYKc6YVc70r6TyjsJHYUshsqdcfr/lwldl1xstkXs/gFJT6ZzxnHMhLH3NEqR494D2pExLz+WmkTsNeWyfZDPnwPoMrIL6gDQWEnoTGU5oN8u6GSyhToYXho26QltNg+N9Q/7jjZ7Atd/Hzr58Nb/kG59lZvKOPsuVA3PwO89rJ1QCUFNFxkmmpyjz8IfFRCRc9UStv5SxWA0IReSqC/H+rTj0F8ZmbkzWF+UNhJKETJfFDsjtdTYgDI+7nBXsAqbZ8tFvviF5jCbw9Pzfab4Y+W6aUf6pMPQ7l0KZTvrnCEMgJwVW4E0HCRmasguzuBY8cGGqXWmMkd/TJjFvDWXiCbA1IKAEU392hdmOSuv0BWhV9uYrChsJNQiJL5oNgdr6fEAIIzUdXNz+i36+qg1Jcx5xCx+PZCsDs3PXHqY6otUZUSGFOtJRrZM4D37ILa9Vco193umL+68VFP5Ual3gppVdtanfb1kacCn3wM07k+cQrwH/utLGMhrON1dzrHNW2GyyZvLURxuh4UdhIKUTIfFLvjdVcRFA1NvgWj5MkTwDOPOaMuYvaDDxtPr9SgOPVjR6FvmTVh1RONlGtvhbr+XktM9W5Ydgerp+ibEFoZCWNH77a9f3LS+rdQgPf2A7a6QKkp05F79+/aY1XLaDXKSqBuntfZqi9EpX4vlfj7oLCTUIia+aCYHW++tHH16T8Br2137hgNbDW/k4whUn3nNJidh4znPb1SXXHqRr0YjKkGqqp847oxrgaATUxdQi0ydRANTZAv6IlGqupwcmLOWcGDHz0G+OhDx1OpydOQe2+/LbZdQJx/EcS5jf5lCPSFqFCiYJaksJPQiKP5wMA3bdweM+35gHDU/DY+E5WFLSzsIvXBnx+B8sPVLru5q1eqLU4d0zJWfLkh/MeOamLf1gqxrVXLEq4e5zyp+zE0e73ZIcvu5Mxlgerx+oLiUzrAJeoQAqO/eQX6ThkFbH4GZvTNhE9pr7dt8h5DKEXt2KNglqSwE+KDx6xgIIS2U1/6PUfNbzOSQ1Ur7jwOE2dsudP05KiSCUC+9TeI2XNNp6a7XowRH+62w4slyyFTaUDNAUrKYT83sN9VeRpeT89ohpY3u4B339I/IHTri8tZOqEWAKDUN0Ld2qKdUyg2U5JtITfs8XqMu5wyPbTM1MGGwk4STyk7aY9dF9AE4IKveuzvWkRNn/W+iMeeF4NDpNJV3loti5dB/v7ftR30rr9A7Wi3nJ8F1osxdvrusgHu62W/U7Dv/rW7AldBOCn9d/G97+ODn66E+Kf/BqvxtbDmqijWZ6S0FocirmkUzJIUdpI47EIOFNYqz42Z/m537AkAEz7lNL+YETXGBwentnulsIvUuHMacNRd3fHYUaf/Qc1ZNdx9BE7u36fXVtdFWK8X4ygbEHC97K9Bf826o/Ap4nXRZRAnj0O++Jy2MzfI9mvZqoaA57JmqKVY+j0ttFKPoTdj3Iu8ppU2S1LYSaLwOPTqG0u2d4pMHZTFy7TjFdIsQ3HWNyll7FG00RsiNaK21tErFDDmbyuy5Sp/645xl+vv1URdKMD8cyDmnm3GtItM/prlvq+Z2ap98PDefijfvwHScIq2bdIEPl1l1W83Qy1fhbp3j7aQ/OgOZwx9BK/JQFDYSaLw/PglBmyVlw/T3KD30nT/uMO67Y5CJEUpiEwdlFW3aT4GibyLmqd4W/U4bzMOV+ip0cA6yLRjhkw++TCwZ5fzhK+9bDa/FjWfApYsh9zXjVNOPQV9U6ZDmJ97FUa8utrWitSVKzwmIHsp4Dhcl1CEfdeuXVi3bh1UVcWiRYtw2WWXhXFYQorGUyCsvhGobyxZeM1dpq2Xpp+4l/tjj0IkRakUOn+PMAtA+mWeuhtYt1mFwPwWUZGpg3LpUqhdf3X6RFRVC1fd/Yr23erx6h+rOWDTn4H6RdrO3ficlNo56/NUnYzJolu2sKuqivvvvx833XQTampq8OMf/xhf+MIXMHXq1DDGR0hR+Np19V1bPoLMIA7B7e+H3NZacpmBfG3hohBJMRSI+kZzZw9AC2E0zCGva+YQsXgZ0POetbu3917NeJtpmN+jX67Bq9utrNSc3RegVY+UVVXA3C8Au17Sn7d8BGb9oLe7rXITMVl0yxb2N954A6eddhomTpwIAKivr8f27dsp7KRieGqX5HHG5evPCdiiY7IqAAm5pQUyXx1xW4ErwzzR9/VvQn74Yd62cPmSpOJo43XjdnxC3xWbZpTXX7XqyBjOS8Nk41ro/L5vT1kB880+z1kv6nHw44CqEc7oHb88hhg5xssW9sOHD6OmpsZ8XFNTg717vfWPW1pa0NLSAgBobm5GbW1tuacOlXQ6HbkxlUJS5gGEM5fjL7yJYznLxDHq3Tcx+osN6OvswAf/9hNttyicCS/GewAAtQ04uug/4+SzT5gREo7Xfc4x4pU2fNz6Z32XB3yw5VmMWNiAPntkhvs8+rlge2wfo0xXYfwtd2FEnSUqJ555HJ+89DxGfukrGHXx4Js/S7kexhgx8hT0+VwH1Dagb+xYfPDTld5rASD92c/hH67+H455+13T7Kmn4ONSJpVKYfzXvwl8/Zvo3/0KquachRF1c3H8T7/DMdUZKjli/kKM/qd/cYwlqpQt7NK3mbC3m3BTUxOamprMxz0BzX0rRW1tbeTGVApJmQcQzlzk1Bla2Bq03diJqTNwsqcH6stbNOGVqp6IogBCOt5jHuOseqB1o7mjOy5SOP67f7dMPa5zfHzyY1PUAQCqir7tW7RzSF0s/M7j2p07xpjtx5GXt5jNpNXNz0A++GsAQN+ul3Hs+HEzYWqwKPZ62McIQIsaUhTv3GsnQfnhak9rOwDIvvk3fPjhhxD278nnmmLqDM1uHtQCzw8hgHMXWSGcX56EkwDkX7ZAffvv2liNkMh0FbJf+6/aeyv4+5o8eXJB7ytb2GtqatDb22s+7u3txfjx48s9LCGhEOhwc9u09YQXP5OHb9ajy2zj6cvZtskpMlJCXPBVbSN69Ignbd669c9BplJakk8eu7ts3+r8fPtWYJCFvVjcY8T0GcDpGfjs+xzms9y+bivd36cXaeA1XbIc8qG7/WPaR40CTpywnVAA6SpPlqvDZKQowJlfhBg7PvS2i4NN2cKeyWRw8OBBHDp0CBMmTEBbWxtWrlwZxtgICQW/qI2iK0Dqx3CnybudeoAmDuK8RZAH3tbLwkJLxDm3EQJWwlRg2zejwuGVKwLHKM4+D9IW3ifOPs/8d1Ts8u4xYvY8YNNTkNmsI9LF07WqvhHqttb8zaL9InGOHfVUEDD52GaoEQKYfw6Ur33b6yt58mHLUSoB8ZlZUC75TulfQoUoW9hTqRSuvvpq3HbbbVBVFRdeeCFOP/30MMZGyKCSL8IiSBzzxVkbn7M7CcUV12C0zOHE1Bn5Fwbh1CRjVxsUSqhccDFUwIyvt9etiUponnLBxci9/x6wsw1YUA9x6igrvFGfO+CfGWxf0OT+fVCffNgxTz/E7LlalItxp2TfudsjZqTUQiC/9m3rKfN700U9Ro5SP0KJY1+wYAEWLFgQxqEIGVLUzc+Y3XmMTjly/z79Oc0sYi/4ZUZyBMRZ+1U8HP3P15j25CDziji3EXLrJmtBONdbCMuNcsHFHvNLlOLhZXcnsOkpbTybngIWL/MkiwWN17xDstnp5Z5dUIFAcfc1mRkLMOCMcHEVNHOUJhAC+Nx8rcRAjMwvdph5SoYtsrtTr+utO8j0zENsedba4eVykL+/B3LKdO0zxi6+5lOQuZxXQMdU69ttvXHDmGoc/9PvIPUde94km1W3DVgIa8D59L4PpFKAiorvOH3L+vqFdOaJ3y/Wl+AwiU2Zbp5LbWsFNj9tvdFVS9+T2BZjUQco7GQYI7s6nLfoQgGOfuBNdJGqJgzbWp1NJdy7TyNLVepdeRb9I+T6e7XQvJSzObafaBRaCMs3kcrl9BPnf7XiDr+gEgCFOELN1/P4EvLhsdsDmt0+22+WXS5mHHGDwk6GLY4CVooCXHQZ0Pma943pKmf6e8Du07KfS+2/t7utrEpZeBu9IPNEXsG3f0bCrEJZSUeqxzRiK/blfl/QwhXkS8hH0PeUt0olgv0ZcYTCToYtvjZZd5XAM78IRXeymR18XLtPw/kqT55wmmGmZazCVFJCvvU3s7xBvuqBQXb4fPZz39rnEXCkGucrZByBgnzBxaZZxe6oDmIgu71flcqkQWEnwxpPGKPjRQHxmVmmkAT2RV17kxUiBwGkFDMu3hF9pzeigBCa7d5oopzLOYSs4Nh7V3nc4DuIoXekeppcFzCOoLo8xS5Qw6XuTj4o7ITA1ebNiIzw6xjkEhRHKVrtGe3zx45qO3g3juqD0grJc+/AS4i999ivKyRwvk2ugxKt7KaXgLo8xUb6DIW9PCq5AkFQ2AmB1ywTlIXq+ZzfgqCkIN/6G/Da9gFOqgApRYvKKVB4i7ED20Mz/bI9w8SxQ9/Wapm08kXD+BVPq1/kyDo1x64o2tqppCB7389rksnnYJZdHeg7pwHwsbEXM9dindtDDYWdDFvcP0I/Jx7gdfZ5Pqs34sC0jNaKbUsLsOtlBKdBAoAA5i0Eqsf5lvENbY779wFbnoNUVUe8fajnsAtdKuWMKtI7KuW927HtxB1Zp0LRcgVyOS2Ec+5CoKMdePFZR9auMYZ8JR/sY/zgz49A+eHqkr+HUpzbQw2FnQxLCu2t6fcD9QiZlJqY7d0DeW6jLmw2UTeaJttJp4GOHVpiVJ6EpHJ2gH5x+oNhZ3cIXdY17/pFwedzxfzbOyLJrg4tq3fzs9pxVQB9n+iVH4MFVfuurUqdZnbr/7U5xsv8Hkpxbg81FHYyLCm6t6Zd2H2FTBMSIWCZZoQA5p2DESNHoO8vm62TnzFTi5h58dm8IlDuDtAvTn8w7Oyevq+2BtB+RbYcO2s95l8sXubxL8juTmckkr1PKYS2MMB1PSB0s41WqRNjqqGuudFZkM32PZSycJbi3B5qKOxk2OBx1AVFmAzwA80nZOLcRmBaxixTgN07MfJf/if6drRpO2clBWXxMgCwTA4+dmNPQaoSdoDuOH13Uk5YuIUO8IZwAgjeWevOZuM9jjh2l98Di/4ReO5x7U5n/b1aq0L39fhPC8yKjJrou0r56s03/MpJFGorL8W5PZQI6VdQfQg4cOBAJU4bSFLqmCdlHkC4c3EX5zK7HAX8CAf6Ydtfdx8n9+DdlgNQCIw453z0jRzlsaXL7k6ty9KWFm0R0Mcl9+/TTCg5WxehdNWAO3a/MYfpzCv3eqgbH4V8/Pcwa+DbdtbG9XBfI9PUYuy6FUWLJoLUFqv/cgWUS74T+F0C8NmxC+D8i4Gtz1n11oWAuOy7pvPXfZyoRL4MWT12QuKAr6Puku8E/mAHij7xS40HNCFF2yaYNnYpNTOMEJDpKoctXWTqLHOJ0dD56T8Br+2w7OIBBancgu0bYqhH9kSl7KxfDXy5r9uM2Al0ShqlAADtuxKKp/qi57u0XWNl1W1Qn34MeFV3aFeN0G8YXOUkxlTrOQm2JLWY9Dh1Q2Enw4Khsn/Krg5rF+h4wd+cImbP1RKVcrpZ4rXt3hKz0zLBzltjl++w+/frTlNZUDieuduVgJieKTjUs1g8mb77us2oF9nWGhjv7rEpzF+oJY4FlVR22eBFpg6p799gzn3cOQ348MMPrUbaupkKx446k9RiXLqXwk6GBWHaP/OZNzxx7eYL/iIhMnUQDU2QL+iVB1VVC+2zLw5P/wnqp06zaq777WztCxdsfUOz/VCffBjKpUu1w+tZslJRIK74HsSU6Q5ThdTHKgsw/ZSCo8SAma0LT7y7vbaMUt8IdWuL5aNwNchwMOcsbXGUmg1eBRwLlVFSQPT0eP4eZHen01ZfvwhKfbw6JxlQ2MmwoZjkniAGilSxJwUZLfD+4fPzcOw/z+E1AAALjElEQVQ/DgYvKNMyzsdfOB/4aztw/CPrvFueNeulYEy1Z2fr2Q2vv9daXPa8CnXvHoj6RktM1Zy2q2+4yMe5WJqztlC82brOcEfAW1tGue72AX0ensUi2685R6X3zgXwN6dFxflZLhR2kmjCzgQsOFZZL/GLdBrpr38TytkNwQc9dhSAgClIL7/gtT+8tRfqnT/WHIeGDf3tbudNQcZZi1x98mFgz6swo2okdOejfjegqtpz6Sqvc3EQTRDuLlSiocnpVC7SH+L4jH2xMExcJUYVxRkKO0ksg5EJWIit3p1S37/7FeDLwSnsjpZuRvifH4Z5JpfVRL2t1dM71Txmpg7KpUuh7t1jNY+ob4ScntHt76rVzLm+cUhs7I6x5at5U4I/xG+xwLSM1UWpgONEKXO0XCjsJLEMRibgQKIkuzu1UDkDJYWqOWfhZNAYjTsKoyyBvfmzord0czf+EApk118di4cjwcp2l+Ieq8jUOToL2ZOChpJ8ZrFSTCJBn/GbaxBRyhwtFwo7SSyDFQmTT5Qc2Z56Sv2IurlAT4839t0VLy0WL7MyK4Xm3ARg7bBtlSDx3ru2ASlWByfbMU37tB7uaJxfnjyhNQEZU205DQdILhpqSvGHBCUNFZxRGqHM0XKhsJPEUglnmFscjJR6x22+UYfdbtfOVwXRbS93M/dsS9DtZYH96qnYYrTlnl3Ivf+e2XA6qD58kpDdnTj+wptmD1o7dJ4SEhPCiIQp9ny+4my/zVd97Oj6rjswVd2wl7vDKAGgox3SLuiAxwFqORdd7GzzH1fMTRF+GIubuwetnaH+exksKOyEhIyvOLvjzKUriWmAgun2BUO+tRfY9ZL1olvU01WeGGzz/O7Wf5+ZDXzQa5XJFaKo+vBxwlksLHkLlx0KOyFDgCPO/PD7kC/otWQM1FxBnYFMs8runZYYA5a4CwViyXJP02fHwnDgHSukcmebo/wAEA0b+2CQJBv6QFDYCRki7MJsprMbXZcKFBp7FA2OHdWSkV58Fvj7Xv0kMCslBp1f3fgopBE3b9j2bfVk4iDopZbbFYuXoapjO/rnLozFPEuFwk7IEOPOEnXHjedr7WavfigWL7MyTIGCF4jARhERaes2EKXGm8vuTsj196IvlwV27/KUG0gSFHZCKkCQky5vZyd3nHX7Vv2xDKwC6Xd8+44/qDpklCNiSo03d3xO5i83EHeUSg+AEGJhZq3aRMtAzJ6r1YhRFBgdhczH6aqCRF1dexPk47+HXH9vcNSO67xRw/M9FGgrd3xOUTwlfpMEd+yEhEAYZgy/rFVHZyefUMqwMivj5FgsNd7c+Nyod9/EcZEqqtxA3KCwE1ImYZkx/LJWB4qzzpsF61ps8rYDdFeH1HewUTVPlBpvLjJ1GP3FBpzs6SlqUYwbFHZCyiSsGiNBWasljSlgsclbfMsomavXZ5fpKiirbkuc6BkkJRnJDwo7IWUSlhkjzJT2oMVmIDFztKHL9kNuazWdq0Gp+CR6UNgJKZMwBTmsXWSpi427WoGUhaXik2hBYSckBKJ2W1/oYuO2w3va0NU3DqtU/EKIQ7w/hZ2QCjHYAuEoQbDxUWfruc3PaBmr77zlKPErMnW+bejiEjETBvmuS1zi/csS9gcffBDt7e1Ip9OYOHEiVqxYgdGjR4c1NkISy1AJhN955P59kA/+2vnGPHZ4e5jgiYTb2Ae6LnFpxlFWgtK8efOwdu1arFmzBpMmTcKGDRvCGhchiWaoEoL8ziPbt3rfOMBOXGTqMPrb/xxJETPuSGR3Z97nCjrWANel1OSooaasHfv8+fPNf8+aNQsvvfRSnncTQgyGKiHI9zxjqiHtLfjO/BKUr30rkqI9EH47bAAl3w0NdF3i0owjNBt7a2sr6uvrwzocIYlmqATC7zwiUwcVgGzfCnH2eZ4Sv3EicIddormkkOsSNUe5H0LKoJboGqtXr8aRI0c8zy9evBgLFy4EADz22GPo7u7GqlWrIAIaBrS0tKClRUuXbm5uRl9fn+/7KkU6nUbWr8NMzEjKPIDkzIXzGDz6OjvwwU9XarH36SqMv+UuAPA8N6LO2nkb8+jr7ED/7ldQNecsx+tRZsSIEQW9b0BhH4jnn38ezz33HG6++WaMHDmy4M8dOHCgnNOGTm1tLXp6eio9jLJJyjyA5MwlCfOQ3Z2RdZ76RbHki2ypra3F+3/Z4iiBHNXoFjeTJ08u6H1lmWJ27dqFJ554ArfccktRok4IiQ9RT1AK6hObN8M2JtEtpVKWsN9///3IZrNYvXo1AGDmzJlYvnx5KAMjhESDJCYoxamaZSmUJey/+tWvwhoHISSiJFEE4xLdUirMPCWE5CVuCUqFZvTGIbqlVCjshMSMfD1RB2sHaq9jHmXikvI/2FDYCYkRQcJFQdNIulO0UNjzlJAYEZSQE6eepYNJXFL+Bxvu2AmJEUGOzLg6OMM2HyXdKVooFHZCYkSQcMVR0AbLfJRkp2ihUNgJiRlBwhU3QRsMe3gcmmAMBRR2QkhFCNt8RAeyBYWdEFIRwjYfMSLGgsJOCKkYYZqP4upAHgwo7ISQRBBHB/JgQWEnhCSGuDmQBwsmKBFCSMKgsBNCSMKgsBNCSMKgsBNCSMKgsBNCSMKgsBNCSMIQUkpZ6UEQQggJD+7Yda6//vpKDyEUkjIPIDlz4TyiRVLmkQ8KOyGEJAwKOyGEJIzUz372s59VehBRYcaMGZUeQigkZR5AcubCeUSLpMwjCDpPCSEkYdAUQwghCYPVHX148skn8dBDD+G+++5DdXV1pYdTNOvXr8eOHTsghMDYsWOxYsUKTJgwodLDKpoHH3wQ7e3tSKfTmDhxIlasWIHRo0dXelglsW3bNjz66KPYv38/br/9dmQymUoPqWB27dqFdevWQVVVLFq0CJdddlmlh1QSd999N3bu3ImxY8di7dq1lR7OoMIdu4uenh50dHSgtra20kMpmUsvvRRr1qzBnXfeiQULFuCPf/xjpYdUEvPmzcPatWuxZs0aTJo0CRs2bKj0kErm9NNPx6pVq/C5z32u0kMpClVVcf/99+OGG27AL37xC2zduhXvvvtupYdVEl/5yldwww03VHoYQwKF3cUDDzyAK664AkKISg+lZEaNGmX++5NPPontXObPn49UKgUAmDVrFg4fPlzhEZXO1KlTMXny5EoPo2jeeOMNnHbaaZg4cSLS6TTq6+uxffv2Sg+rJD7/+c9jzJgxlR7GkEBTjI0dO3ZgwoQJOOOMMyo9lLL5wx/+gM2bN2PUqFH46U9/WunhlE1rayvq6+srPYxhx+HDh1FTU2M+rqmpwd69eys4IlIIw07YV69ejSNHjnieX7x4MTZs2ICbbrqpAqMqnnzzWLhwIZYsWYIlS5Zgw4YNePrpp3H55ZdXYJQDM9A8AOCxxx5DKpXC+eefP9TDK4pC5hI3/ILm4noHOJwYdsL+k5/8xPf5t99+G4cOHcJ1110HAOjt7cW//uu/4o477sC4ceOGcogFETQPNw0NDWhubo6ssA80j+effx7t7e24+eabIy8ohV6TOFFTU4Pe3l7zcW9vL8aPH1/BEZFCGHbCHsS0adNw3333mY+///3v44477ohlVMzBgwcxadIkAJp5KY62XUCLxnjiiSdwyy23YOTIkZUezrAkk8ng4MGDOHToECZMmIC2tjasXLmy0sMiA8AEpQDiLOxr1qzBwYMHIYRAbW0tli9fHstwxx/84AfIZrOmw2vmzJlYvnx5hUdVGi+//DJ+85vf4OjRoxg9ejTOOOMM3HjjjZUeVkHs3LkTDzzwAFRVxYUXXohvfetblR5SSfzyl7/Enj178NFHH2Hs2LG4/PLL0djYWOlhDQoUdkIISRgMdySEkIRBYSeEkIRBYSeEkIRBYSeEkIRBYSeEkIRBYSeEkIRBYSeEkIRBYSeEkITx/wE5Wtj89zVuYAAAAABJRU5ErkJggg==\n", "text/plain": ["<Figure size 432x288 with 1 Axes>"]}, "metadata": {}, "output_type": "display_data"}], "source": ["plt.plot(fact[:,0], fact[:,1], \".\");"]}, {"cell_type": "code", "execution_count": 30, "metadata": {}, "outputs": [], "source": []}, {"cell_type": "code", "execution_count": 31, "metadata": {}, "outputs": [], "source": []}], "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.7.0"}}, "nbformat": 4, "nbformat_minor": 2}