FactorAnalysis.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "FactorAnalysis.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "toc_visible": true,
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/rkn2/factorAnalysisExample/blob/master/FactorAnalysis.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "60YyXAwxgeo6",
        "colab_type": "text"
      },
      "source": [
        "# Background on the problem\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "QZkl03-xlDMc",
        "colab_type": "text"
      },
      "source": [
        "A wall in a building next to a river is deteriorating. What is causing it?\n",
        "\n",
        "Below find an image of how the wall is situated in the building, where doors are, and where current sensors are located. \n",
        "\n",
        "**Your job is to look at the data coming from the sensors, and figure out what they are telling you.**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "sXxYuxEzgk5T",
        "colab_type": "text"
      },
      "source": [
        "![alt text](https://drive.google.com/uc?id=1TWhPuLFF7pNcXlycMDmdDWlFxWQJq3DQ)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "nYeZzHCmllHU",
        "colab_type": "text"
      },
      "source": [
        "# Pre-processing our data\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dkYopb1bx4G7",
        "colab_type": "text"
      },
      "source": [
        "## Loading the data"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "bANTW6kuUClT",
        "colab_type": "text"
      },
      "source": [
        "This step just installs any packages we will need to run the analysis code. You don't have to worry about any specifics here."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "tP4cR7ZrSuLc",
        "colab_type": "code",
        "outputId": "a7511d7a-e5ea-4921-e015-2963b423a920",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 207
        }
      },
      "source": [
        "# Import required libraries\n",
        "import pandas as pd\n",
        "from sklearn.datasets import load_iris\n",
        "#!pip install factor_analyzer\n",
        "!pip install factor_analyzer==0.2.3\n",
        "from factor_analyzer import FactorAnalyzer\n",
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "\n",
        "\n",
        "!pip install -U -q PyDrive\n",
        "from pydrive.auth import GoogleAuth\n",
        "from pydrive.drive import GoogleDrive\n",
        "from google.colab import auth\n",
        "from oauth2client.client import GoogleCredentials\n",
        "# Authenticate and create the PyDrive client.\n",
        "auth.authenticate_user()\n",
        "gauth = GoogleAuth()\n",
        "gauth.credentials = GoogleCredentials.get_application_default()\n",
        "drive = GoogleDrive(gauth)"
      ],
      "execution_count": 1,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Collecting factor_analyzer==0.2.3\n",
            "  Downloading https://files.pythonhosted.org/packages/79/1b/84808bbeee0f3a8753c3d8034baf0aa0013cf08957eff750f366ce83f04a/factor_analyzer-0.2.3-py2.py3-none-any.whl\n",
            "Requirement already satisfied: scipy in /usr/local/lib/python3.6/dist-packages (from factor_analyzer==0.2.3) (1.3.1)\n",
            "Requirement already satisfied: numpy in /usr/local/lib/python3.6/dist-packages (from factor_analyzer==0.2.3) (1.16.5)\n",
            "Requirement already satisfied: pandas in /usr/local/lib/python3.6/dist-packages (from factor_analyzer==0.2.3) (0.24.2)\n",
            "Requirement already satisfied: pytz>=2011k in /usr/local/lib/python3.6/dist-packages (from pandas->factor_analyzer==0.2.3) (2018.9)\n",
            "Requirement already satisfied: python-dateutil>=2.5.0 in /usr/local/lib/python3.6/dist-packages (from pandas->factor_analyzer==0.2.3) (2.5.3)\n",
            "Requirement already satisfied: six>=1.5 in /usr/local/lib/python3.6/dist-packages (from python-dateutil>=2.5.0->pandas->factor_analyzer==0.2.3) (1.12.0)\n",
            "Installing collected packages: factor-analyzer\n",
            "Successfully installed factor-analyzer-0.2.3\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "i5os7BeoT8zV",
        "colab_type": "text"
      },
      "source": [
        "Now we will load the data. This line reads in the comma separated value sheet that I made in excel.\n",
        "\n",
        "If you are given sensor data in excel and want to export it to csv in the future, see this link:\n",
        "https://www.ablebits.com/office-addins-blog/2014/04/24/convert-excel-csv/"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab_type": "code",
        "id": "kX8bxB9mpiYp",
        "colab": {}
      },
      "source": [
        "url = 'https://raw.githubusercontent.com/rkn2/factorAnalysisExample/master/testDataSmall.csv' #you can change this link depending on your data\n",
        "df = pd.read_csv(url)\n",
        "#df= pd.read_csv(\"testDataSmall.csv\") # you can also just upload data using the interface on the left and use this reader here."
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "r_MFE9oIs5KZ",
        "colab_type": "text"
      },
      "source": [
        "## Cleaning the data"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ed6U4YqqVfvn",
        "colab_type": "text"
      },
      "source": [
        "Now we can start to process our data and remove any columns we don't want. We don't need the 'true ground water' 'wind driven rain' 'corSet', or 'notCor' columns since those were only needed to generate the synthetic data.\n",
        "\n",
        "Additionally, we want to remove anywhere that there is missing data."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "cbE21AiDU_5L",
        "colab_type": "code",
        "outputId": "e09a408f-a71c-4a5b-9467-6f991729ce61",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 527
        }
      },
      "source": [
        "#index columns\n",
        "df.columns\n",
        "\n",
        "#drop unnecessary columns\n",
        "unnecessaryColumns = ['corSet', 'notCor', 'true ground water', 'wind driven rain'] #this list is the only thing you would need to change\n",
        "df.drop(unnecessaryColumns,axis=1,inplace=True) \n",
        "\n",
        "#drop any missing value rows\n",
        "df.dropna(inplace=True)\n",
        "\n",
        "#view information if you want\n",
        "df.info()\n",
        "#df.head()\n",
        "\n",
        "#calcualte the number of variables\n",
        "numVars = df.shape[1]-len(unnecessaryColumns)\n"
      ],
      "execution_count": 26,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "<class 'pandas.core.frame.DataFrame'>\n",
            "Int64Index: 125 entries, 0 to 124\n",
            "Data columns (total 25 columns):\n",
            "water pH 1                  125 non-null float64\n",
            "water pH 2                  125 non-null float64\n",
            "ground water level 1        125 non-null float64\n",
            "ground water level 2        125 non-null float64\n",
            "soil moisture 1             125 non-null float64\n",
            "soil moisture 2             125 non-null float64\n",
            "soil moisture 3             125 non-null float64\n",
            "motor monitor 1             125 non-null float64\n",
            "motor monitor 2             125 non-null float64\n",
            "salinity I                  125 non-null float64\n",
            "salinity E                  125 non-null float64\n",
            "air RH I                    125 non-null float64\n",
            "air RH E                    125 non-null float64\n",
            "air temp I                  125 non-null float64\n",
            "air temp E                  125 non-null float64\n",
            "wall RH 1 I                 125 non-null float64\n",
            "wall RH 2 I                 125 non-null float64\n",
            "wall RH 3 I                 125 non-null float64\n",
            "wall RH 1 E                 125 non-null float64\n",
            "wall RH 2 E                 125 non-null float64\n",
            "wall RH 3 E                 125 non-null float64\n",
            "rain intensity              125 non-null float64\n",
            "rain duration               125 non-null float64\n",
            "wind speed                  125 non-null float64\n",
            "wind direction (radians)    125 non-null float64\n",
            "dtypes: float64(25)\n",
            "memory usage: 25.4 KB\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "gn7I2_5StNZ1",
        "colab_type": "text"
      },
      "source": [
        "# Check that factor analysis applies\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "dDVz8CyFWwh5",
        "colab_type": "text"
      },
      "source": [
        "Factor analysis code based on: https://www.datacamp.com/community/tutorials/introduction-factor-analysis"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "fPO4G7tHyBbg",
        "colab_type": "text"
      },
      "source": [
        "## Bartlett's test"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "NrqV1ZNGWIc7",
        "colab_type": "text"
      },
      "source": [
        "\n",
        "\n",
        "The first way we can do this is using *bartlett's test* which checks whether the data is *statistically significant* by checking whether or not the observed variabels intercorrelate at all using the observed correlation matrix against the identity matrix. If the test is found statistically insignificant, you should not use factor analysis. "
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "aXtTewA1WHwR",
        "colab_type": "code",
        "outputId": "2a3504c5-27b0-4eb5-90af-ffa07f6fa586",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        }
      },
      "source": [
        "#this is not something you need to change\n",
        "from factor_analyzer.factor_analyzer import calculate_bartlett_sphericity\n",
        "chi_square_value,p_value=calculate_bartlett_sphericity(df)\n",
        "chi_square_value, p_value"
      ],
      "execution_count": 27,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "(2408.015065762507, 0.0)"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 27
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "HdtI64_7WuSW",
        "colab_type": "text"
      },
      "source": [
        "In this Bartlett ’s test, the p-value is 0. The test was statistically significant, indicating that the observed correlation matrix is not an identity matrix."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "p0IHvEg4yE9V",
        "colab_type": "text"
      },
      "source": [
        "## Kaiser-Meyer-Olkin test"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "jOOzLoLEWohJ",
        "colab_type": "text"
      },
      "source": [
        "The second step we can do is the *Kaiser-Meyer-Olkin (KMO*) Test which measures the *suitability of data for factor analysis*. It determines the adequacy for each observed variable and for the complete model. KMO estimates the proportion of variance among all the observed variable. Lower proportion id more suitable for factor analysis. KMO values range between 0 and 1. Value of KMO less than 0.6 is considered inadequate."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "gJbu2WwCXBa0",
        "colab_type": "code",
        "outputId": "2acae92c-439c-4cf8-cc0f-46d8dd0fd297",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        }
      },
      "source": [
        "#this is not something you need to change\n",
        "from factor_analyzer.factor_analyzer import calculate_kmo\n",
        "kmo_all,kmo_model=calculate_kmo(df)\n",
        "\n",
        "kmo_model"
      ],
      "execution_count": 28,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "0.794519644398728"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 28
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ZK-I2yJrXEgB",
        "colab_type": "text"
      },
      "source": [
        "The overall KMO for our data is 0.79, which is excellent. This value indicates that you can proceed with your planned factor analysis."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "QmEB220wyL85",
        "colab_type": "text"
      },
      "source": [
        "# Applying factor analysis"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "KuLv3yjyyPWr",
        "colab_type": "text"
      },
      "source": [
        "## Choosing the number of factors"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "mtR0FtQGXHRA",
        "colab_type": "text"
      },
      "source": [
        "Now we need to choose the number of factors. For this, we can use the Kaiser criterion and scree plot. Both are based on eigen values. "
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Z6zPLCYTXNgL",
        "colab_type": "code",
        "outputId": "8e49488c-93c1-4648-e5cc-f4bac068d63e",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 824
        }
      },
      "source": [
        "# Create factor analysis object and perform factor analysis\n",
        "fa = FactorAnalyzer()\n",
        "fa.analyze(df, numVars, rotation=None) \n",
        "# Check Eigenvalues\n",
        "ev, v = fa.get_eigenvalues()\n",
        "ev"
      ],
      "execution_count": 29,
      "outputs": [
        {
          "output_type": "execute_result",
          "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>Original_Eigenvalues</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>0</th>\n",
              "      <td>5.952378</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1</th>\n",
              "      <td>3.157398</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2</th>\n",
              "      <td>1.804729</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>3</th>\n",
              "      <td>1.637225</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>4</th>\n",
              "      <td>1.390228</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>5</th>\n",
              "      <td>1.266333</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>6</th>\n",
              "      <td>1.162993</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>7</th>\n",
              "      <td>1.092643</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>8</th>\n",
              "      <td>1.043082</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>9</th>\n",
              "      <td>1.004171</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>10</th>\n",
              "      <td>0.920566</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>11</th>\n",
              "      <td>0.827604</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>12</th>\n",
              "      <td>0.725411</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>13</th>\n",
              "      <td>0.632156</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>14</th>\n",
              "      <td>0.595860</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>15</th>\n",
              "      <td>0.481387</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>16</th>\n",
              "      <td>0.412096</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>17</th>\n",
              "      <td>0.288582</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>18</th>\n",
              "      <td>0.222748</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>19</th>\n",
              "      <td>0.191830</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>20</th>\n",
              "      <td>0.118418</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>21</th>\n",
              "      <td>0.021652</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>22</th>\n",
              "      <td>0.019398</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>23</th>\n",
              "      <td>0.017750</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>24</th>\n",
              "      <td>0.013361</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "    Original_Eigenvalues\n",
              "0               5.952378\n",
              "1               3.157398\n",
              "2               1.804729\n",
              "3               1.637225\n",
              "4               1.390228\n",
              "5               1.266333\n",
              "6               1.162993\n",
              "7               1.092643\n",
              "8               1.043082\n",
              "9               1.004171\n",
              "10              0.920566\n",
              "11              0.827604\n",
              "12              0.725411\n",
              "13              0.632156\n",
              "14              0.595860\n",
              "15              0.481387\n",
              "16              0.412096\n",
              "17              0.288582\n",
              "18              0.222748\n",
              "19              0.191830\n",
              "20              0.118418\n",
              "21              0.021652\n",
              "22              0.019398\n",
              "23              0.017750\n",
              "24              0.013361"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 29
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "AYqhbUM5XQLR",
        "colab_type": "text"
      },
      "source": [
        "Here, you can see only for 10-factors eigenvalues are greater than one. It means we need to choose only 10 factors (or unobserved variables). But this is only a first step metric, next we can generate a Scree plot to see which of those 10 are really distinct from the pack."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "8EcInPmvXS4E",
        "colab_type": "code",
        "outputId": "fcf80af1-9ab0-42f6-f3ea-b017d8d7d045",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 295
        }
      },
      "source": [
        "# Create scree plot using matplotlib\n",
        "plt.scatter(range(1,df.shape[1]+1),ev)\n",
        "plt.plot(range(1,df.shape[1]+1),ev)\n",
        "plt.title('Scree Plot')\n",
        "plt.xlabel('Factors')\n",
        "plt.ylabel('Eigenvalue')\n",
        "plt.grid()\n",
        "plt.show()"
      ],
      "execution_count": 30,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAEWCAYAAABliCz2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl8nGW5//HPNdnXpknTtA1twtKG\nVrbSAmrBpigWhSPFnwv+jp6jR+V4jnhcEFnccEOOiNvveBbEBY5KcYGKKCBqWwoq0A260VLomu5r\nkjZps1y/P2aSpu1MMmnyzCTzfN+v17wy88yzXDdDr3nmeu7nvs3dERGRzBdJdwAiIpIaSvgiIiGh\nhC8iEhJK+CIiIaGELyISEkr4IiIhoYQvEjAzqzUzN7PsdMci4aaEL8OOmV1qZn8xs4Nmts/MnjGz\ni9IcU72ZdZpZs5k1mdlaM/vAKezndjP7aRAxiuiMQ4YVMysFHgX+BfgFkAtcBhzp536y3b19kMPb\n5u6nmZkB1wC/MrNngcODfByRU6IzfBluJgG4+wPu3uHuLe7+B3d/sWsFM/uwma2JnWmvNrMLY8s3\nmtnNZvYicMjMss1snJn92sx2m9kGM/u3HvuJmNktZvaKme01s1+YWXlfAXrUPGA/MOXE92PHfCT2\n62S9mX04tvxK4Dbg3bFfCi8M8L+VyHGU8GW4WQd0mNl9ZvYWMxvZ800zeydwO/APQCnwNmBvj1Xe\nA1wFlAGdwG+BF4Bq4I3AJ8xsdmzdjwFzgJnAOKIJ/Pt9BRj7org2dowVcVaZC2yN7fMdwB1mdrm7\nPw7cATzo7sXufn5fxxLpDyV8GVbcvRG4FHDgB8Du2NlyVWyVDwHfcPfnY2fa6919U49dfM/dt7h7\nC3ARUOnuX3b3o+7+amyf18XW/QjwWXff6u5HiH6RvKOXi6/jzOwAsAf4IvA+d1/bcwUzGw/MAG52\n91Z3Xw7cS/QLSiRQquHLsOPua4D3A5jZ2cBPge8QPXsfD7zSy+Zbejyv4ViS7pIFLOrx/sNm1tnj\n/Q6gCmiIs+9t7n5aH+GPA/a5e1OPZZuA6X1sJzJgSvgyrLn7S2b2E+CfY4u2AGf2tkmP51uADe4+\nMcG6W4B/cvdnBhzoMduAcjMr6ZH0J3DsC0TD10pgVNKRYcXMzjazG83stNjr8UTP7P8WW+Ve4NNm\nNs2izjKzmgS7ew5oil3ILTCzLDM7p0cXz/8Gvta1vZlVmtk1A4nf3bcAfwG+bmb5ZnYe8EGiv1IA\ndgK1ZqZ/mzLo9D+VDDdNwCXAs2Z2iGiiXwncCODuvwS+Bvw8tu48IG7PGnfvAK4GLgA2EK293wuM\niK3yXeAR4A9m1hQ71iWD0Ib3ALVEz/YfBr7o7n+MvffL2N+9ZrZ0EI4l0s00AYqISDjoDF9EJCSU\n8EVEQkIJX0QkJJTwRURCYkj1wx81apTX1tZy6NAhioqK0h1O2oS5/Wp7ONsO4W7/QNq+ZMmSPe5e\nmcy6Qyrh19bWsnjxYhYsWEB9fX26w0mbMLdfba9PdxhpE+b2D6TtZrap77WiVNIREQkJJXwRkZBQ\nwhcRCQklfBGRkAg04ZtZmZn9ysxeis1A9LogjyciIokF3Uvnu8Dj7v4OM8sFCgM+noiIJBBYwjez\nEcAbiE1U4e5HgaODfZx5yxq464m1bDvQwriyAm6aXcecqdWDfRgRkWEvsNEyzewC4B5gNXA+sAT4\nuLsfOmG964HrAaqqqqbNnTuX5uZmiouL+zzGgZY2Gva30NbptHVCfhZEzKgeWUBZQc6gtylVkm1/\nJlLbw9l2CHf7B9L2WbNmLXH3pGZMCzLhTyc6fvgMd3/WzL4LNLr75xNtM336dO/PjVcz7vwzDQda\nTlpeXVbAM7dcfurBp5luQKlPdxhpEea2Q7jbP8Abr5JO+EFetN0KbHX3Z2OvfwVcOJgH2BYn2fe2\nXEQkzAJL+O6+A9hiZnWxRW8kWt4ZNOPKCvq1XEQkzILuh/8x4Gdm9iLRaeTuGMyd3zS7joKcrOOW\nFeRkcdPsugRbiIiEV6DdMt19OZBUbelUdPXG+fJvV7Pv8FFGl+Rx21snq5eOiEgcw/5O2zlTq7n/\ngxcD8OVrXqNkLyKSwLBP+AATKqL3c23aezjNkYiIDF0ZkfBL83MYWZjDpn1K+CIiiWREwgeYUFHE\nZp3hi4gklDEJv6a8kE37DvW9oohISGVOwq8oZNuBVto6OtMdiojIkJQxCX9CeSEdnU7Dft1lKyIS\nT8Yk/JqK6IzvunArIhJfBiX8aNfMzXtVxxcRiSdjEv7okjzycyLqiy8ikkDGJHwzY0J5oUo6IiIJ\nZEzCB5hQrr74IiKJZFTCr6koZPO+wwQ1qYuIyHCWcQm/pa2D3U1H0h2KiMiQk1EJf0J5bBA11fFF\nRE6SUQm/uy++6vgiIifJqIRfXVZAxNQXX0QknoxK+LnZEcaVFaikIyISR0YlfIheuN2oko6IyEky\nLuFH++KrpCMicqKMS/g1FYXsP9xGY2tbukMRERlSMi/hl3cNoqayjohITxmX8DWhuYhIfBmX8I+N\ni686vohIT9lB7tzMNgJNQAfQ7u7TgzweQHFeNhVFuSrpiIicINCEHzPL3fek4DjdJlQUqqQjInKC\njCvpANRWFLFZN1+JiBzHghxK2Mw2APsBB/7H3e+Js871wPUAVVVV0+bOnUtzczPFxcWnfNyHXz7K\nI6+0cc+bC8mJ2CnvJ10G2v7hTG0PZ9sh3O0fSNtnzZq1JOlyubsH9gCqY39HAy8Ab+ht/WnTprm7\n+/z5830gfr1ki9fc/Kiv39U0oP2ky0DbP5yp7eEV5vYPpO3AYk8yJwda0nH3htjfXcDDwMVBHq/L\nsQnNVdYREekSWMI3syIzK+l6DrwZWBnU8XqaUN41TLK6ZoqIdAmyl04V8LCZdR3n5+7+eIDH6zaq\nOJfC3CyNmiki0kNgCd/dXwXOD2r/vTEzJpQXqqQjItJDRnbLhGgdX2f4IiLHZHDCj/bF7+wMrtup\niMhwkrEJf0J5IUfbO9nZ1JruUEREhoSMTfg1GjVTROQ4mZvwY10zdeFWRCQqYxP+uLJ8siOmYZJF\nRGIyNuFnZ0WoHlmgko6ISEzGJnyIXrjVqJkiIlEZnfBrNC6+iEi3zE745UUcbGnjwOGj6Q5FRCTt\nMjrha0JzEZFjMjrhd/fFVx1fRCSzE/6E8q5x8dU1U0QkoxN+YW42lSV5KumIiJDhCR+gplyjZoqI\nQAgS/oQKjYsvIgIhSPg15UXsaGylta0j3aGIiKRV5if8WE+dLSrriEjIZXzCV198EZGojE/4NeXq\niy8iAiFI+OVFuRTnZasvvoiEXsYnfDNjgrpmiohkfsKH6IVbdc0UkbALScIvYsv+w3R0erpDERFJ\nm8ATvpllmdkyM3s06GMlUlNRSFuHs/1gS7pCEBFJu1Sc4X8cWJOC4yRU0z2Imso6IhJegSZ8MzsN\nuAq4N8jj9GWChkkWEcHcg6trm9mvgK8DJcCn3f3qOOtcD1wPUFVVNW3u3Lk0NzdTXFw8aHF0uvPh\nPxxmdm0O76rLHbT9BmWw2z+cqO3hbDuEu/0DafusWbOWuPv0ZNbNPqUjJMHMrgZ2ufsSM6tPtJ67\n3wPcAzB9+nSvr69nwYIF1Ncn3OSU1CxZAMUl1NdPG9T9BiGI9g8Xant9usNImzC3P1VtD7KkMwN4\nm5ltBOYCl5vZTwM8Xq8maEJzEQm5wBK+u9/q7qe5ey1wHfBnd39vUMfrS015tC9+kCUsEZGhLBT9\n8AEmVBTRdKSd/Yfb0h2KiEhapCThu/uCeBdsU6l7EDWNqSMiIRWaM/waDZMsIiEXmoQ/vlwJX0TC\nLTQJPz8nizGl+Wzap5KOiIRTaBI+aEJzEQm3pBK+mVWZ2Q/N7LHY6ylm9sFgQxt8NRoXX0RCLNkz\n/J8ATwDjYq/XAZ8IIqAg1VQUsrvpCIePtqc7FBGRlEs24Y9y918AnQDu3g50BBZVQCZUFAGwWWf5\nIhJCySb8Q2ZWATiAmb0WOBhYVAGpUU8dEQmxZAdP+xTwCHCmmT0DVALvCCyqgHT1xdeFWxEJo6QS\nvrsvNbOZQB1gwFp3H3ZjFJQV5lKan62umSISSkklfDP7hxMWXWhmuPv9AcQUqJqKIpV0RCSUki3p\nXNTjeT7wRmApMOwS/oSKQlY2DLvLDyIiA5ZsSedjPV+bWRnRMe6HnZryQp5YuYP2jk6ys0J135mI\nhNypZrxDwOmDGUiq1FQU0t7pbDvQmu5QRERSKtka/m+Jdckk+iUxBfhFUEEFaUJ5tC/+pn2Huic3\nFxEJg2Rr+N/s8bwd2OTuWwOIJ3A9h0m+bGKagxERSaFka/gLgw4kVcaU5pObHdHdtiISOr0mfDNr\n4lgp57i3AHf30kCiClAkYowfWaCZr0QkdHpN+O5ekqpAUkl98UUkjJKt4QNgZqOJ9sMHwN03D3pE\nKdDe0cnaHU3U3vI7qssKuGl2HXOmVqc7LBGRQCU7Hv7bzOxlYAOwENgIPBZgXIGZt6yBv76yt7tO\n1XCghVsfWsG8ZQ1pjUtEJGjJ9sP/CvBaYJ27n070Ttu/BRZVgO56Yi1tncdflmhp6+CuJ9amKSIR\nkdRINuG3ufteIGJmEXefD0wPMK7AbDvQ0q/lIiKZItka/gEzKwaeAn5mZruI3m077IwrK6AhTnIf\nV1aQhmhERFIn2TP8a4DDwCeBx4FXgL/rbQMzyzez58zsBTNbZWZfGliog+Om2XUU5GQdt6wgJ4ub\nZtelKSIRkdRI9gz/n4EH3b0BuC/JbY4Al7t7s5nlAE+b2WPuntbaf1dvnC8/upp9h45SUZTL56+e\nol46IpLxkj3DLwH+YGaLzOwGM6vqawOPao69zIk94t3ElXJzplaz6DOzyM2K8PYLq5XsRSQUkkr4\n7v4ld38N8FFgLLDQzP7Y13ZmlmVmy4FdwJPu/uyAoh1ERXnZXHT6SBau253uUEREUsLckz/pNrMx\nwDuB64ASdz8vye3KgIeBj7n7yhPeux64HqCqqmra3LlzaW5upri4OOm4TtVjG9p4cO1R7p5ZQEXB\n0BkbP1XtH4rU9nC2HcLd/oG0fdasWUvcPblek+7e5wP4V2ABsAq4HZiSzHYn7OMLwKd7W2fatGnu\n7j5//nxPhbU7Gr3m5kf9589uSsnxkpWq9g9Fant4hbn9A2k7sNiTzMPJntaOBz7h7q9x99vdfXVf\nG5hZZezMHjMrAK4AXkryeCkxcXQxY0fks3CtyjoikvmSHR751lg9flzPbbz3sXTGAveZWRbRawW/\ncPdHBxTtIDMz6usqefSF7bR1dJKjKQ9FJIMlO+PVDURLOTuBzthiBxLW8N39RWDqAOML3MxJlTzw\n3BaWbtrPJWdUpDscEZHAJNsP/xNAnUeHV8gorz9rFNkRY8G63Ur4IpLRkq1hbAEOBhlIupTm53Bh\nzUjV8UUk4yV7hv8qsMDMfkf0DloA3P1bgUSVYjMnVXLXE2vZ1djK6NL8vjcQERmGkj3D3ww8CeQS\nveu265ER6usqAXQTlohktGR76XwJwMwK3T3j5gacMraUypI8Fq7bzTunj093OCIigUh2xqvXmdlq\nYv3ozex8M/vPQCNLITNj5qRKFr28h47OITHcj4jIoEu2pPMdYDawF8DdXwDeEFRQ6TBzUiUHW9pY\nvuVAukMREQlE0ncaufuWExZ1DHIsaXXZxFFETHV8EclcSXfLNLPXA25mOWb2aWBNgHGlXFlhLheM\nL2Ph2l3pDkVEJBDJJvyPEB0auRpoAC6Ivc4oMyeN5sWGg+xtPtL3yiIiw0yy4+Hvcfe/d/cqdx/t\n7u/NxLtu6+sqcYen1+9JdygiIoMu2bF0vhdn8UGiw3L+ZnBDSp9zq0dQXpTLgrW7ueYCzYIlIpkl\n2ZJOPtEyzsuxx3nAacAHzew7AcWWcpGI8YaJo3hq3W461T1TRDJMskMrnAfMcPcOADP7L2ARcCmw\nIqDY0mJmXSXzlm9j1bZGzj1tRLrDEREZNMme4Y8Ees6/VQSUx74AMuoK52UTo8MsLFBvHRHJMMkm\n/G8Ay83sx2b2E2AZcJeZFQF9TmY+nIwqzuO800aoP76IZJxke+n8EHg9MI/oZOSXuvu97n7I3W8K\nMsB0mDmpkqWb93PwcFu6QxERGTS9JnwzOzv290KiUxZuiT3GxJZlpPq6SjrVPVNEMkxfF21vBD4M\n3B3nPQcuH/SIhoDzTyujND+bhet2cdV5Y9MdjojIoOg14bv7h2N/Z6UmnKEhOyvCZRMrWbhuN+6O\nmaU7JBGRAeurpPOZHs/fecJ7dwQV1FAws66SnY1HeGlHU7pDEREZFH1dtL2ux/NbT3jvykGOZUiZ\nOUmzYIlIZukr4VuC5/FeZ5Sq0nwmjy1Vf3wRyRh9JXxP8Dze64wzc1Ilizfup/lIe7pDEREZsL4S\n/vlm1mhmTcB5seddr8/tbUMzG29m881stZmtMrOPD1rUKTJzUiXtnc4z6p4pIhmg14Tv7lnuXuru\nJe6eHXve9Tqnj323Aze6+xTgtcBHzWzKYAWeCtNqRlKcl606vohkhKSnOOwvd9/u7ktjz5uIzpA1\nrMYczs2O8PozK1i4Nto9U0RkOAss4fdkZrXAVODZVBxvMNXXjabhQAuv7G5OdygiIgNiQZ+5mlkx\nsBD4mrs/FOf964HrAaqqqqbNnTuX5uZmiouLT1w1Lfa2dHLjwhbec3Yus2v7qmINjqHU/lRT28PZ\ndgh3+wfS9lmzZi1x9+lJrezugT2AHOAJ4FPJrD9t2jR3d58/f74PJW+6e4G/996/pex4Q639qaS2\nh1eY2z+QthOdeTCpnBxYScei4xH8EFjj7t8K6jipMHNSJc9u2EfL0Y50hyIicsqCrOHPAN4HXG5m\ny2OPtwZ4vMDkZEU42t7J5C88zow7/8y8ZQ3pDklEpN+SneKw39z9aTLgbtx5yxr48TMbul83HGjh\n1oeiszrOmTqsOh2JSMilpJfOcHbXE2tpbe88bllLWwd3PbE2TRGJiJwaJfw+bDvQEnd5w4EWOjvV\nN19Ehg8l/D6MKytI+N6bvr2Q//3bJg4f1Vg7IjL0KeH34abZdRTkZB23LD87wvteW0NxXjafn7eS\n197xJ77+2BoaEvwaEBEZCgK7aJspui7M3vXEWrYdaGFcWQE3za5jztRq3J2lm/fzw6c38IOnXuXe\nRRu48pwx/NOM09m89xDf/MO6k7YREUkXJfwkzJlaHTdZmxnTasqZVlPO1v2Huf+vm3jguc387sXt\nmEHXTczq2SMiQ4FKOoPktJGF3PbWyfzt1jcyoiCnO9l3Uc8eEUk3JfxBVpSXTWNLW9z3VOMXkXRS\nwg9Aop49WWaaMlFE0kYJPwDxevbkZkWoKM7l/T9+nk89uJz9h46mKToRCSsl/ADMmVrN199+LtVl\nBRhQXVbAN95xHk99ZhYfu/wsHnlhG2/61kIeeWGbJlYRkZRRL52AJOrZc+Ob63jruWO55dcv8m8P\nLOM3yxr46rXnMHZE4hu8REQGg87w02Dy2FIe+tcZfO6qyTzzyh6u+NZT/PRvm3ho6VZm3PlnVjQc\n1KicIjLodIafJlkR40OXncEVU6q49aEVfG7eSiIGnQ6MV999ERl8OsNPs5qKIn72oUsoK8ihayy2\nRzdHPxb13ReRwaSEPwSYGQd79N3f3HxsGoGGAy08vGwrja3x+/aLiCRLJZ0hYlxZQfeNWR+Z3MG3\nV0Y/mojBJx98gZwsY8ZZo3jLOWO4YsoYyotymbesIe4YPyIi8SjhDxE3za7j1odW0NLWQSR2gl+Q\nk8XXrj2HmooiHl+5ncdW7uDmX6/gtodXcsaoIjbuOURbrA6kmr+I9EUJf4joOSonNFF9whn7tJqR\n3PbWyaza1shjK7fzPwtfpf2ECVha2jr4xuMvKeGLSFxK+ENIV9/9BQsW8LG/rz/pfTPjnOoRnFM9\ngv+c/0rcfWw72MqH71/MZRNHcelZozh9VBFmx64JqAwkEl5K+MNUz5p/T4W5WazZ3siTq3cC0bt8\nL5s4iksnjuLg4Ta++rs1tLR1ACoDiYSNEv4w1bPm36UgJ4s7rj2Xay4Yx6a9h1m0fg9Pv7yb363Y\nztznt8TdT1fXTyV8kcynhD9M9TYTF0DtqCJqRxXxvtfW0N7RyQtbD/J//usvcffVcKCFrzy6mslj\nS5k8toSJo0vIzT7WY1dlIJHMoIQ/jCUar+dE2VkRptWMpDpBGSgny/jp3zZxpL0zun7EOGt0MZPH\nltLR6Ty+cgdHO6LvJVMG0heEyNCkhB8iicpAX3/7uVx93lg27j3E6u1NrNneyJrtjfzllT3sbDxy\n0n5a2jr47LwVHGxpY3x5AeNHFnLayEIKcrOYt6zhuGPoOoHI0BFYwjezHwFXA7vc/ZygjiPJ66sM\ndNboEs4aXcLbzh/XvU3tLb+Lu69DRzr44iOrjltWWZLHwcNt3b8GurS0dXDnYy/xtvPHEYkY8XT9\nKrhufBOfvfPP+lUgEoAgz/B/AvwHcH+Ax5B+SrYM1CVRGWjciHzm3TCDLfta2Lr/MFv2HWbLvhYe\nXBz/4vCOxlbqPv8Yo0vyGV2aR1VJPlWleYwuzadhfwu/WrKVox2dtFXrV4FIUAJL+O7+lJnVBrV/\nSY1EZaDPXHl2NHmX5DOtZmT3e0+v3xP3C2JEQQ5/f8kEdjS2sqvxCK/sbuYvr+yhsbX9uPX+Y1V0\nprCWtg6+9NtVXDC+jJqKwuPuJeiiawUi/WNBzrgUS/iP9lbSMbPrgesBqqqqps2dO5fm5maKi4sD\ni2uoG2rtP9DSxs6DrRzt6CQ3K0LViHzKCnISrtuwv4XOHv9fRcyoHlkQd5sjHc5zmxppbofmNqOl\nw1h/ELa3GG2d0SRfnANnlGVx5ogIZ5ZlccaICEfb2/t1nOFgqH3uqRbm9g+k7bNmzVri7tOTWTft\nCb+n6dOn++LFi1mwYAH19fWBxTXUDff29/fMe8adf+7+VXDjue3cvSL6w3N0SR6fumISyzYfYNmW\n/by8qxl3MItOCH/i0BIQLUE9c8vlwTQsYMP9cx+oMLd/IG03s6QTvnrpyKDr73WCRGWj2946mTlT\nq7nu4gkANLa28eKWgyzbvJ+7n1wXd18NB1pYte0gZ48pJeuEC8QqAUnYKeFL2vU1cFyX0vwcLo0N\nEzH3+S1xrxUAXPW9pynNz+ai2nIuPj36eHX3IT43b6W6i0qoBdkt8wGgHhhlZluBL7r7D4M6ngxv\nfQ0cd6KEF5Nn11FWlMNzG/bx7IZ9/OmlXQAYcGIBSMNKSNgE2UvnPUHtW6SvewqunXoaALubjvDc\nhn189OdL4+6n4UALa3c0MamqOG5PIJFMopKODFvJXCuoLMnjqvPGcsfv499PADD7O08xqjiPS8+q\n4NKJlcw4q4KxIwoA1f0lsyjhSyj0VgIqysvm6fV7WPTyHuYt3wbAmZVFjCsr4NlX9/VrHCGRoUwJ\nX0KhrxLQuy4aT2en89KOJp5Zv4dF6/ewaN3uuHX/f9esYjJMKeFLaPRVAopEjCnjSpkyrpQPv+GM\nhOMIbT/YytX/bxEX11Z09wIqL8rtfl9lIBmqlPBFEkg0jlBJXjbFedn87NlN/OiZDQBMqirmktMr\nMINfPL+F1vb+DyetgeMkaEr4Igkkqvt/Zc45zJlazZH2DlZsPcizsS6gDy3dyqGjHSftp6UtOrJo\ndpZRWRwdMK6yJI+i3Cx+s3zbsWOM13UCCZYSvkgCfdX987KzmF5bzvTacj46C9o7Ojnrs4/F3dfB\nljZu+Pmy45YV5GRxtKOTjtgQEQu2R2cZ0/0BEhQlfJFe9GeYiOysSMIy0JjSfO77p4vZ1dTK7qYj\n3Y97n97Qvc7yvcfuA2g40MIvF29h1tmjGVWcN/CGiKCELzKoEpWBbnnL2dSNKaFuTMlx6z+2ckf3\nF8RHp3TwvVXRf5IRg5t+9SJmcMH4Mt40uYo3Th5NXVUJZqYLw3JKlPBFBlFfZaAT9fyCyInNG1+Q\nk8Ud157DxKoS/rhmJ39+aRd3PbGWu55YS3VZAaePKuK5Dbo/QPpPCV9kkPWnDNTXwHHnVI/gE2+a\nxM7GVua/tIs/rtnFn9bsjHt/wB2/X8NV540lJysS91j6VSBK+CJplszAcVWl+Vx38QSuu3hCwvsD\ndjUdYcoXHues0SVMHlPC5LGlnD22hLPHlPLM+j2aXF6U8EWGm0QXhkcW5vDuiybw0o5GnnllDw8t\na+h+L2Jw4nwx6g0UPkr4IsNMogvDX/y71xyXvPcdOspLOxpZs72Jrzy6Ou6+Gg600NrWQX5OVuBx\nS/op4YsMM8leGC4vyuX1Z47i9WeO4kdPb0g4Wui0rzzJrLNHc+U5Y5hVN5qivGhaUM0/8yjhiwxD\ngzGNZH52hPfPqOVgSztPrt7Boy9uJy87whsmVVJZnMdDS7f2a4gIGfqU8EVCoK9fBV+dcw6LN+7j\nsZU7eGLVDrYfbD1pH6r5D39K+CIh0duvgqyIcckZFVxyRgVfuHoKZ9z2+7jrNRxo4QM/fo5JY0qo\nq4reSHZmZXH3NQCVgYY2JXwROU4kYgl7AhXkZLH9YCtPr99DW0e020/EoHZUESV52aza1kh7rDtQ\nw4EWbvn1i0DvZSCNFpo6SvgicpJEPYG+/vZzmTO1mraOTjbuOcTanU2s29HE2p1N/HHNru6B4Lq0\ntnfyyQeX880/rKWiKJfyolxGFuXGnuexae8hHlrawNGOTlrH6VpB0JTwReQkfdX8c7IiTKwqYWJV\nCZwX3eb0BDeEOXBRbTl7Dx1lT/NR1u1sZu+hI7S2dR633vdXR9NRS1sHn//NSkrys3nNuBFUleYd\nN8G8ykanTglfROLqb0+gcQnKQNVlBXz73RectLzlaAeTv/B49+vLxnSwaEf0WkBTazsfvG8xABVF\nubymegSvGVfK4SPtzH1+C0fUe+iUxB90Q0Skn26aXUfBCTdwFeRkcdPsurjrF+RmUV1W0P364spj\n5aCxI/L55Udex+1/N4XLzx7NnqYj3LvoVe7766buZN+lpa2Dbzz+0iC2JHPpDF9EBkV/RwqFxNcK\nbr7ybC6qLeei2vLu5UfaO6jNT0/PAAAHnUlEQVT73OPxdsO2g628+3/+yuvOrOB1Z1RwwYQy8rKP\nffmcShmov9sM5BipumAdaMI3syuB7wJZwL3ufmeQxxOR9OpvGaiv0UJ7ysvOSth7qDgvm0NH2/nu\nn17mO398mbzsCNNqRvK6Mypo6+jknqde7fc8w/0ZbK6/65+0TYqmtwws4ZtZFvB94ApgK/C8mT3i\n7vEH9RCRUEpmtNAuiX4RfDU2z/DBw208t3Eff31lL399dS93P7ku7n5a2jr43LyVrN7eiAFmhlm0\ni6lh3PeXjccdo2ubL/xmZfcXTiS2jQHfn78+7vq3P7KKlrYOssyIRIysSHS7rIjx5d+uPvYFcejY\nNkHe3BbkGf7FwHp3fxXAzOYC1wBK+CJySvoqG40ozOGKKVVcMaUKgP2HjjL1K0/G3VfzkXb+96+b\n6HTHAXfHHTrdTxpZtEtja3vs10hyDrS0dZ+192b+9mPlp20JxjwaDOaeoGUD3bHZO4Ar3f1Dsdfv\nAy5x9xtOWO964HqAqqqqaXPnzqW5uZni4uJA4hoOwtx+tT2cbYfg2r92R1P37GA95WZFTppy8sRt\n3OmebMaB3EiEiWOK6UqbXe+v29ncfQzvsTw7EuH0yiI6neMe7vDqnkMciR2jqgC6ep72Flc8s2bN\nWuLu05NZN+0Xbd39HuAegOnTp3t9fT0LFiygvr4+vYGlUZjbr7bXpzuMtAmq/QdOqK/DsZvI6hOU\nTnrb5oo427T2sn4yNfwbz23n7hXZfcY1UEEm/AZgfI/Xp8WWiYikzKn0HurvNgM9Rl8XrAdLkAn/\neWCimZ1ONNFfB/zfAI8nIhJXf3sPnco2AzlGMhesB0NgCd/d283sBuAJot0yf+Tuq4I6noiI9C7Q\nGr67/x6IP86qiIiklIZWEBEJCSV8EZGQUMIXEQkJJXwRkZAI7E7bU2Fmu4FNwChgT5rDSacwt19t\nD68wt38gba9x98pkVhxSCb+LmS1O9lbhTBTm9qvt4Ww7hLv9qWq7SjoiIiGhhC8iEhJDNeHfk+4A\n0izM7VfbwyvM7U9J24dkDV9ERAbfUD3DFxGRQaaELyISEkMu4ZvZlWa21szWm9kt6Y4nlcxso5mt\nMLPlZrY43fEEzcx+ZGa7zGxlj2XlZvakmb0c+zsynTEGJUHbbzezhtjnv9zM3prOGINiZuPNbL6Z\nrTazVWb28djyjP/se2l7Sj77IVXDj018vo4eE58D7wnLxOdmthGY7u6huPnEzN4ANAP3u/s5sWXf\nAPa5+52xL/yR7n5zOuMMQoK23w40u/s30xlb0MxsLDDW3ZeaWQmwBJgDvJ8M/+x7afu7SMFnP9TO\n8LsnPnf3o0DXxOeSgdz9KWDfCYuvAe6LPb+P6D+GjJOg7aHg7tvdfWnseROwBqgmBJ99L21PiaGW\n8KuBLT1ebyWF/zGGAAf+YGZLYpO7h1GVu2+PPd8BVKUzmDS4wcxejJV8Mq6kcSIzqwWmAs8Sss/+\nhLZDCj77oZbww+5Sd78QeAvw0djP/tDyaL1x6NQcg/dfwJnABcB24O70hhMsMysGfg18wt0be76X\n6Z99nLan5LMfagk/1BOfu3tD7O8u4GGiJa6w2Rmrc3bVO3elOZ6Ucfed7t7h7p3AD8jgz9/Mcogm\nvJ+5+0OxxaH47OO1PVWf/VBL+N0Tn5tZLtGJzx9Jc0wpYWZFsYs4mFkR8GZgZe9bZaRHgH+MPf9H\n4DdpjCWlupJdzLVk6OdvZgb8EFjj7t/q8VbGf/aJ2p6qz35I9dIBiHVH+g7HJj7/WppDSgkzO4Po\nWT1E5xr+eaa33cweAOqJDg27E/giMA/4BTCB6FDZ73L3jLu4maDt9UR/0juwEfjnHjXtjGFmlwKL\ngBVAZ2zxbURr2Rn92ffS9veQgs9+yCV8EREJxlAr6YiISECU8EVEQkIJX0QkJJTwRURCQglfRCQk\nlPAl45lZR49RCJfHbmnv7z7eb2bjBj86kdTJTncAIinQ4u4XDHAf7yd6M8y2ZDcwsyx37xjgcUUG\njc7wJZTMrNbMFpnZ0tjj9T3euzk2L8ELZnanmb0DmA78LPYLocDM3mhmy2Lr/cjM8mLbbjSzfzez\npcA7zezfYmOfv2hmc9PUXBFAZ/gSDgVmtjz2fIO7X0t0nJYr3L3VzCYCDwDTzewtRIfpvcTdD5tZ\nubvvM7MbgE+7+2Izywd+ArzR3deZ2f3AvxC9Qxxgb2wQPMxsG3C6ux8xs7KUtVgkDp3hSxi0uPsF\nsce1sWU5wA/MbAXwS2BKbPmbgB+7+2GABLf21xH94lgXe30f0HNk0wd7PH+R6C+D9wLtg9MckVOj\nhC9h9UmiY9icT7RckzuI+z7U4/lVwPeBC4HnzUy/qiVtlPAlrEYA22PD0b6P6GB9AE8CHzCzQojO\nsxpb3gSUxJ6vBWrN7KzY6/cBC088gJlFgPHuPh+4OXbM4gDaIpIUJXwJq/8E/tHMXgDOJnZW7u6P\nEx2md3Gs7v/p2Po/Af47tsyADwC/jJWEOoH/jnOMLOCnsXWWAd9z9wPBNUmkdxotU0QkJHSGLyIS\nEkr4IiIhoYQvIhISSvgiIiGhhC8iEhJK+CIiIaGELyISEv8fxty8PUlfpD4AAAAASUVORK5CYII=\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": []
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "9AmVXMkJXTdi",
        "colab_type": "text"
      },
      "source": [
        "The scree plot method draws a straight line for each factor and its eigenvalues. Number eigenvalues greater than one considered as the number of factors.\n",
        "\n",
        "Here, you can see again 10-factors eigenvalues are greater than one. It means we need to choose only 10 factors (or unobserved variables).\n",
        "\n",
        "However something we can additionally, see is that after the first two, there really isnt a diffence in the eigenvalues. So really only the first two are significant. Therefore we will only choose **two factors**."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "rP0mIgpSyt6j",
        "colab_type": "text"
      },
      "source": [
        "## Running factor analysis\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "es4vNDC1XYdW",
        "colab_type": "code",
        "outputId": "5dd2ef4c-8ce1-40e4-802f-3bc46651a37f",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 824
        }
      },
      "source": [
        "# Create factor analysis object and perform factor analysis\n",
        "fa = FactorAnalyzer()\n",
        "numFactors = 2 # this will vary based on what you saw above\n",
        "fa.analyze(df, numFactors, rotation=\"varimax\")\n",
        "fa.loadings"
      ],
      "execution_count": 31,
      "outputs": [
        {
          "output_type": "execute_result",
          "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>Factor1</th>\n",
              "      <th>Factor2</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>water pH 1</th>\n",
              "      <td>0.119985</td>\n",
              "      <td>-0.021944</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>water pH 2</th>\n",
              "      <td>0.044566</td>\n",
              "      <td>0.000270</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>ground water level 1</th>\n",
              "      <td>0.799025</td>\n",
              "      <td>-0.094661</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>ground water level 2</th>\n",
              "      <td>-0.044711</td>\n",
              "      <td>0.022975</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>soil moisture 1</th>\n",
              "      <td>0.975765</td>\n",
              "      <td>-0.099547</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>soil moisture 2</th>\n",
              "      <td>0.987725</td>\n",
              "      <td>-0.085610</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>soil moisture 3</th>\n",
              "      <td>0.993876</td>\n",
              "      <td>-0.086575</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>motor monitor 1</th>\n",
              "      <td>0.155475</td>\n",
              "      <td>0.133223</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>motor monitor 2</th>\n",
              "      <td>0.052529</td>\n",
              "      <td>-0.083869</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>salinity I</th>\n",
              "      <td>-0.041657</td>\n",
              "      <td>0.030068</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>salinity E</th>\n",
              "      <td>0.101357</td>\n",
              "      <td>-0.130272</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>air RH I</th>\n",
              "      <td>0.103788</td>\n",
              "      <td>-0.033804</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>air RH E</th>\n",
              "      <td>-0.041323</td>\n",
              "      <td>0.001375</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>air temp I</th>\n",
              "      <td>0.033716</td>\n",
              "      <td>-0.124489</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>air temp E</th>\n",
              "      <td>-0.142752</td>\n",
              "      <td>0.004564</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>wall RH 1 I</th>\n",
              "      <td>0.887812</td>\n",
              "      <td>0.014285</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>wall RH 2 I</th>\n",
              "      <td>0.863680</td>\n",
              "      <td>-0.062129</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>wall RH 3 I</th>\n",
              "      <td>-0.028089</td>\n",
              "      <td>0.107985</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>wall RH 1 E</th>\n",
              "      <td>0.507402</td>\n",
              "      <td>0.828640</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>wall RH 2 E</th>\n",
              "      <td>0.515331</td>\n",
              "      <td>0.827161</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>wall RH 3 E</th>\n",
              "      <td>-0.024609</td>\n",
              "      <td>0.967112</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>rain intensity</th>\n",
              "      <td>-0.005517</td>\n",
              "      <td>0.449248</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>rain duration</th>\n",
              "      <td>-0.012829</td>\n",
              "      <td>0.475564</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>wind speed</th>\n",
              "      <td>-0.064279</td>\n",
              "      <td>0.141067</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>wind direction (radians)</th>\n",
              "      <td>0.247983</td>\n",
              "      <td>-0.257226</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "                           Factor1   Factor2\n",
              "water pH 1                0.119985 -0.021944\n",
              "water pH 2                0.044566  0.000270\n",
              "ground water level 1      0.799025 -0.094661\n",
              "ground water level 2     -0.044711  0.022975\n",
              "soil moisture 1           0.975765 -0.099547\n",
              "soil moisture 2           0.987725 -0.085610\n",
              "soil moisture 3           0.993876 -0.086575\n",
              "motor monitor 1           0.155475  0.133223\n",
              "motor monitor 2           0.052529 -0.083869\n",
              "salinity I               -0.041657  0.030068\n",
              "salinity E                0.101357 -0.130272\n",
              "air RH I                  0.103788 -0.033804\n",
              "air RH E                 -0.041323  0.001375\n",
              "air temp I                0.033716 -0.124489\n",
              "air temp E               -0.142752  0.004564\n",
              "wall RH 1 I               0.887812  0.014285\n",
              "wall RH 2 I               0.863680 -0.062129\n",
              "wall RH 3 I              -0.028089  0.107985\n",
              "wall RH 1 E               0.507402  0.828640\n",
              "wall RH 2 E               0.515331  0.827161\n",
              "wall RH 3 E              -0.024609  0.967112\n",
              "rain intensity           -0.005517  0.449248\n",
              "rain duration            -0.012829  0.475564\n",
              "wind speed               -0.064279  0.141067\n",
              "wind direction (radians)  0.247983 -0.257226"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 31
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "y1ylue-xy25x",
        "colab_type": "text"
      },
      "source": [
        "## Translating this into something we can understand a bit better"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "jXtE5F5VaLvW",
        "colab_type": "code",
        "outputId": "48e97d98-1567-4202-907a-be2f5b4ecbc7",
        "cellView": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 71
        }
      },
      "source": [
        "\n",
        "L = np.array(fa.loadings)\n",
        "headings = list( fa.loadings.transpose().keys() )\n",
        "factor_threshold = 0.25\n",
        "for i, factor in enumerate(L.transpose()):\n",
        "  descending = np.argsort(np.abs(factor))[::-1]\n",
        "  contributions = [(np.round(factor[x],2),headings[x]) for x in descending if np.abs(factor[x])>factor_threshold]\n",
        "  print('Factor %d:'%(i+1),contributions)"
      ],
      "execution_count": 33,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Factor 1: [(0.99, 'soil moisture 3'), (0.99, 'soil moisture 2'), (0.98, 'soil moisture 1'), (0.89, 'wall RH 1 I'), (0.86, 'wall RH 2 I'), (0.8, 'ground water level 1'), (0.52, 'wall RH 2 E'), (0.51, 'wall RH 1 E')]\n",
            "Factor 2: [(0.97, 'wall RH 3 E'), (0.83, 'wall RH 1 E'), (0.83, 'wall RH 2 E'), (0.48, 'rain duration'), (0.45, 'rain intensity'), (-0.26, 'wind direction (radians)')]\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "4LMjFfoJxaDZ",
        "colab_type": "text"
      },
      "source": [
        "**So what does this mean?**"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "U2ubfqwlXkym",
        "colab_type": "text"
      },
      "source": [
        "## Solution\n",
        "\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "hMuZ2xQJxXdO",
        "colab_type": "text"
      },
      "source": [
        "Factor 1 has high factor loadings for soil moisture, RH interior and exterior up to a certain height, as well as ground water only for the first sensor.\n",
        "\n",
        "\n",
        "\n",
        "Factor 2 has high factor loadings for only external RH values all the way up and down the wall as well as rain duration, intensity, and wind direciton."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "RgwoZseOxQ_n",
        "colab_type": "text"
      },
      "source": [
        "# Computing the variance"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "wGC0sDVcXthY",
        "colab_type": "code",
        "cellView": "both",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 142
        },
        "outputId": "45cc89de-6194-4f3f-c5fb-77654f9f1411"
      },
      "source": [
        "# Get variance of each factors\n",
        "fa.get_factor_variance()"
      ],
      "execution_count": 34,
      "outputs": [
        {
          "output_type": "execute_result",
          "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>Factor1</th>\n",
              "      <th>Factor2</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>SS Loadings</th>\n",
              "      <td>5.769659</td>\n",
              "      <td>2.929944</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Proportion Var</th>\n",
              "      <td>0.230786</td>\n",
              "      <td>0.117198</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Cumulative Var</th>\n",
              "      <td>0.230786</td>\n",
              "      <td>0.347984</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "                 Factor1   Factor2\n",
              "SS Loadings     5.769659  2.929944\n",
              "Proportion Var  0.230786  0.117198\n",
              "Cumulative Var  0.230786  0.347984"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 34
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "0Qn7XAVxXyFS",
        "colab_type": "text"
      },
      "source": [
        "Total 35% cumulative Variance explained by the 2 factors."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "AenH31yxX3E5",
        "colab_type": "text"
      },
      "source": [
        "**Pros and Cons of Factor Analysis**\n",
        "\n",
        "Factor analysis explores large dataset and finds interlinked associations. It reduces the observed variables into a few unobserved variables or identifies the groups of inter-related variables, which help the market researchers to compress the market situations and find the hidden relationship among consumer taste, preference, and cultural influence. Also, It helps in improve questionnaire in for future surveys. Factors make for more natural data interpretation.\n",
        "\n",
        "Results of factor analysis are controversial. Its interpretations can be debatable because more than one interpretation can be made of the same data factors. After factor identification and naming of factors requires domain knowledge."
      ]
    }
  ]
}