{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# A Statistical Analysis of the Virginia Smart Scale Scoring Methodology\n",
    "\n",
    "#### September 17, 2018\n",
    "\n",
    "#### Christopher Garcia, Ph.D., Mukesh Srivastava, D.B.A.\n",
    "#### University of Mary Washington, Center for Business Research\n",
    "\n",
    "## Executive Summary\n",
    "\n",
    "The Smart Scale methodology is intended to provide an objective, transparent, and consistent way to evaluate and compare the benefits of competing transportation projects in Virginia. This study provides a statistical analysis of a set of 404 Smart Scale Round 2 projects scored using the Smart Scale methodology. The study's first goal was to provide a basic characterization of the Smart Scale score distribution across different project features, while the second goal was to answer the following two questions:\n",
    "\n",
    "1. Was the Smart Scale methodology applied consistently across all projects?\n",
    "\n",
    "2. If Accessibility (a key concern for the Fredericksburg area) was prioritized in Smart Scale score calculation for area types A and D (the area types where Fredericksburg projects occur), how would this impact the competitiveness of Fredericksburg-based projects?\n",
    "\n",
    "To answer the first question, the Smart Scale scores were independently calculated for each project using the methodology prescribed in the November 2017 Smart Scale Technical Guide. These scores were then compared to the project scores assigned by the state. There was no discrepancy between the calculated scores and scores assigned by the state, indicating that the Smart Scale methodology was consistently applied accross all projects. \n",
    "\n",
    "To answer the second question, two alternative scenarios were used which both increased the importance of Accessibility in the  scoring calculation. All projects were scored according to the methodology defined in each scenario, and the percentile rankings of the Fredericksburg scores under each alternative scenario were compared to that obtained by the current methodology. There were small increases in the percentile rankings of some Fredericksburg project scores under one of the alternative scenarios when compared to the current scoring methodology. This suggests that the full impact of some Fredericksburg-based projects may be underestimated under the current Smart Scale methodology and, depending on the cut score, these projects may be edged out of funding accordingly.\n",
    "\n",
    "\n",
    "## I. Introduction\n",
    "\n",
    "The Virginia Smart Scale methodology was finalized in November 2017 in order to provide an objective scoring system for the state's transportation projects competing for funding. The complete technical specification of the Smart Scale methodology can be found in the November 2017 technical guide, found here: [http://vasmartscale.org/documents/20171115/ss_technical_guide_nov13_2017.pdf](http://vasmartscale.org/documents/20171115/ss_technical_guide_nov13_2017.pdf).\n",
    "\n",
    "The purpose of this analysis is 1) to provide a descriptive statistical characterization of recently funded projects, 2) to independently assess whether the Smart Scale scores (abbreviated by SSS hereafter) assigned to these projects were consistent with the Smart Scale methodology as specified in the technical guide, and 3) to determine how increasing the emphasis on Accessibility (which is a major concern to the Fredericksburg area) in the SSS calculation methodology would impact the competitiveness of Fredericksburg-based projects. This analysis uses a dataset provided by the Fredericksburg Area Metropolitan Planning Organization (FAMPO), containing data on 404 Smart Scale Round 2 projects. This analysis also uses the R statistical computing environment within a Jupyter notebook. This allows the analysis to be both explained as well as executed. It furthermore enables complete transparancy of the analytical methodology and code used in its execution. This removes any ambiguity and also guarantees reproducibility of the results.\n",
    "\n",
    "## II. Data Inspection and Preparation\n",
    "\n",
    "The data was provided in a Microsoft Excel spreadsheet. This spreadsheet was first converted to a CSV file (ss-data.csv) and then a number of transformation processes were applied to clean and scale the data as well as to impute missing values. To begin, we first import several libraries which will be used later on and also define several functions which will be used in the data preparation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Turn off verbose output messages.\n",
    "options(warn=0)\n",
    "\n",
    "# Include all libraries that will be used in this notebook's code\n",
    "library(ggplot2)\n",
    "library(gridExtra)\n",
    "\n",
    "# Given a vector, a vector of original values, and a corresponding vector of \n",
    "# replacement values, apply the replacements and return the transformed vector.\n",
    "recode <- function(vec, from.vals, to.vals) {\n",
    "    f = function(v) {\n",
    "        for(i in 1:length(from.vals)) {\n",
    "            if(v == from.vals[i]) { return(to.vals[i]) }\n",
    "        }\n",
    "        return(v)\n",
    "    }\n",
    "    return(sapply(as.vector(vec), f))\n",
    "}\n",
    "\n",
    "# Function for imputing mean of a vector to its missing values.\n",
    "impute.mean <- function(vec) {\n",
    "    vec[is.na(vec)] = mean(vec, na.rm=TRUE)\n",
    "    vec\n",
    "}\n",
    "\n",
    "# Function for imputing value of 0 to each missing value in the vector.\n",
    "impute.0 <- function(vec) {\n",
    "    vec[is.na(vec)] = 0\n",
    "    vec\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next we read in the dataset and inspect the first few rows of data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>App.Id</th><th scope=col>Area.Type</th><th scope=col>District</th><th scope=col>Organization.Name</th><th scope=col>Project.Title</th><th scope=col>Statewide.High.Priority</th><th scope=col>District.Grant</th><th scope=col>Throughput.Score</th><th scope=col>Delay.Score</th><th scope=col>Crash.Frequency.Score</th><th scope=col>...</th><th scope=col>Travel.Time.Reliability.Score</th><th scope=col>Land.Use.Score</th><th scope=col>Project.Benefit.Score</th><th scope=col>Total.Project</th><th scope=col>Score.Divided.by.Total.Cost</th><th scope=col>SMART.SCALE.Request</th><th scope=col>SMART.SCALE.Score</th><th scope=col>Benefit.Rank</th><th scope=col>State.Rank</th><th scope=col>District.Rank</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>1414                                                        </td><td>A                                                           </td><td>NOVA                                                        </td><td>Northern Virginia Transportation Commission                 </td><td>VRE Fredericksburg Line Capacity Expansion                  </td><td>x                                                           </td><td>                                                            </td><td>17.69                                                       </td><td>87.66                                                       </td><td>100                                                         </td><td>...                                                         </td><td> -                                                          </td><td>64.016                                                      </td><td>64.236                                                      </td><td>216034920                                                   </td><td>2.973                                                       </td><td> 92636120                                                   </td><td>6.934                                                       </td><td>1                                                           </td><td>109                                                         </td><td>16                                                          </td></tr>\n",
       "\t<tr><td>1057                                                        </td><td>A                                                           </td><td>Hampton Roads                                               </td><td>Hampton Roads Transportation Planning Organization          </td><td>I-64 Southside Widening and High Rise Bridge - Phase 1      </td><td>x                                                           </td><td>                                                            </td><td>100                                                         </td><td>88.96                                                       </td><td>45.16                                                       </td><td>...                                                         </td><td>29.044                                                      </td><td>10.464                                                      </td><td>62.042                                                      </td><td>600000000                                                   </td><td>1.034                                                       </td><td>100000000                                                   </td><td>6.204                                                       </td><td>2                                                           </td><td>120                                                         </td><td>22                                                          </td></tr>\n",
       "\t<tr><td>1090                                                        </td><td>A                                                           </td><td>Hampton Roads                                               </td><td>Hampton Roads Transportation Planning Organization          </td><td>I-64/I-264 Interchange Improvements                         </td><td>x                                                           </td><td>                                                            </td><td>60.61                                                       </td><td>53.86                                                       </td><td> -                                                          </td><td>...                                                         </td><td>24.245                                                      </td><td>31.671                                                      </td><td>48.747                                                      </td><td>350091800                                                   </td><td>1.392                                                       </td><td> 50000000                                                   </td><td>9.749                                                       </td><td>3                                                           </td><td> 79                                                         </td><td>18                                                          </td></tr>\n",
       "\t<tr><td>1293                                                        </td><td>A                                                           </td><td>NOVA                                                        </td><td>Prince William County                                       </td><td>Route 234 At Balls Ford Intrchng and Rel/Widen Balls Ford Rd</td><td>x                                                           </td><td>x                                                           </td><td>65.48                                                       </td><td>100                                                         </td><td>15.61                                                       </td><td>...                                                         </td><td>1.05                                                        </td><td> -                                                          </td><td>41.289                                                      </td><td>126027000                                                   </td><td>3.276                                                       </td><td>124027000                                                   </td><td>3.329                                                       </td><td>4                                                           </td><td>170                                                         </td><td>30                                                          </td></tr>\n",
       "\t<tr><td>1249                                                        </td><td>A                                                           </td><td>NOVA                                                        </td><td>Fairfax County                                              </td><td>VA 286 - Popes Head Road Interchange                        </td><td>x                                                           </td><td>x                                                           </td><td>61.09                                                       </td><td>54.17                                                       </td><td>8.22                                                        </td><td>...                                                         </td><td>0.118                                                       </td><td> -                                                          </td><td>37.194                                                      </td><td> 64303070                                                   </td><td>5.784                                                       </td><td> 50558370                                                   </td><td>7.357                                                       </td><td>5                                                           </td><td>101                                                         </td><td>14                                                          </td></tr>\n",
       "\t<tr><td>1240                                                        </td><td>A                                                           </td><td>NOVA                                                        </td><td>Loudoun County                                              </td><td>Loudoun County Parkway (Shellhorn Road to US Route 50)      </td><td>x                                                           </td><td>x                                                           </td><td>45.27                                                       </td><td>69.78                                                       </td><td>5.43                                                        </td><td>...                                                         </td><td>1.6                                                         </td><td>16.364                                                      </td><td>33.732                                                      </td><td>112053000                                                   </td><td>3.01                                                        </td><td>112053000                                                   </td><td>3.01                                                        </td><td>6                                                           </td><td>179                                                         </td><td>32                                                          </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|llllllllllllllllllllllllllll}\n",
       " App.Id & Area.Type & District & Organization.Name & Project.Title & Statewide.High.Priority & District.Grant & Throughput.Score & Delay.Score & Crash.Frequency.Score & ... & Travel.Time.Reliability.Score & Land.Use.Score & Project.Benefit.Score & Total.Project & Score.Divided.by.Total.Cost & SMART.SCALE.Request & SMART.SCALE.Score & Benefit.Rank & State.Rank & District.Rank\\\\\n",
       "\\hline\n",
       "\t 1414                                                         & A                                                            & NOVA                                                         & Northern Virginia Transportation Commission                  & VRE Fredericksburg Line Capacity Expansion                   & x                                                            &                                                              & 17.69                                                        & 87.66                                                        & 100                                                          & ...                                                          &  -                                                           & 64.016                                                       & 64.236                                                       & 216034920                                                    & 2.973                                                        &  92636120                                                    & 6.934                                                        & 1                                                            & 109                                                          & 16                                                          \\\\\n",
       "\t 1057                                                         & A                                                            & Hampton Roads                                                & Hampton Roads Transportation Planning Organization           & I-64 Southside Widening and High Rise Bridge - Phase 1       & x                                                            &                                                              & 100                                                          & 88.96                                                        & 45.16                                                        & ...                                                          & 29.044                                                       & 10.464                                                       & 62.042                                                       & 600000000                                                    & 1.034                                                        & 100000000                                                    & 6.204                                                        & 2                                                            & 120                                                          & 22                                                          \\\\\n",
       "\t 1090                                                         & A                                                            & Hampton Roads                                                & Hampton Roads Transportation Planning Organization           & I-64/I-264 Interchange Improvements                          & x                                                            &                                                              & 60.61                                                        & 53.86                                                        &  -                                                           & ...                                                          & 24.245                                                       & 31.671                                                       & 48.747                                                       & 350091800                                                    & 1.392                                                        &  50000000                                                    & 9.749                                                        & 3                                                            &  79                                                          & 18                                                          \\\\\n",
       "\t 1293                                                         & A                                                            & NOVA                                                         & Prince William County                                        & Route 234 At Balls Ford Intrchng and Rel/Widen Balls Ford Rd & x                                                            & x                                                            & 65.48                                                        & 100                                                          & 15.61                                                        & ...                                                          & 1.05                                                         &  -                                                           & 41.289                                                       & 126027000                                                    & 3.276                                                        & 124027000                                                    & 3.329                                                        & 4                                                            & 170                                                          & 30                                                          \\\\\n",
       "\t 1249                                                         & A                                                            & NOVA                                                         & Fairfax County                                               & VA 286 - Popes Head Road Interchange                         & x                                                            & x                                                            & 61.09                                                        & 54.17                                                        & 8.22                                                         & ...                                                          & 0.118                                                        &  -                                                           & 37.194                                                       &  64303070                                                    & 5.784                                                        &  50558370                                                    & 7.357                                                        & 5                                                            & 101                                                          & 14                                                          \\\\\n",
       "\t 1240                                                         & A                                                            & NOVA                                                         & Loudoun County                                               & Loudoun County Parkway (Shellhorn Road to US Route 50)       & x                                                            & x                                                            & 45.27                                                        & 69.78                                                        & 5.43                                                         & ...                                                          & 1.6                                                          & 16.364                                                       & 33.732                                                       & 112053000                                                    & 3.01                                                         & 112053000                                                    & 3.01                                                         & 6                                                            & 179                                                          & 32                                                          \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "App.Id | Area.Type | District | Organization.Name | Project.Title | Statewide.High.Priority | District.Grant | Throughput.Score | Delay.Score | Crash.Frequency.Score | ... | Travel.Time.Reliability.Score | Land.Use.Score | Project.Benefit.Score | Total.Project | Score.Divided.by.Total.Cost | SMART.SCALE.Request | SMART.SCALE.Score | Benefit.Rank | State.Rank | District.Rank | \n",
       "|---|---|---|---|---|---|\n",
       "| 1414                                                         | A                                                            | NOVA                                                         | Northern Virginia Transportation Commission                  | VRE Fredericksburg Line Capacity Expansion                   | x                                                            |                                                              | 17.69                                                        | 87.66                                                        | 100                                                          | ...                                                          |  -                                                           | 64.016                                                       | 64.236                                                       | 216034920                                                    | 2.973                                                        |  92636120                                                    | 6.934                                                        | 1                                                            | 109                                                          | 16                                                           | \n",
       "| 1057                                                         | A                                                            | Hampton Roads                                                | Hampton Roads Transportation Planning Organization           | I-64 Southside Widening and High Rise Bridge - Phase 1       | x                                                            |                                                              | 100                                                          | 88.96                                                        | 45.16                                                        | ...                                                          | 29.044                                                       | 10.464                                                       | 62.042                                                       | 600000000                                                    | 1.034                                                        | 100000000                                                    | 6.204                                                        | 2                                                            | 120                                                          | 22                                                           | \n",
       "| 1090                                                         | A                                                            | Hampton Roads                                                | Hampton Roads Transportation Planning Organization           | I-64/I-264 Interchange Improvements                          | x                                                            |                                                              | 60.61                                                        | 53.86                                                        |  -                                                           | ...                                                          | 24.245                                                       | 31.671                                                       | 48.747                                                       | 350091800                                                    | 1.392                                                        |  50000000                                                    | 9.749                                                        | 3                                                            |  79                                                          | 18                                                           | \n",
       "| 1293                                                         | A                                                            | NOVA                                                         | Prince William County                                        | Route 234 At Balls Ford Intrchng and Rel/Widen Balls Ford Rd | x                                                            | x                                                            | 65.48                                                        | 100                                                          | 15.61                                                        | ...                                                          | 1.05                                                         |  -                                                           | 41.289                                                       | 126027000                                                    | 3.276                                                        | 124027000                                                    | 3.329                                                        | 4                                                            | 170                                                          | 30                                                           | \n",
       "| 1249                                                         | A                                                            | NOVA                                                         | Fairfax County                                               | VA 286 - Popes Head Road Interchange                         | x                                                            | x                                                            | 61.09                                                        | 54.17                                                        | 8.22                                                         | ...                                                          | 0.118                                                        |  -                                                           | 37.194                                                       |  64303070                                                    | 5.784                                                        |  50558370                                                    | 7.357                                                        | 5                                                            | 101                                                          | 14                                                           | \n",
       "| 1240                                                         | A                                                            | NOVA                                                         | Loudoun County                                               | Loudoun County Parkway (Shellhorn Road to US Route 50)       | x                                                            | x                                                            | 45.27                                                        | 69.78                                                        | 5.43                                                         | ...                                                          | 1.6                                                          | 16.364                                                       | 33.732                                                       | 112053000                                                    | 3.01                                                         | 112053000                                                    | 3.01                                                         | 6                                                            | 179                                                          | 32                                                           | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "  App.Id Area.Type District     \n",
       "1 1414   A         NOVA         \n",
       "2 1057   A         Hampton Roads\n",
       "3 1090   A         Hampton Roads\n",
       "4 1293   A         NOVA         \n",
       "5 1249   A         NOVA         \n",
       "6 1240   A         NOVA         \n",
       "  Organization.Name                                 \n",
       "1 Northern Virginia Transportation Commission       \n",
       "2 Hampton Roads Transportation Planning Organization\n",
       "3 Hampton Roads Transportation Planning Organization\n",
       "4 Prince William County                             \n",
       "5 Fairfax County                                    \n",
       "6 Loudoun County                                    \n",
       "  Project.Title                                               \n",
       "1 VRE Fredericksburg Line Capacity Expansion                  \n",
       "2 I-64 Southside Widening and High Rise Bridge - Phase 1      \n",
       "3 I-64/I-264 Interchange Improvements                         \n",
       "4 Route 234 At Balls Ford Intrchng and Rel/Widen Balls Ford Rd\n",
       "5 VA 286 - Popes Head Road Interchange                        \n",
       "6 Loudoun County Parkway (Shellhorn Road to US Route 50)      \n",
       "  Statewide.High.Priority District.Grant Throughput.Score Delay.Score\n",
       "1 x                                      17.69            87.66      \n",
       "2 x                                      100              88.96      \n",
       "3 x                                      60.61            53.86      \n",
       "4 x                       x              65.48            100        \n",
       "5 x                       x              61.09            54.17      \n",
       "6 x                       x              45.27            69.78      \n",
       "  Crash.Frequency.Score ... Travel.Time.Reliability.Score Land.Use.Score\n",
       "1 100                   ...  -                            64.016        \n",
       "2 45.16                 ... 29.044                        10.464        \n",
       "3  -                    ... 24.245                        31.671        \n",
       "4 15.61                 ... 1.05                           -            \n",
       "5 8.22                  ... 0.118                          -            \n",
       "6 5.43                  ... 1.6                           16.364        \n",
       "  Project.Benefit.Score Total.Project Score.Divided.by.Total.Cost\n",
       "1 64.236                216034920     2.973                      \n",
       "2 62.042                600000000     1.034                      \n",
       "3 48.747                350091800     1.392                      \n",
       "4 41.289                126027000     3.276                      \n",
       "5 37.194                 64303070     5.784                      \n",
       "6 33.732                112053000     3.01                       \n",
       "  SMART.SCALE.Request SMART.SCALE.Score Benefit.Rank State.Rank District.Rank\n",
       "1  92636120           6.934             1            109        16           \n",
       "2 100000000           6.204             2            120        22           \n",
       "3  50000000           9.749             3             79        18           \n",
       "4 124027000           3.329             4            170        30           \n",
       "5  50558370           7.357             5            101        14           \n",
       "6 112053000           3.01              6            179        32           "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "data <- read.csv('ss-data.csv')\n",
    "head(data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the technical guide (Table 4.2, p. 36) there are six major factors used to compute the Smart Scale score, which are themselves calculated from more basic measurements:\n",
    "\n",
    "1. Congestion Mitigation\n",
    "2. Economic Development\n",
    "3. Accessibility\n",
    "4. Safety\n",
    "5. Environmental Quality\n",
    "6. Land Use\n",
    "\n",
    "Of these six factors only Land Use is directly given in the dataset. All other major metrics are calculated from their basic constituent measurements as specified on pages 27-31. Within the data there are four area types: A, B, C, and D. Each of these area types has a unique factor weighting used in calculating in its Smart Scale Score. These weightings are specified for each area in Table 4.2 on page 36 within the technical guide. It is noted that for areas C and D, the Land Use factor is not used (or equivalently, carries 0% weight). Upon inspection of the whole dataset it was apparent that this was represented by empty cells for projects falling within these area types. Accordingly, all such Land Use values are set to 0. By this precedent it also is taken that all empty cells in numeric columns correspond to zero. Accordingly, a value of zero is imputed to all empty cells.\n",
    "\n",
    "Based upon this, the following basic transformations are made to the data prior to analysis:\n",
    "\n",
    "* All columns with present/absent markers (x and no x) are recoded to 1 and 0, respectively. This converts them into a numeric equivalent (called binarization).\n",
    "\n",
    "* All empty cells in basic constituent measurement columns (which are used in calculating the factors) and other numeric values are changed to zero.\n",
    "\n",
    "* All empty cells in the Land Use column are changed to zero.\n",
    "\n",
    "Additionally, prior to imputing zeros to empty cells we inspect the percentages of empty cells in the numeric colums. All this is done as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Column_Name</th><th scope=col>Percent_Missing_Values</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>Throughput.Score             </td><td>30.94                        </td></tr>\n",
       "\t<tr><td>Delay.Score                  </td><td>25.99                        </td></tr>\n",
       "\t<tr><td>Crash.Frequency.Score        </td><td>17.08                        </td></tr>\n",
       "\t<tr><td>Crash.Rate.Score             </td><td>19.80                        </td></tr>\n",
       "\t<tr><td>Access.to.Jobs               </td><td>33.91                        </td></tr>\n",
       "\t<tr><td>Disadvantaged.Access.to.Jobs </td><td>33.91                        </td></tr>\n",
       "\t<tr><td>Multimodal.Access.Score      </td><td>41.34                        </td></tr>\n",
       "\t<tr><td>Air.Quality.Score            </td><td>34.16                        </td></tr>\n",
       "\t<tr><td>Enviro.Impact.Score          </td><td> 0.99                        </td></tr>\n",
       "\t<tr><td>Econ.Dev.Support.Score       </td><td>37.87                        </td></tr>\n",
       "\t<tr><td>Intermodal.Access.Score      </td><td>60.15                        </td></tr>\n",
       "\t<tr><td>Travel.Time.Reliability.Score</td><td>31.19                        </td></tr>\n",
       "\t<tr><td>Land.Use.Score               </td><td>56.44                        </td></tr>\n",
       "\t<tr><td>Project.Benefit.Score        </td><td> 0.74                        </td></tr>\n",
       "\t<tr><td>Total.Project                </td><td> 0.00                        </td></tr>\n",
       "\t<tr><td>Score.Divided.by.Total.Cost  </td><td> 0.74                        </td></tr>\n",
       "\t<tr><td>SMART.SCALE.Request          </td><td> 0.00                        </td></tr>\n",
       "\t<tr><td>SMART.SCALE.Score            </td><td> 0.74                        </td></tr>\n",
       "\t<tr><td>Benefit.Rank                 </td><td> 0.00                        </td></tr>\n",
       "\t<tr><td>State.Rank                   </td><td> 0.00                        </td></tr>\n",
       "\t<tr><td>District.Rank                </td><td> 0.00                        </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|ll}\n",
       " Column\\_Name & Percent\\_Missing\\_Values\\\\\n",
       "\\hline\n",
       "\t Throughput.Score              & 30.94                        \\\\\n",
       "\t Delay.Score                   & 25.99                        \\\\\n",
       "\t Crash.Frequency.Score         & 17.08                        \\\\\n",
       "\t Crash.Rate.Score              & 19.80                        \\\\\n",
       "\t Access.to.Jobs                & 33.91                        \\\\\n",
       "\t Disadvantaged.Access.to.Jobs  & 33.91                        \\\\\n",
       "\t Multimodal.Access.Score       & 41.34                        \\\\\n",
       "\t Air.Quality.Score             & 34.16                        \\\\\n",
       "\t Enviro.Impact.Score           &  0.99                        \\\\\n",
       "\t Econ.Dev.Support.Score        & 37.87                        \\\\\n",
       "\t Intermodal.Access.Score       & 60.15                        \\\\\n",
       "\t Travel.Time.Reliability.Score & 31.19                        \\\\\n",
       "\t Land.Use.Score                & 56.44                        \\\\\n",
       "\t Project.Benefit.Score         &  0.74                        \\\\\n",
       "\t Total.Project                 &  0.00                        \\\\\n",
       "\t Score.Divided.by.Total.Cost   &  0.74                        \\\\\n",
       "\t SMART.SCALE.Request           &  0.00                        \\\\\n",
       "\t SMART.SCALE.Score             &  0.74                        \\\\\n",
       "\t Benefit.Rank                  &  0.00                        \\\\\n",
       "\t State.Rank                    &  0.00                        \\\\\n",
       "\t District.Rank                 &  0.00                        \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "Column_Name | Percent_Missing_Values | \n",
       "|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|\n",
       "| Throughput.Score              | 30.94                         | \n",
       "| Delay.Score                   | 25.99                         | \n",
       "| Crash.Frequency.Score         | 17.08                         | \n",
       "| Crash.Rate.Score              | 19.80                         | \n",
       "| Access.to.Jobs                | 33.91                         | \n",
       "| Disadvantaged.Access.to.Jobs  | 33.91                         | \n",
       "| Multimodal.Access.Score       | 41.34                         | \n",
       "| Air.Quality.Score             | 34.16                         | \n",
       "| Enviro.Impact.Score           |  0.99                         | \n",
       "| Econ.Dev.Support.Score        | 37.87                         | \n",
       "| Intermodal.Access.Score       | 60.15                         | \n",
       "| Travel.Time.Reliability.Score | 31.19                         | \n",
       "| Land.Use.Score                | 56.44                         | \n",
       "| Project.Benefit.Score         |  0.74                         | \n",
       "| Total.Project                 |  0.00                         | \n",
       "| Score.Divided.by.Total.Cost   |  0.74                         | \n",
       "| SMART.SCALE.Request           |  0.00                         | \n",
       "| SMART.SCALE.Score             |  0.74                         | \n",
       "| Benefit.Rank                  |  0.00                         | \n",
       "| State.Rank                    |  0.00                         | \n",
       "| District.Rank                 |  0.00                         | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "   Column_Name                   Percent_Missing_Values\n",
       "1  Throughput.Score              30.94                 \n",
       "2  Delay.Score                   25.99                 \n",
       "3  Crash.Frequency.Score         17.08                 \n",
       "4  Crash.Rate.Score              19.80                 \n",
       "5  Access.to.Jobs                33.91                 \n",
       "6  Disadvantaged.Access.to.Jobs  33.91                 \n",
       "7  Multimodal.Access.Score       41.34                 \n",
       "8  Air.Quality.Score             34.16                 \n",
       "9  Enviro.Impact.Score            0.99                 \n",
       "10 Econ.Dev.Support.Score        37.87                 \n",
       "11 Intermodal.Access.Score       60.15                 \n",
       "12 Travel.Time.Reliability.Score 31.19                 \n",
       "13 Land.Use.Score                56.44                 \n",
       "14 Project.Benefit.Score          0.74                 \n",
       "15 Total.Project                  0.00                 \n",
       "16 Score.Divided.by.Total.Cost    0.74                 \n",
       "17 SMART.SCALE.Request            0.00                 \n",
       "18 SMART.SCALE.Score              0.74                 \n",
       "19 Benefit.Rank                   0.00                 \n",
       "20 State.Rank                     0.00                 \n",
       "21 District.Rank                  0.00                 "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Clean and recode numeric columns.\n",
    "numeric.columns <- colnames(data)[8:ncol(data)]\n",
    "for(col in numeric.columns) {\n",
    "    data[[col]] <- suppressWarnings(as.numeric(as.character(data[[col]]))) \n",
    "}\n",
    "\n",
    "imputer.f <- impute.0 # Missing value imputation function - can change if needed\n",
    "\n",
    "# Impute missing values to component scores and print out the percent missing in each column\n",
    "missing.ind <- 1\n",
    "missing.perc <- c()\n",
    "missing.titles <- c()\n",
    "for(i in 8:ncol(data)) {\n",
    "    missing.titles[missing.ind] <- colnames(data)[i]\n",
    "    missing.perc[missing.ind] <- round(100*(1 - (length(sort(data[[i]]))/nrow(data))), 2) \n",
    "    data[[i]] <- imputer.f(data[[i]])\n",
    "    missing.ind <- missing.ind + 1\n",
    "}\n",
    "\n",
    "# Properly binarize binary columns.\n",
    "data$Statewide.High.Priority <- as.numeric(sapply(data$Statewide.High.Priority, function(x){if(x == 'x') return(1); return(0);}))\n",
    "data$District.Grant <- as.numeric(sapply(data$District.Grant, function(x){if(x == 'x') return(1); return(0);}))\n",
    "\n",
    "# Display summary of missing values\n",
    "data.frame(Column_Name=missing.titles, Percent_Missing_Values=missing.perc)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the output above, the percentage of missing values in the key numeric columns ranges from 0% to 74%. We note that 74% of rows are missing a SMART.SCALE.Score and Project.Benefit.Score, meaning that only 26% of the projects have a Project Benefit Score and Smart Scale Score assigned - which is interpreted as only 26% of the projects having a Project Benefit Score and Smart Scale Score greater than zero."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "## III. Data Exploration\n",
    "\n",
    "In this section we explore the data further using descriptive statistics to provide several characterizations. We begin by looking at the overall distribution of Smart Scale Scores (SSS) ordered from least to greatest:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {},
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0gAAANICAMAAADKOT/pAAAARVBMVEUAAAAAv8QzMzNNTU1o\naGh8fHx8rgCMjIyampqnp6eysrK9vb3HfP/Hx8fQ0NDZ2dnh4eHp6enr6+vw8PDy8vL4dm3/\n//+EhyYwAAAACXBIWXMAABJ0AAASdAHeZh94AAAgAElEQVR4nO2di2KiyrZF6eMjMTFPt/7/\npx55KSoiWJO1qipj3ns6nbSOzKxibBDRFAdCSHAK7wKE5BBEIkQQRCJEEEQiRBBEIkQQRCJE\nEEQiRBBEIkQQRCJEEEQiRBCBSMUp5d8f3fqt9xafr4uiWG52o7/p7nXq/e9Wm/y97+fzpf74\n8XKcxsvHtB79N335FPQis8dapM/eW7y0hNFbzSVmzP3vVZv+ve/mZ/FTf2iIy58JPe7ctGGS\nyCMRafxt34u+W78Xi3Ir/jn+6/cz33TU/e/UfOJ7383Lpvzz6NGqRH6uij4JJop02LyE1iIG\nMRXpe1Us+269KJoN7r14vf3Xx9901P3v1Hzie9/LZ63istg0X9gUy9E9+lLd9FuwqySzRyxS\n8/fyP8bvt5tMcdzE2i/23Kv52/F/b4tiedx83o8f3quvb44Crr7qW+yWxUtzJNl7/+rbF6+7\nmzuWf34d/2nV3TIH7vt9fPC0eP0+f9fL+3+/Fsd//jqjVpU2n8V5F7IsJTjdt5zK202PPvTF\nTZerA4k+s4j0Vm3nmxuRyo2yT6TV6b/h9b/U9/96rT6UJi3Oj2KK0qLNpUiX9z/U/7o4XN/x\n+MdH/Xnn5vfv+1lcf9eL+39dP7T6Kd6qH7KzB/ks93LtfWv063WPHvTFTY/T5FFS/JlDpO+i\n+DhuRIveg5i+L/6U/+X/+D7donzUsikW9YdluSmVm+hbsar+dbW7xlze/7u8/64S5PKO5T+V\ne6fjEebn4/sef4q33WF33KZ/Tt+1e/9l+UMeN//T0dtbTV0W5/N/u/JfT/c9T6XL6UFfDvCz\nFpREHeVZu0P78Ljai5T/ie29dc8Xf+pjtZfmWKf8sGv+W1/evtk2m+/wdYu5uP9r9e2rbfjy\njmWz6vNd5+jr7n03za7qtfzYfNfu/W9+kJfi+6ZZfaja3Leaym5RXHNu0N2bllpxuiH+zCHS\nst0mRot0vPHn2/EhQrPr2HVu2Hz4rv795st99+/uEq7uuOx2fXDfZXNE9VPvVZqvne//Uu7J\nugddi55mnacEWvTLVY8edPemJWTROzISU+Y4tOt7GNR76+scH26/XatSfXhfXKp6B3Nx/ypX\ndyw6W+uj+57+2tnau/evny1qToZ07rDosHelA1eNV1c9etDdmz4YGYkkMYh0/lLPhld+eD/u\nLjYfP/dE6r9/mes73n7z+/ftF6l738/qbMjb1R16TjZc3Pe6x8APNPQfDRJX5hBp6qFd57H/\n7UZbH+qcHxf1bFxX9190D88u77i4edJ16L63h3Y39/9+PR93Nbuiz+J8uro5/d1tvLvu0YPu\n3vTAoV0SmUOk5rHy+1iRzme+PsqN8Fak9uzVHZGu7v/anKxY3N7xtX7S9fu8rd+/76Z5grY5\n2dD8/fr+nSLNyYabJ2SbWzT//H7dowfdvSknG9LIHCLVZ28/+h6J9B+nrIrFx3FP8LMpbv8L\nXu+R3k9ng09f/rl3/89i8d2cwr6+43f1PM33onP0dfe+xwPCTX36+/zsV/f+9envzXkrf2+o\n5SVC5Y7wq7lEqLnv+/EbtVPpcnrQ3ZuWnc6Pw0ismUOk5gnZi5MDN7fufv1n1d7h7dAn0nv7\nz18XJ7nu3L95OrPcadzcsXmSdXPve3fue/WEbH3rzv2bJ2TPl9P9tNif9gxcfdHq+dRB9V2u\ne/SgL256rMQTsvFnFpGq62xevkaLVL7+oHopQ/eJmO6H9/JZ06/PztM338uLRw7d+x9vfvzH\neqO+uePPZtE90zZ038tLhJobd+7/VV0i1N0xnp6bvXgZxem+H53rfs6cPvTFTblEKIXMeEJo\n99eO7T+Lr8c3mpwvLlpNIXOIVD1EKh8i/LVj+5fQ68f78vrH/nOUaOYQqX2I9OcOSX4Kxcts\nL7PjEVISmeXQ7rN8iLDqf6V11vnU7z14qXka4UlzQgRBJEIEQSRCBEEkQgRBJEIEQSRCBEEk\nQgRBJEIEQSRCBAkW6Xcgg/84ISpOfCAKKUAKEUKDSK4gCilAChFCg0iuIAopQAoRQoNIriAK\nKUAKEUKDSK4gCilAChFCg0iuIAopQAoRQoNIriAKKUAKEUKDSK4gCilAChFCg0iuIAopQAoR\nQoNIriAKKUAKEUKDSK4gCilAChFCg0iuIAopQAoRQoNIriAKKUAKEUKDSK4gCilAChFCg0iu\nIAopQAoRQoNIriAKKUAKEUKDSK4gCilAChFCg0iuIAopQAoRQoNIriAKKUAKEUKDSK4gCilA\nChFCg0iuIAopQAoRQoNIriAKKUAKEUKDSK4gCilAChFCg0iuIAopQAoRQoNIriAKKUAKEUKD\nSK4gCilAChFCg0iuIAopQAoRQoNIriAKKUAKEUKDSK4gCilAChFCg0iuIAopQAoRQoNIriAK\nKUAKEUKDSK4gCo0C7ff7wZtEEERyBVFoDGi/HzZJIUJoEMkVRKERoP3+gUkKEUKDSK4gCo0A\n5SLS4pi+j3Ushhnj6sbFybdQPiItmj+uPzaxGGaMqxsXJ99C+TxGQqQZQRQaBcrnrB0izQOi\nkAI0gxeTEyLS/8rMVYyQlDJGpPrkAnukGUAUUoBmMWNiOLRzBVFIAZrBi8lBJFcQhRSgGbyY\nHM7auYIopADNo8a0IJIriEIK0DxqTAtXNriCKKQAzaPGtHCtnSuIQgqQQoTQIJIriEIKkEKE\n0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIr\niEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQ\nQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAg\nkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hC\nCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE\n0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIr\niEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCY0D7/X74JhEEkVxB\nFBoB2u8fmKQQITSI5Aqi0GPQfv/IJIUIoUEkVxCFHoMQKefVjYyTbyFEEg4zwtWNjJNvIR4j\nCYcZ4+rGxcm3EGfthMOMcnWj4uRbiOeRhMNMcnVNOfkWQiThMJNcXVNOvoUQSTjMJFfXlJNv\nIUQSDjPJ1TXl5FsIkYTDTHJ1TTn5FkIk4TCTXF1TTr6FEEk4zCRX15STbyFEEg4zydU15eRb\nCJGEw0xydU05+RZCJOEwk1xdU06+hRBJOMwkV9eUk28hRBIOM8nVNeXkWwiRhMNMcnVNOfkW\nQiThMJNcXVNOvoUQSTjMJFfXlJNvIUQSDjPJ1TXl5FsIkYTDTHJ1TTn5FkIk4TCTXF1TTr6F\nEEk4zCRX15STbyFEEg4zydU15eRbCJGEw0xydU05+RZCJOEwk1xdU06+hRBJOMwkV9eUk28h\nRBIOM8nVNeXkWwiRhMNMcnVNOfkWQiThMJNcXVNOvoUQSTjMJFfXlJNvIUQSDjPJ1TXl5FsI\nkYTDTHJ1TTn5FkIk4TCTXF1TTr6FEEk4zCRX15STbyFEEg4zydU15eRbCJGEw0xydU05+RZC\nJOEwk1xdU06+hRBJOMwkV9eUk28hRBIOM8nVNeXkWwiRhMNMcnVNOfkWQiThMJNcXVNOvoUQ\nSTjMJFfXlJNvIUQSDjPJ1TXl5FsIkYTDTHJ1TTn5FkIk4TCTXF1TTr6FEEk4zCRX15STbyFE\nEg4zydU15eRbCJGEw0xydU05+RZCJOEwk1xdU06+hRBJOMwkV9eUk28hRBIOM8nVNeXkWwiR\nhMNMcnVNOfkWQiThMJNcXVNOvoUQSTjMJFfXlJNvIUQSDjPJ1TXl5FsIkYTDTHJ1TTn5FkIk\n4TCTXF1TTr6FEEk4zCRX15STbyFEEg4zydU15eRbCJGEw0xydU05+RZCJOEwk1xdU06+hRBJ\nOMwkV9eUk2+hPyMSIYQ9kjOIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQ\nIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuI\nQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBC\nhNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCS\nK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIK\nkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQ\nIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuI\nQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBC\nhNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCS\nK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIK\nkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQ\nIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuI\nQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBChNAgkiuIQgqQQoTQIJIriEIKkEKE0CCSK4hCCpBC\nhNAgkiuIQgqQQoTQIJIriEKPQfv9/tFNIggiuYIo9DD7/UOTFCKEBpFcQRR6lP3+sUkKEUKD\nSK4gCj0KIgmHGd3q5vuTRVcIkYTDjG518/3JoivEYyThMONb3Wx/sugKcdZOOMwIVzcyTr6F\neB5JOMwkV9eUk28hRBIOM8nVNeXkWwiRhMNMcnVNOfkWQiThMJNcXVNOvoUQSTjMJFfXlJNv\nIUQSDjPJ1TXl5FsIkYTDTHJ1TTn5FkIk4TCTXF1TTr6FEEk4zCRX15STbyFEEg4zydU15eRb\nCJGEw0xydU05+RZCJOEwk1xdU06+hRBJOMwkV9eUk28hRBIOM8nVNeXkWwiRhMNMcnVNOfkW\nQiThMJNcXVNOvoUQSTjMJFfXlJNvIUQSDjPJ1TXl5FsIkYTDTHJ1TTn5FkIk4TCTXF1TTr6F\nEEk4zCRX15STbyFEEg4zydU15eRbCJGEw0xydU05+RZCJOEwk1xdU06+hRBJOMwkV9eUk28h\nRBIOM8nVNeXkWwiRhMNMcnVNOfkWQiThMJNcXVNOvoUQSTjMJFfXlJNvIUQSDjPJ1TXl5Fso\nI5EWx/R9rGMxzCRX15STb6F8RFo0f1x/bGIxzCRX15STbyFEEg4zydU15eRbKB+RqiDSPCAK\nKUAzeDE5ISL9r8xcxQhJKeNEWhzYI80CopACNIMXk4NIriAKKUAzeDE5o0RadP9AJCGIQgrQ\nDF5MzhiRFuc/EUkLopACNIMXkzPmCdnOB0TSgiikAM3gxeSMeB5p0VzKwJUNehCFFKC55JgS\nrrVzBVFIAVKIEBpEcgVRSAFSiBAaRHIFUUgBUogQGkRyBVFIAVKIEBpEcgVRSAFSiBAaRHIF\nUUgBUogQGkRyBVFIAVKIEBpEcgVRSAFSiBAaRHIFUUgBUogQGkRyBVFIAVKIEBpEcgVRSAFS\niBAaRHIFUUgBUogQGkRyBVFIAVKIEBpEcgVRSAFSiBAaRHIFUUgBUogQGkRyBVFIAVKIEBpE\ncgVRSAFSiBAaRHIFUUgBUogQGkRyBVFIAVKIEBpEcgVRSAFSiBAaRHIFUUgBUogQGkRyBVFI\nAVKIEBpEcgVRSAFSiBAaRHIFUUgBUogQGkRyBVFIAVKIEBpEcgVRSAFSiBAaRHIFUUgBUogQ\nGkRyBVFIAVKIEJobkd5fiuKw+h4NsBhmkqtrysm3UKIi7ZbFMYei+BoLsBhmkqtrysm3UKIi\nvRabo0WHj2I1FmAxzCRX15STb6FERTpKdPrfuFgMM8nVNeXkWwiRhMNMcnVNOfkWSlSk5tBu\nU7yOBVgMM8nVNeXkWyhRkXaLosriZyzAYphJrq4pJ99CiYp0OLwti2K52Y0GWAwzydU15eRb\nKFmRpsZimEmurikn30IKkYpxD/mLU57w4Oo+q9GPjdpYDDPJ1TXl5FtIINLnUY3PMS4IRVpM\nZlgMM8nVNeXkW0gg0mvxMvrk2VMSVXe8/PR7tRl9mqGOxTCTXF1TTr6FBCIVxa72oyi+F6vD\nYfdaFK/VWYCvl6JYbC5vW/65K5bth+MXXopVJcX5fo9Fmr53sxhmkqtrysm3ULhIn8fd0Wt1\nbFcUq3LXVJ2ZXh7qY75juiY1G/6muv1H8Xb8wmt5Frv053Q/RIoORCEFaHgLLSX6rI7tamfe\nyj83xfvhsCw+jgdhF9t688l3dY3cS/FVyrc7rMq7nO83QqTpsRhmkqtrysm3ULhIlRvNH+Uh\n2rI+zHsp//z5fFv1iXR06Pv4yaI6HDzerNwRde/X811CJCpjMcwkV9eUk2+hYJGa47fyWK19\npHQ65lrdHH21n3wfhal2Y+19Lu43QqTdhidkDUEUUoAGN9DXZvt/vRXptVi+f/70inTc//xs\nuvJNFOmHS4RMQRRSgAY30EVR7hV29WFarciFNbt+kT6LzeJ0OPhTPmRaDh693Vy0Wp7p+1lx\n0aoNiEIK0ND2+dVsyq/ViYPyb5vypEH1irvy9au7/sdIpTfVCYfyz+Nt3rr3GyFSy+GsnQ2I\nQgrQ0Pa5aV7s/Vm/ruFwujD7u/y308Ha7YZ/fGz1UYtUHqFd3A+RogNRSAEa2j4Xi/Nfmq36\n5/iwaVXpVf3ljkjNOb7jn8fjs5/L+z0WiUM7WxCFFKCx2+qkfNXPvY7dpXCywRVEIQVovB0T\nsqovdH1SJE5/24IopACN3ljHpyiaswrPijQ5FsNMcnVNOfkWchJp0V7BgEhJgCikAIVuw4pc\ni/RSP5G75DGSCYhCCpDUiCdzJdKmvSKCs3YmIAopQFolnsvNK2Sr8+TfPI9kA6KQAqRV4rnw\nhKwriEIK0OAG+t9AJutyP1fCvFSvpd1teO9vGxCFFKDBDdRHpNMTsqN/r4vFMJNcXVNOvoUS\nFal9Qnb8O6BYDDPJ1TXl5FsoVZEmx2KYSa6uKSffQogkHGaSq2vKybdQkiLtNuWnH4vihUM7\nGxCFFKDBDdRDpEV51vuraN/Ha1Qshpnk6ppy8i2Uokjv5Vt4HZar5lW142IxzCRX15STb6EU\nRVrV7/Pw2rxVxLhYDDPJ1TXl5FvIQqTgEwVXjOpyho9qZ8SVDTYgCilAgxtoV5z9ft8jUqEw\nqYuo3n5oU727AyLZgCikAA1uoJceXZrUOiAW6aV8B7Bl+Ur1Ly4RsgFRSAEa3ECvPLowqXVA\nLNL78eHRZ/kOXrvVvbcKv43FMJNcXVNOvoXmF6k4KB4ldQnVG3eVJ77v//KK21gMM8nVNeXk\nWyhFkQ7fy/qp2PEnvxEpDk6+hWZ/jFTcevBMuETIFUQhBWhwA3101u75Xxt7EURyBVFIARrc\nQP8byNkAREobRCEFaHADfSBSK0CoCIjkCqKQAjS4gSKSnhMfiEIK0OAG+ujQThREcgVRSAEa\n3EARSc+JD0QhBWhwA/US6f2lKA6r0e99gkhRcPItlKhIu2V1Tr0o7v0+pZtYDDPJ1TXl5Fso\nUZFe618QePc3Zd7GYphJrq4pJ99CiYpU/xp0XkZhBaKQAjS4gSKSnhMfiEIK0OAG6npot+G3\nUdiAKKQACX14OtcnG/gdsqYgCilAgxuozx7pcHjjd8gagiikAA1uoF4iTY3FMJNcXVNOvoUQ\nSTjMJFfXlJNvoRRFKroZC7AYZpKra8rJtxAiCYeZ5OqacvItlKJIT8VimEmurikn30KIJBxm\nkqtrysm3UKoibTi0swRRSAEa3EB9RNrwGMkURCEFaHADvTBnfSvSPO8itCi+V8XPbsXLKGxA\nFFKABjfQC4/W61uRekWYnNuLVt+Kz8OOl1HYgCikAA1uoG4ifZbv+82hnQ2IQgrQ4AZ66dGl\nSR0DxCK9FB8/xbL89ZdjARbDTHJ1TTn5FpKLtO4VSf0YqTRoVT724mUUJiAKKUCDG+i1Rzci\nzfOWxZ/L8kVJE95F32KYSa6uKSffQnqR1lci9YswNTwh6wqikAI0uIH+N5COAYiUNIhCCtDg\nBuok0kv1hWLJK2RNQBRSgAY3UB+RNvWjrgknGwiJOg9FmunKhuqShm9Of9uAKKQADW6gj0QS\npeftuLofH8dimEmurikn30KJivRSvO4Oh92GS4RsQBRSgAY3UB+Rftq34xr9LvoWw0xydU05\n+RZKVKTjzqh6O67RJ+0QKQpOvoVSFWlyLIaZ5OqacvIthEjCYSa5uqacfAslKdJuU376sShe\nOLSzAVFIARrcQD1EWpRnvb+qkw2j37PYYphJrq4pJ99CwSIZpSvSe7E6+rNclRc4jL7822KY\nSa6uKSffQinukVbF8Yjup7w6aFcsxgIshpnk6ppy8i2UokjV5Qwf1c6IKxtsQBRSgAY3UAeR\nFuUnm6J8LhaRbEAUUoAGN1AHkV6K8iHS8lCecOASIRMQhRSgwQ3UQaT348Ojz+Lt+BBpVb6T\n0LhYDDPJ1TXl5FsoRZGq33tZnvguiuVogMUwk1xdU06+hVIU6fC9rJ+KnfDeJ4gUBSffQkmK\n9Ewshpnk6ppy8i2ESMJhJrm6ppx8C1mINMf72k2OxTCTXF1TTr6FhCL9q3MtUvM2JYEeIJIr\niEIK0OAG+lgkiQmI5AqikAI0uIE+EinYgDqI5AqikAI0uIEikp4TH4hCCtDgBuojUnv2YsHV\n3yYgCilAgxuog0iLopOxAIthJrm6ppx8C81++nsGkd47HnGtnQmIQgrQ4Ab6SKR5ztpNf2LK\nYphJrq4pJ99CBiLN8TzSavKb51sMM8nVNeXkWyjRKxsW7JFMQRRSgAY30BEiKXIlzvdqwpus\nVrEYZpKra8rJt1CiInHWzhZEIQVocANFJD0nPhCFFKDBDdRHpOmxGGaSq2vKybcQIgmHmeTq\nmnLyLZSqSBsO7SxBFFKABjdQH5E2PEYyBVFIARL68HRunkf6XhU/u1X9O5nHxGKYSa6uKSff\nQonukY57orfi87DjDSJtQBRSgAY3UDeRPssLVjm0swFRSAEa3EB9RHopPn6KZfk7ksYCLIaZ\n5OqacvItlKhIpUGr8lzD6ItXLYaZ5OqacvItlKhIh8/l4fBaTHirVYthJrm6ppx8C6Uq0uRY\nDDPJ1TXl5FsIkYTDTHJ1TTn5FkIk4TCTXF1TTr6FEhVp98qVDZYgCilAgxvohTlbK5FeuETI\nFEQhBWhwA30k0sTt/V5unpD9mAiwGGaSq2vKybfQ/CL1ijA5V/df8p4NpiAKKUCDG+hIkcRv\nx/Wz3OymASyGmeTqmnLyLfS73+8ffq+h+Ih0+OAxkiWIQo+y3z82aXADPUt0ioFInGywBVHo\nQfb7ESYNbqA+eyRONtiCKPQgqYr0wskGUxCFHiRVkQ4vr7xBpCGIQo+ie4xkfGjHYyRLEIUe\nRnrW7r+7IomfR0IkWxCFFKDBDfShSLNc2TA9FsNMcnVNOfkWSvSi1emxGGaSq2vKybdQkiLt\nNuWnH4viZfwZB4thJrm6ppx8CyUp0qI8VPwqjxgXoy8Ushhmkqtrysm3UIoivReroz/LVfmG\nq6PftMFimEmurikn30IpirQqjkd0P+UbCO2KxViAxTCTXF1TTr6FUhSpOgf4Ue2MOP1tA6KQ\nAjS4gTqIVP0C2U3xfUAkKxCFFCChD0+nK8xLUT5EWh7KEw6897cJiEIKkF6L6bk82fB6+Cze\njg+RVuX7f4+LxTCTXF1TTr6FUhRpt2hOfBfFcjTAYphJrq4pJ99CKYp0+F7WT8VOeMdiRIqC\nk2+hJEV6JhbDTHJ1TTn5FkIk4TCTXF1TTr6FEEk4zCRX15STbyFEEg4zydU15eRbCJGEw0xy\ndU05+RZCJOEwk1xdU06+hRBJOMwkV9eUk28hRBIOM8nVNeXkWwiRhMNMcnVNOfkWQiThMJNc\nXVNOvoUQSTjMJFfXlJNvIUQSDjPJ1TXl5FsIkYTDTHJ1TTn5FkIk4TCTXF1TTr6FEEk4zCRX\n15STbyFEEg4zydU15eRbCJGEw0xydU05+RZCJOEwk1xdU06+hRBJOMwkV9eUk28hRBIOM8nV\nNeXkWwiRhMNMcnVNOfkWQiThMJNcXVNOvoUQSTjMJFfXlJNvIUQSDjPJ1TXl5FsIkYTDTHJ1\nTTn5FkIk4TCTXF1TTr6FEEk4zCRX15STbyFEEg4zydU15eRbCJGEw0xydU05+RZCJOEwk1xd\nU06+hRBJOMwkV9eUk28hRBIOM8nVNeXkWwiRhMNMcnVNOfkW+l2P+F4RBJFcQRR6lPX6sUkK\nEUKDSK4gCj3IGpF0w4xudTP+yWIrtB5jkkKE0CCSK4hCw1nnJdKi/vOY7sc68w8zutUVgig0\nmPV6lEkzeDE5Y0Rq/Gn+OH9SZfZhCjnxgSikAM0hxtSMEGlxQKS5QBRSgGYxY2LGH9oh0gwg\nCilAM3gxOSEi/a/MXMUISSnskXeC56QAABSjSURBVFxBFFKAZvBichDJFUQhBWgGLyYHkVxB\nFFKAZvBichDJFUQhBWgGLyYHkVxBFFKAZvBicriywRVEIQVoBi8mh2vtXEEUUoAUIoQGkVxB\nFFKAFCKEBpFcQRRSgBQihAaRXEEUUoAUIoQGkVxBFFKAFCKEBpFcQRRSgBQihAaRXEEUUoAU\nIoQGkVxBFFKAFCKEBpFcQRQazhaRdMOMbnWFIAoN5V+VEd8rgiCSK4hCQ0Ek4TCVnPhAFBoK\nIgmHqeTEB6LQUBBJOEwlJz4Qhe5mu93WIm0ff68IgkiuIAoNhT2ScJhKTnwgCg0FkYTDVHLi\nA1FoODyPJBxmdKsrBFFoOIgkHGZ0qysEUUgBUogQGkRyBVFIAVKIEBpEcgVRSAFSiBAaRHIF\nUUgBUogQGkRyBVFIAVKIEBpEcgVRSAFSiBAaRHIFUUgBUogQGkRyBVFIAVKIEBpEcgVRSAFS\niBAaRHIFUUgBUogQGkRyBVFIAVKIEBpEcgVRaCjbkSCFCKFBJFcQhYaCSMJhKjnxgSg0FEQS\nDlPJiQ9EoaEgknCYSk58IAoNBZGEw1Ry4gNR6F62p4z4XhEEkVxBFBrIlj2SbphSTnwgCg0E\nkYTDlHLiA1FoIIgkHKaUEx+IQvey3+8RSTVMNSc+EIX60znVgEjZre4MIAr1Z19mux8JUogQ\nGkRyBVGoP4gkHKaeEx+IQv1BJOEw9Zz4QBS6k0qksSCFCKFBJFcQhfqzRyTdMPWc+EAU6s2+\nzliQQoTQIJIriEK9QSThMGfgxAeiUG8QSTjMGTjxgSjUn7NHiBQ8TD0nPhCF7uTkESKFD1PO\niQ9EIQVIIUJoEMkVRKH+dF7Nh0jZre4MIAr1B5GEw9Rz4gNRqD+IJBymnhMfiEL9QSThMPWc\n+EAU6g8iCYep58QHotBtrt4/CJGyWt2ZQBTqD3sk4TD1nPhAFOoPIgmHqefEB6JQfxBJOEw9\nJz4QhfqDSMJh6jnxgSikAClECA0iuYIopAApRAgNIrmCKKQAKUQIDSK5giikAClECA0iuYIo\npAApRAgNIrmCKKQAKUQIDSK5giikAClECA0iuYIo1Jvub+lDpNxWdw4QhXqDSMJhzsCJD0Sh\n3iCScJgzcOIDUag3iCQc5gyc+EAU6g0iCYc5Ayc+EIWus716XR8iBQxzLk58IAr1hj2ScJgz\ncOIDUag3iCQc5gyc+EAU6g0iCYc5Ayc+EIV6g0jCYc7AiQ9EoducflffaJBChNAgkiuIQjfp\n/EaXsSCFCKFBJFcQha7T/R1jY0EKEUKDSK4gCl0HkYTDnIsTH4hC10Ek4TDn4sQHotBNeIwk\nHOZMnPhAFLpIdWVQedZuez5vh0i5rO6cIArdZnv5KSJltbozgSh0G0QSDnMeTnwgCt0GkYTD\nnIcTH4hCt0Ek4TDn4cQHopACpBAhNIjkCqKQAqQQITSI5AqikAKkECE0iOQKopACpBAhNIjk\nCqKQAqQQITSI5AqikAKkECE0iOQKopACpBAhNIjkCqKQAqQQITSI5AqikAKkECE0iOQKopAC\npBAhNIjkCqKQAqQQITSI5AqikAKkECE0iOQKopACpBAhNIjkCqLQTdbTQQoRQoNIriAKXWeN\nSLphzsWJD0Sh66yvTUKkjFZ3NhCFrrL+syIRosu6ineLZ8IeyRVEocvUIk0FKUQIDSK5gih0\nkfX61iREymV15wRR6Crbm68gUj6rOx+IQuf8K7P99wRIIUJoEMkVRKFzbn6f+WiQQoTQIJIr\niEKn7Nkj6YY5Myc+EIVOQSThMGfmxAei0CmIJBzmzJz4QBRqs0Yk3TDn5sQHolCbPSLphjk3\nJz4QhZqsb37l5QSQQoTQIJIriEJNEEk4zNk58YEoVGe9Xt/xCJEmD3N+TnwgCv2WlzTs7+6P\nEGniMIcHJeLEB6LQb0ekZ0EKEUKDSK4gCv0iUh3RMIcHJeLEB6LQb3W16n2PEGnaMIcHJeLE\nB6JQc9X3XY8QadIwHwxKxIkPRKFGpGOeBylECA0iuYIohEhNNMN8MCgRJz4QhVqRQkAKEUKD\nSK4gCiFSE80wHwxKxIkPRCFEaqIZ5oNBiTjxgSikAClECA0iuYIo9Nv3zkETQQoRQoNIriAK\n/SJSHdEwQ2eZKIhCv2tEqiIZ5qNBiTjxgSi03u8RqYximA8HJeLEB/rzhcqXIW3vXR00FqQQ\nITSI5Ar684X63xRyIkghQmgQyRX01wtVr4vd3r1gdSRIIUJoEMkV9McLrRHplOBhjkh0m1u+\nP5l5IR4jNVEM8+GgRJz4QH+9UPWCPs7alQkf5ohBiTjxgf5woX/DL4wdD0Kk0Yluc8v3JzMV\naeitGqY0UogQGkRyBf3hQoh0mcAZjEp0m1u+P5lZof12u0WkTgJnMCrRbW75/mQ2hcq9EY+R\nLhM4g1GJbnPL9yezFWngzYOmNFKIEBpEcgX90ULt+538Gz7zPbaRQoTQIJIr6I8Wqi6uK59B\nQqRTAmcwKtFtbvn+ZCaF1qdfKYZIpwTOYFSi29zy/cksCq3PIoWBTjeJIIjkCvqLhdZ1RuyM\nHoDON4kgiOQK+nuFjg+OGpFGmoRICa2uG+jPFVrv16fIGilECA0iuYL+WKF949F2rEXjGilE\nCA0iuYL+VqH99vEry59ppBAhNIjkCvpThfb79gWx4zVCpJEzUM0yUdCfKoRIQwmcgWqWiYL+\nVKHy8dHDV5aPAd3cJIIgkivobxVaVyZtx59ouAe6vkkEQSRX0N8qtK0u9p5yXHcHdH2TCIJI\nrqC/VWiiQvdB1zeJIIjkCvo7hZ448T26kUKE0CCSK+hvFWKPNJDAGahmmSjo7xSacsH3IKjv\nJhEEkVxBf6cQIj1I4AxUs0wU9GcKTXoJ0hCo9yYRBJFcQX+j0Ha93j/43eVBjRQihAaRXEF/\no1D5NCwiPUjgDFSzTBSUfaFKH0QakcAZqGaZKCj7QrVHiPQ4gTNQzTJRUPaFGpHWY95T9elG\nChFCg0iuoNwL7U9v0fCsR4g0bgaqWSYKyr3QvtkhTXiPhicaKUQIDSK5gvIu1D40eurSoC7o\n0U0iCCK5gnIuVL0eFpHGJnAGqlkmCsq40P6U5y5WndBIIUJoEMkVlG+h5g0anj5ZdwoiRbi6\n0YHyLXQWKRCESBGubnSgfAvtg5486gSRIlzd6EA5FqpOc+/3zZNHT570ntRIIUJoEMkVlF2h\ndfWCiSm/SEzQSCFCaBDJFZRboeNOaN38juVt+FHd2EYKEUKDSK6g3AodPdrXHiHSxATOQDXL\nREGZFVrvEenZBM5ANctEQZkVKn9ry1kks0YKEUKDSK6gvAqtq19/1D5GsmukECE0iOQKyqjQ\nun3q6F9lkuVPphAhNIjkCsqh0La6lu5fsyNqn4JFpIkJnIFqlomCMii0LkUq3yXo8pIgRJqY\nwBmoZpkoKP1C66NI++rZo8tLghBpYgJnoJploqB0C5W/5Whd/t+6vIThKNLV68kRaWICZ6Ca\nZaKgFAtVp+PWzbsDVfaUH/5dvaAckSYmcAaqWSYKSqtQuf/5V51VWNci7etr6tb15xdvt4VI\nExM4A9UsEwWlUWhbn5qrdjiNR81jovVpr3T9gglEmpjAGahmmSgo+kLlTuh3++/f6f3pqrPc\nrUiVWtveF/Ah0sQEzkA1y0RBURdq3Pn3r+NRKVL1XnX/2v3Qfvvb9/I9RJqYwBmoZpkoKN5C\n2/Z8woVI2+POab8+fbURqe9qIESamMAZqGaZKCiSQrUL5Rm3kzZHS5r3dqx/8+u/dbtH6l5K\nV79HECIhkjPIvdC6emV46cK+Pm9Q72iaPyqLOjuk5vccjXmLLUSaGIthum9us4E8Cm2bUwi/\n7S+KaB7m1M8I1Z7UIq3X3UO7+vR2uy9CpMsgkivIrNC2uiiuUmjbeZDTnDDYnjxqTiZUV871\nijT2d5Mj0sRYDDO67T+9n2y7XrcnD5pnfzoeVTmLtK+P6s6/0mj7230ghEh9QSRXkEGhf+0z\nP/v2UO501qASqXnDn/oRUvvWJe05uRpxKVJooYlBpKQ2NyfQTIXW9Vm46v9rM2p11pcirTt7\npPX5ZMLNC1yfebErIk2MxTCj2/5j+8m2v4dt/bFSqH5cU+5fyp1LpUW51zl+qfy47op0foxU\n7n/OH8Le+P4XkSbHYpjRbf8uP9m2vfC6Om2wbc5bV69naEXaNwbVH+o9TfVujbU66072p6eL\nLi5O6Dwha/eThYIUIoQGkVxBDzjNbqZ+tudoTPUioPrytu3pFxzX2Tba1AdsJ4+OX2mP21qR\ntrU03cvjzo+AopsQIgmHmeTqjuQ0W3B1QFbvaurnSNuLDMrjrvZUdecRTnU01hWpeaDzr3q/\n4LNI7Zv6lN+mMWdwXxPdhBBJOMyUVnd78bF+vN/8fd19NqeR5fil7pUD2/3pRQntRQbNfuby\nTEGd6quNR6U+9RNA9R6peYTUXtOzPVdDpDmCSFPSPkyp/txu289LL9bl5/UxUulIZxs/idL5\na/P39cW/N7++ofVofbGfac8MNJe+1a+pK/9ylKc+gV3fp/Mk6umMXFvYYEIuIIUIoUGk80bW\n/1zJuvPnWaTTq3TOZ46rRyr/9h0zRom0Pn2t61HzTov1C1DXrTTto53mt3h1qB35rg73pr1z\nMCI9mVRFah+Gb0+7iO325r+7zZe22/aZlUPHhTNoWz84Ob9BW2cz/Hed87Mu2/bTrkide6wv\n7txcWdA6c3brdEDWyNLskKqvdvToPEI6vSi1q+dp7zPuUrgJo44fpBAhNPOK1Gze2/NhT3NM\n1F3ubbMRb7ufnwxo7vzbHjM1B1R3RWo325pUfbJtH65fqHH+W3sqqzoJfHFa+Fak/aVIHY/a\nhzRnj/bV9647VGI0n5ysq69pO3VoRDkxr/Yz23+/7fM927bs6T8eh217SgGRXPKESItjzp8N\n/oADIl1toNvfbedLnQ24unNnmz7dqr2Mpf7H03HWv+Y9BaovNr+dvnm4fqHG/lKk+uaX/94j\nUsejdXNtzf60M6mPw86GrDsPgmoxzvuW5mvtJW2nZ3euziW0tcppVtOrD9W6/yVqRWo+DXwC\nCJGezHSRFqc/qgz+gO0e4bdWoHkxS7m9ba+2z1akdvO5EKm7TW9Pd+ls2A2vvtX5xNfpmrHz\nA/cW05Gku+3WgPL/r0VqRdmfpGn2SKenOCvs6VBrfSFSfYP1rUgX+a13Ne0FBo04nfOAN389\ni/T0hnq1ZCIOIj3M8yKdTsVe/6e++tr2pNDVBvx7sU1vT3fperRvj7P29bZ9tuJ0Wvjf/rT5\nt/+pP6e9dPNapOaxTLfpyZ1tK+mdx0j/rndnzQ9wuZVsf7tX5Fw9gjvtZnr+bfLmNiqI9GQM\nRTpJcn3M1B6/1P/ZP23xF1t5e61LK0PnYGp//mv5T/UO6eL+66E90rp/j3TxDM5l12YvsT1/\nofppt+2xV/3Z7+/F309/uRWpc5NBWQZH/dS95uMg0sOcRfpfmbu32/bktP0ft9dWh7MU52Oz\niz3SobtzGLFHun7K8t/2fCbgTL4U6fYx0s1Zu+pnOjRlttUPWH6+7fzA5Wft59vuKLp/6f7D\n7ack0cS1R9qf/OiI9Ns8pL8QafAx0v7yu2zbA66Lb1zf4PQoY9+eCDvn9ge684Rs84+dk4hj\nE90OILpCf2GPVGXwB3zuMdLwWbv29Pe+eVzenLWrtuN9q0Pz1/o2p3Pl955H+h16QvbqxqNX\nd1Si226jK4RIv7env5vNvHo00+4OGilqXeqb3H0eqbn14PNId56Q1Yh0/snG3tCIk28hRPrV\nPyHb/OsdkTSLMi7ZbrfRFUKk37sizTHLREEUUoB0Ojyfea9ssJtloiAKKUA6HZ5PqhetZgKi\nkAKkECE0iOQKopACpBAhNIjkCqKQAqQQITSI5AqikAKkECE0iOQKopACpBAhNIjkCqKQAqQQ\nITSI5AqikAKkECE0iOQKopACpBAhNIjkCqKQAqQQITSI5AqikAKkECE0iOQKopACpBAhNIjk\nCqKQAqQQITSI5AqikAKkECE0iOQKopACpBAhNIjkCqKQAqQQITSI5AqikAKkECE0iOQKopAC\npBAhNIjkCqKQAqQQITSI5AqikAKkECE0iOQKopACpBAhNIjkCqKQAqQQITSI5AqikAKkECE0\niOQKopACpBAhNIjkCqKQAqQQITSI5AqikAKkECE0iOQKopACpBAhNIjkCqKQAqQQITTBIg3l\n/q/zcwqFHiW6QhE26g0iuYZCDxNfo94gkmso9DDxNeoNIrmGQg8TX6PeIJJrKPQw8TXqzawi\nEfJXgkiECIJIhAiCSIQIgkiECDKjSJe//dw1dY+2kHux6yIU6kt8jYYyn0iL0x/uWZy7LCIo\ndl2EQn2pFi2qRoP5CyItDoiUVqFDu2gxNRrOXxDpEJlIh6siFLrbJqJGD4JITolruz3EVgiR\nTonpZ49PpMi2kkVkj0gWh+hG9CCI5JP4tpKY9kinAtE0ehhEckl0hSITqU5EjR4GkTyyOP8Z\nQ6HrIu6FDm2BuBoNBZG86kRUCJHCw5UNDm2a45ZoCnFlQ3i41o4QQRCJEEEQiRBBEIkQQRCJ\nEEEQiRBBEIkQQRCJEEEQiRBBEIkQQRDJNEWTl6+7NzDtQ1Rh3UxTnHLHJERKNKybaRpPdpti\nOXgDklpYN9OcPKn+8vVSFItN9dnPS/u34x+b4tOtIXkuiGSaiz3SZ32Qtym/vDj9rfRo41yT\nTA4imebiMdKy+Dgcvkt1imK1O7wXi0qkTfHmXZNMDiKZptXo9bv69OfzbVWL9HOo90bH3RLH\ndSkGkUxTHdp9Fqv6s1VtVXvEV/+tKN4dC5Ing0imqY3Z1Ca9Fsv3z59rkd4W1e6JpBVEMk1z\nsmFZPQyqPtldi3T4KF4cG5LngkimaUT6Lorv8pOvw251I9LxiO/DtSR5IohkmvZ5pLfy9Pem\n6HuMVGq22HmWJE8EkUxzekK2Orh7LYrV161IR81e/SqSp4JIhAiCSIQIgkiECIJIhAiCSIQI\ngkiECIJIhAiCSIQIgkiECPJ/QpdQ28EIFjcAAAAASUVORK5CYII=",
      "text/plain": [
       "plot without title"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sssr <- data[order(data$SMART.SCALE.Score),]\n",
    "sssr <- cbind(sssr, data.frame(ind=1:nrow(data)))\n",
    "ggplot(sssr, aes(x=ind, y=SMART.SCALE.Score, shape=Area.Type, color=Area.Type)) +\n",
    "    geom_point() +\n",
    "    labs(x='Rank', y='Smart Scale Score', \n",
    "         title='Fig. 1: Smart Scale Scores (Ordered)')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In Fig. 1 above, the SSS's appear to follow a power law, with a long tail on the left and a relatively small number of projects having extreme scores. The distribution of scores clearly does not follow a linear pattern. Below we inspect the percentile thresholds from this sample by SSS:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Percentile</th><th scope=col>Score</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td> 5    </td><td> 0.149</td></tr>\n",
       "\t<tr><td>10    </td><td> 0.283</td></tr>\n",
       "\t<tr><td>25    </td><td> 0.737</td></tr>\n",
       "\t<tr><td>50    </td><td> 2.301</td></tr>\n",
       "\t<tr><td>75    </td><td> 7.347</td></tr>\n",
       "\t<tr><td>90    </td><td>25.141</td></tr>\n",
       "\t<tr><td>95    </td><td>44.623</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|ll}\n",
       " Percentile & Score\\\\\n",
       "\\hline\n",
       "\t  5     &  0.149\\\\\n",
       "\t 10     &  0.283\\\\\n",
       "\t 25     &  0.737\\\\\n",
       "\t 50     &  2.301\\\\\n",
       "\t 75     &  7.347\\\\\n",
       "\t 90     & 25.141\\\\\n",
       "\t 95     & 44.623\\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "Percentile | Score | \n",
       "|---|---|---|---|---|---|---|\n",
       "|  5     |  0.149 | \n",
       "| 10     |  0.283 | \n",
       "| 25     |  0.737 | \n",
       "| 50     |  2.301 | \n",
       "| 75     |  7.347 | \n",
       "| 90     | 25.141 | \n",
       "| 95     | 44.623 | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "  Percentile Score \n",
       "1  5          0.149\n",
       "2 10          0.283\n",
       "3 25          0.737\n",
       "4 50          2.301\n",
       "5 75          7.347\n",
       "6 90         25.141\n",
       "7 95         44.623"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sss.ord <- sort(data$SMART.SCALE.Score)\n",
    "percentiles <- c(0.05, 0.1, 0.25, 0.5, 0.75, 0.9, 0.95)\n",
    "cut.scores <- sapply(percentiles, function(p){sss.ord[ceiling(p * nrow(data))]})\n",
    "percentile.cutoffs <- data.frame(Percentile=100 * percentiles, Score=cut.scores)\n",
    "percentile.cutoffs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Because only a small number of projects exceed a SSS of 30, we will truncate all subsequent exploratory visualizations at 30 to best enable the majority variation to be seen clearly.\n",
    "\n",
    "### 3.1 Descriptive Statistics for Majority (Lower-Scoring) Projects (SSS < 30)\n",
    "\n",
    "To drill down further, we look at the distribution of Smart Scale Scores by area type and district. We begin with area type, looking at the number of projects within each of the four area types together with their average and median smart scale scores:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Area_Type</th><th scope=col>Num_Projects</th><th scope=col>Mean_Smart_Scale_Score</th><th scope=col>Median_Smart_Scale_Score</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>A        </td><td>114      </td><td>16.287246</td><td>2.9795   </td></tr>\n",
       "\t<tr><td>B        </td><td> 88      </td><td> 8.776795</td><td>4.3685   </td></tr>\n",
       "\t<tr><td>C        </td><td> 86      </td><td> 6.257140</td><td>1.9760   </td></tr>\n",
       "\t<tr><td>D        </td><td>116      </td><td> 7.978250</td><td>1.4635   </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|llll}\n",
       " Area\\_Type & Num\\_Projects & Mean\\_Smart\\_Scale\\_Score & Median\\_Smart\\_Scale\\_Score\\\\\n",
       "\\hline\n",
       "\t A         & 114       & 16.287246 & 2.9795   \\\\\n",
       "\t B         &  88       &  8.776795 & 4.3685   \\\\\n",
       "\t C         &  86       &  6.257140 & 1.9760   \\\\\n",
       "\t D         & 116       &  7.978250 & 1.4635   \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "Area_Type | Num_Projects | Mean_Smart_Scale_Score | Median_Smart_Scale_Score | \n",
       "|---|---|---|---|\n",
       "| A         | 114       | 16.287246 | 2.9795    | \n",
       "| B         |  88       |  8.776795 | 4.3685    | \n",
       "| C         |  86       |  6.257140 | 1.9760    | \n",
       "| D         | 116       |  7.978250 | 1.4635    | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "  Area_Type Num_Projects Mean_Smart_Scale_Score Median_Smart_Scale_Score\n",
       "1 A         114          16.287246              2.9795                  \n",
       "2 B          88           8.776795              4.3685                  \n",
       "3 C          86           6.257140              1.9760                  \n",
       "4 D         116           7.978250              1.4635                  "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "dt1 <- data.frame(table(data$Area.Type))\n",
    "m1 <- data.frame(aggregate(data$SMART.SCALE.Score ~ data$Area.Type, FUN=mean))\n",
    "mm1 <- data.frame(aggregate(data$SMART.SCALE.Score ~ data$Area.Type, FUN=median))\n",
    "dt1 <- cbind(dt1, m1$data.SMART.SCALE.Score, mm1$data.SMART.SCALE.Score)\n",
    "colnames(dt1) <- c('Area_Type', 'Num_Projects', 'Mean_Smart_Scale_Score', 'Median_Smart_Scale_Score')\n",
    "dt1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By comparing the mean SSS to the medians, it is apparant that the few high-scoring projects significantly inflate the means. Below we look more closely at the spread of the SSS by area type (truncated at 30, as discussed above):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0gAAANICAMAAADKOT/pAAAAM1BMVEUAAABNTU1oaGh8fHyM\njIyampqnp6eysrK9vb3Hx8fQ0NDZ2dnh4eHp6enw8PD//wD////3XPWHAAAACXBIWXMAABJ0\nAAASdAHeZh94AAAgAElEQVR4nO3d6YKiOhCG4bCIisrx/q/2yBINNjpIVbCC7/NjehkgZTqf\nQkB0VwBi7tsFAFtAkAAFBAlQQJAABQQJUECQAAUECVBAkAAFBAlQQJAABQQJUECQAAUECVBA\nkAAFBAlQQJAABQQJUECQAAUECVBAkAAFBAlQQJAABQQJUECQAAUECVBAkAAFBAlQQJAABQQJ\nUECQAAUECVBAkAAFBAlQQJAABQQJUECQAAUECVBAkAAFBAlQQJAABQQJUGA1SC7Q//jvdS5V\n5lxZT/5fvbv9X16dF5SyW7altyULynlqo1NUl3eNPj2CYaF/9enuzRZfljNnseq2XDV7o8HG\nR0PCGIMldRYEqR4Wn/obVX5b+08LOWXjlmdv6V3Jy8v504ZXv2z0+RHMC5JfK0KQlkWBIC3y\neZCa+/KHP/93eGzs9Hkhy7b0pmRBOX/aGG9pqtE/v5sVJP+/+kE6duUeZ2/1sXGC9LnPO2t/\n28M5X5vdbZfpz//lt9ep5rbvV9x2/WSFzN/Sm0cgKGeyjaa+bTH7tJB5QVpQzj/cHvXtFbn4\neOvzW/gGs3WNe8z/eLnlpDhM9mcxelIeL+F/aoIn40Pu8tsKh8wVw8vCsXTtgctlWOKSu+rP\nE+Dzlq5NdRvERT25hWGBzGX+GGbORm4HULff7+r7sm0h7zfUZbN+/KLZ33rDlcfr/Yk83NK4\nE7LdZbSt7uvEvsBzVdf61sZudJQ32uIj2tnob3Fp/+NW7uOoburxhT059YCftj607IbHMtVX\nkaUVpFP/9y2mglSW4b7IeIlbF5fhTtSwjdvfsHrsFhXhTlL7Z/QLhZt63tIl65eoprYQLhCu\n9HYj962UvtS8e/6e2tCjtNsB4u7+C79ou9ooErn/1agTxnuFL4L0p6rhSeY8LuexxZ0/bqvH\nh61Ve2S4v/9u+vGNenLqAT9tvS2sWz67vOr0yNIKkh8gU0EanPq/9XiJ7qAk2x39U9R9K9lj\ndNyWKZru71zelzj+aex5S/eK6okthAtkb8oJN3It78WVj1KP0xt6lNb0e7T9L3bdCk3RHi6O\nItFvyd2DFG7SvQ3SZFWt3bicxxZPPnTleEDfHkfTlpuF6zw9vnFPTj3gp60/Ws5fdXpkhoPk\nDT9eu8PU7DbWjpl7HaTcP1WN+HGQ3ye3bn+odjzn5+5Lv2KwS9Yvcf2T6Kct3dbNzt2QzSe2\ncPVDotmNj67fbKSdezzc1ti7+85aV8jkhoLiwkZdX8goXH8fkut7s876TY6C9PzjRFXturtx\n74y3OHTHZXzQWveHR2Uwz/j38Y17ctzC8N14623Lp77l+lWnx5VUkHzvH18HqXhxGNsejneq\nYXOn0Zdxy/2/LyaUR1saKmry/WViC/0Cw+At521k52cdq/vOWlfI5IZeBKl9otnVzws9PyT/\ncz1+DZ8M0ouqmj9dF2xx2H2rxtOow4Otg53EWY/v+RfjrbshM/0e7otOjyqpIN2PWl8G6Zaj\nrJn+r8tx173iH67Pg+W+tcuxKtx9RDavmgq29Py/4y2MHkc2byP3hi/jQiY39CJI+37B3fip\n4PkhuXFvvg3SRFUTvTPe4rD7lvlVO/d9uuzxuCYfX9CTUy08bT1oOX/T6REZDtLfH5/+9H+0\n88n5u6maSzna2xl/OfpXiVlNlf5PFvxyYgtPTwj/3sjjhz9D9u+GRuMqeGB+tjG7XCcG/bIg\nvaoqeDxPW+xefOrxYVRwCq1/oZp+fKOenGrhaevjlt92eiRJBekfr0hN9mK/LnhSHI2D0Zd2\nfzHfHc5vg/S8pdH/Tm0hm6j07UZc8NyfBe1PbShYczRrd+uKYz/vVVwnBv3jYYevUP8K0mRV\nf4IUbrE7HCrGx6z3gLjh0Gny8Y17cvIBj7c+bnmyryJLKkj/OEYqXl3Ctbs/K473TEZf8sff\n5Pp31L3YUhEe3kxtoZyY+3i7kfLv0Uj349SGguKK0XmkTr178erxeNj+uGI4RmrH4mkySK+q\n+hOkcIvtgK7Hu1cnFwoPUEePb9yTkw94vHX/+ta3PNlXkSUVJD9rN/2avXt5urydc+pO1bWz\nOsE4GH0Zfph6RWpebmk04Ta1hWM/nXQMa3u7kYn5sW6lqQ3d//NUjuaw8/txh3/1eDqO8A+7\nHfdtEe0lf1n3NHTKgk55rPWqqj9BCrc47GGG1xNWj4m0Q/+0N/n4xj05biHc1n3rQcuHF30V\nWVJBejqPNF7m/Ob/HmdB+knTqSD1r2d19rR++3P1ekv3ig4vtnBf4DRvI49zkU97/1MbCrYT\nPoW3ez2X+8SWfwRTQep18zO7Uf89r/WiKvcnSMEWuxfc0VTDNfixmeho//jGPTle//59uPWg\n5Vd9FVlaQRou8C7cRFjGA+Fp/fs4yMJdhtEXv9eR+TN8j80Wr7d0Ci5KmNzC5DXpbzby+M+n\ntzFMbejxkMcz236yoQgfgfsbpDJY9zJsfvjf57Wmq3LPQQq32B8RhTPQx3DmoQzPDo8f37gn\nxy08fgi27vyVFv2VDe/eCBBJWkG6nnfdNWlTYXn3ajW8AciV+6kJ4OFLu+1sd750f53H+uV4\n2ulpS91VXf5NUJNb6K6je36X1JuNDP8ZXtV2fbUh/4CLqhkv3R0fFYfRI5gI0vVQPK5P6/r2\nGB60jNaarOpPkEZb7HbCw4JHMw91uBf39PhGPfnUwl2w9fb3x9xlVfN3W+uwGqS3mjVPEGC5\nQ9Q/VLD15+fN9X27/Y+4/lD1LH77AVZxzmLuXYVbJ0gfeRzxrj+9iU/1fyjpu+nnbZ0gfeT+\n9oA1jyKxUPeHEr6Xfu7WCdJnmn07K5TteD1KQN7NXayzdYIEbAJBAhQQJEABQQIUECRAAUEC\nFBAkQAFBAhQQJEABQQIUECRAAUECFBAkQAFBAhQQJEABQQIUECRAAUECFBAkQAFBAhQQJEAB\nQQIUECRAAUECFBAkQAFBAhQQJEABQQIUECRAAUECFBAkQAFBAhQQJEABQQIUrBAkByRmwSjX\nD84XmoCUe/r62wgSFiJIIYKEpVzw788jSFiqOy5YcnCwRQQJyy06xt4mggQoIEiAAoIEKCBI\ngAKCBCggSIACggQoWDVIp33ZXd9XVqdYTQBfsWKQmjy4VraI0gTwJSsGqXLZ8dx9d6kzV8Vo\nAviSFYOUufP9+7PLYjQBfMmKQRpdl/X+Ii2ChMTwigQoWPcYqb5033GMhK1Zc/q7CGbt8iZK\nE8B3rHseqerOI2XlnvNI2BaubAAU2AmS8N5Gm0RPJOMLQTpkLj/EbWIbuCdCQtYM0rl02eG6\n5xKhmbhLT0JWDNK5S1Dlds31Urq3r0kMniv3jUvLikHateeOqv5MbOPyGE1sCkFKyeqXCLky\n+EG7iU0hSClZPUjHfp+OS4T+jWOkhKy6a7fzlzM0Oy4R+jdm7RKy5hv7svugcO9fkAjSgPNI\nyVj1PFLl45O9fT0iSEiOnSsbVm4C0ESQsF0r7hoTJGzVqpM1BAlbterpA4JkGbN2Auue0CZI\ndnEeSYQgEaQeVzaIECQGT4dr7YQ4RmLstAiSELN2jJ0WQRLjPFL8JhLAMVJCCJJdzNolhCBZ\nxnmkZBAkQAFBAhQQJEABQQIUECRAAUECFBAkQAFBAhQQJEABQQIUECRAAUHCdvE2ivhNYPN4\nYx9BggLeak6QIMfNT9YLEu/32TCCtNpj5x2oW0aQ1n3sBGmrOEZa58Fzl56NY9aOIEEF55Hi\nN0GQoOmHg8QxEvT8dJCYtYOWXw4S55Gg5reDBCghSIACggQoIEiAAoJkGZMhySBIdjE9n5Df\nDpLtZ3xOGCfkl4Nk/BmfS5hS8tNBWrGtBQhSSn44SNYHqvX6ECJIdgeq8VdMhAiS3ZFq/BgO\noR8OUgLP+LZnFRH46SDxjA8tvxwknvGh5reDBCghSIACggQoIEiAAoIEKCBIgAKCBCggSIAC\nggQoIEiAAoIEKCBIgAKCBCggSIACggQoIEiAAoIEKCBIgAKCBCggSIACggQoIEiAAoIEKFg/\nSIfcubKO2gSwthWD1N+MsXCdKkoTwJesHaTKVc31eqncIUYTwJesHaTMNe33jctjNAF8ydpB\n8jfbfn/TbYIEDSve233tIO18kLIYTQAPq37ayKpBKveH2h1v3zbV+9kGggS5VT//atUg9bpv\nsyZGE8Ddup/IuOZ5pPP5cCjLbsqhepsjggS57QbJVBPYul8NkgvFaQI/ZavHSNdm51wxXBzE\n9Ddi2+qsXZN1LzZlvxGChOi2eR6puyyoOWRFtxGChC1ZMUhZv+Ilyy8ECRuz+tXftxeloiBI\n2JgVg5Q7f/IoLwgStmXFIB3cbvju4gqChE1Zc/q7uqen/sd0CkFCYlY9IXsu/XeXHUHClti5\nsmHlJgBNBAlQQJAABQQJUECQAAUECVBAkAAFBAlQQJAABQQJUECQAAUECVBAkAAFBAlQQJAA\nBQQJUECQAAUECVBAkAAFBAlQQJAABQQJUECQAAUECVBAkAAFBAlQQJAABQQJUECQAAUECVBA\nkAAFBAlQQJAABQQJUECQAAUECVBAkLDcPz6c/pcQJCzVpYgo9QgSlnLBvz+PIGEh9/T1txEk\nLESQQgQJCxGkEEHCUhwjBQgSlmLWLkCQsBznke4IEqCAIAEKCBKgYONBclJahWDjNh4kMw1h\n4wgSoIAgWca+ZTIIkl2c8EzIjwfJNC7BSQhBMouLQlNCkMwiSCn58SBZHqQEKSUEyS6OkRJC\nkOxi1i4hBMkyziMlgyABCn48SIAOggQoIEiAgh8PEomFDoIEKCBIljH9nQyCZBcnZBNCkOzi\nEqGE/HiQLHPDrh2dsdiKu8YEySw37NrRGQutumtMkMwahgBBWmrVXeMfD5LlQUqQZNZ9PxdB\nMotdOxmCRJA6TDbIECSC1LM//W375XK7x0infdndULusTrGa+JDlcWD+hCz1hY2tskqnyYOb\n0xdRmtganvFltnkeqXLZ8dx9d6kzV8VoAita9xjEuhWDlLnz/fuzy2I0gRURpNCKQRq9zL5/\nzbX8gowBQQr9+CsSg0DA/jHSitY9Rqov3Xd2jpEYBQLWZ+1Wteb0dxHM2uVNlCY+xSgQsT2r\nuKp1zyNV3XmkrNxzHgnb8uNXNgA67ASJDxNHwtYMUrNzrqiHjdiY/gZ0rHmJUNZfaNdvxEaQ\nSCx0rDr9fbil6ZB1l9kRJGzKqidkuy+XLL8QJGzMFy4RaoqCIGFjVgxS7vxJ2LwgSNiWFYN0\ncLvhu4srjAQJ0LHm9Hd1T0/9j1NFBAmJWfWE7Ln03112BAlbYufKhpWbWLkhbBxBAhQQJECB\nOEh12c4blBeleqaaiIggQYc0SEV/rbbLVJNEkJAYYZAOrmjaID3OEalgfCMxwiBlronx1n2C\nhMQIg9Tt1hEk/DxhkPLhFenscrWSrhwjITk6x0h11r7XSA9BQmKks3blrJvii5qIiSBBh8p5\nJFcelcqZbCIiggQdXNkAKBAGqXx75+HFGN9IjML0dwQECYlRmP6OgCAhMcIgNWXxj9t4L8Ix\nEhIj3rWLcpdhgoTEECRAAdPfgAKCBCgQB+lYpHxlA6BD5R2y6V5rB+gQX/2dtR94lOzV34AO\n8QnZc/eV9yPht2ldIsT0N36a2itSplPP3yaiIkjQ8ePHSAQJOn581o4gQYf8PFLS75AFdPz4\nlQ2ADoIEKJAGqana6bqs0n1/H8dISIwwSJdsuM8qN9Ff1r7Ul+uHJwxS4Xbta1FTufLV4kv8\nTJD+xXp98LiywTTr9cFT+DSKVkOQorBeHzxhkCrX3fzkVDjVG9wxfgZ0RCp+/MoGQIfOlQ2F\n6pV2BAnJ+fETsiQWOgiSadbrgycJUlN1355yl6W6a2d9oFqvD54kSP1VDXXKkw3WB6r1+uAJ\ngtR+7OXtS5adr03hVN9IQZAG1uuDJwhS4drr605u3/2r+pLE+BnQEakQBKm/mKFyp8cPWhg/\nSIw4SLkLftBCkJAYQZDydtfu4nbt9w13EcJPEwSpaicbdq69i9D10OdJC0EaWK8PniBITXaf\n9z644f52SgjSwHp98EQnZHeuv+jbOd2LvwmSZ70+eCqXCLlS+YNkCdLAen3wfvxaO+voCJEV\nb2pBkLBVw315VmpslVUMNoHNc8G/KzUWexWDTazcEFbnnr6u01rcVQw2sXJDC1mvzzKCRJDu\nrNdnGUEiSHfW6zONYySC5Fmvz7S0Zu3qsi21VL31N+PHoyNEEjqPVPR3ck/1JvqADmGQ2reb\nt0FK9epvQIcwSO29vyPsilp+QQYmCIPU7dYRpGis1wdPGKR8eEU6u1ytpCtBurNeHzydY6Q6\nc6q3iCRIA+v1wZPO2pVpfxqF9YFqvT54KueRXKl6e0jGzx0dkYofv7IB0EGQAAWCIEX8oHqO\nkZAYgmSa9frg/fiunfWBar0+eOsH6ZA7V9ZRm5jP+kC1Xh88rSCdyn+v1604fAz6+xtKEqSB\n9frgSYNUzT9G6hapXNVcr5fq/ZUQjJ8BHZEKYZAeOfrHztp1CFJ7ufi1/fSKt9fmMX6QGPHb\nKI7tJ/ddCvfvmxaPrhJ//wpGkJAYhbdR7G+vRucZF9t12dn5IL39PCWOkZAYhSDV7fHOrGOk\ncn+ou09tbqr3sw0EaWC8vhXviWCdMEjlLReX2/HOaVaQ7rMSzmWNclXLWB8Hputb9S491gmD\nVLf92M1oz7hnw/l8OJRlN+VQvc0RQfJM1+eCf3+edPp73x/48EFjcViuzz19/W12LhGKduFe\nyix3BEEK2QnSyk1AiiCFpEFqqnYe+1/HPKImYBTHSAFhkC6Zn4VL806rjAIBZu0CwiAVbte+\nFjWVm3PR6uzDIII0MF4fR7N3Cidkx9+8diBIH7NeHzzxtXb9wVEz56npnM29aRdBGlivD574\n6u+ivVr1VMw6kXSee7qJIA2s1wdPOmtXfHSDyIM7x6pqm+iIVIjPIx3bG0QWqjcsZvwgOZyQ\nBRT8eJBILHRIgtRU3ben3GWp7tpZD5L1+uBJgtRf1VDzaRTxWK8PniBI7Wcj3b5k2fnaFE71\n8ygI0sB6ffAEQSpce33dye27f1VfkgjSwHp98ARB6i9mqPr7ByV672/r6IhUiIOUu+AHLYwf\nJEYQpLzdtbv0N2to3t9ea2kTQCIEQarayYZdf4vVw5ybn3zeRHQkFjoEQWqy+7z3wc28iO7D\nJqKzHiTr9cETnZD1dw/656dLLG4iNusD1Xp98FQuEXLlv+/8LWwiEusD1Xp98LjWzjTr9cH7\n8SBZR0ekgiABCggSoODHg0RioYMgmWa9PngEyTTr9cEjSKZZrw+e1p1WszQvWrU+UK3XB08p\nSBfeRhEFHZEKQZDq0a288y9XBXyT5BUpD3OkerUdQUJitI6RdHGMhMQIg1Tqvn1iqomorAfJ\nen3weEUyzXp98IRByp3uh8dONBGV9YFqvT54wiA1ZaH8nr4/TURlfaBarw+eeNdu1kdZSpr4\naXREKggSoODHr7UDdPx4kEgsdGgF6VRKK/lnEzFYD5L1+uBJg1SlfYxkfaBarw+eMEiPHNVq\nJV0J0p31+uAJg5S5Y/s5SZfCpXnRqvWBar0+eAqXCO1vr0bnRD9ozDo6IhUKQardQfuiO8YP\nEiO9+vu2a3dx+fVEkPDThEGq2wAV7WQDn4+EXyad/t63P+34WJdIrNcHjysbTDNen/Lpw5QR\nJNNM19eliCj1CJJpputzwb8/T3yJUJb0JULWWe4I9/T1t6ldIkSQfg1BColPyB7USnnRBIwi\nSCHuIoSlOEYKiHftuItQTKbrY9YuIJ1sKIqLVimvmojJ+igwXh/nke6kQarTnmywPg6s1wdP\nGKR94rN21geq9frgid/Yx6xdTHREKn581g7QId61S3vWDtAhfhsF9/4GFHbtmGyIyXp98AiS\nadbrg8fbKEyzXp9xK54wJkimWa/PtFUvYZIEqam6b0+5y5TPJjF+BnSEwKoX1UqClHVp7y8S\nUr0/JOMHcuu+zUMQpIMr2pNIWXa+NoU7frkqYCyZIBWuvfD75Pbdv2nespjEblcyQeoP46r+\n9vlMf0dhvT7TUjlG6rOTu+AHLQRpYL0+01KZtcvbXbtLf6/ixmWKRREkz3p9xqVxHqlqJxt2\n/SeMHbj3dxTW64MnCFKT3ee9D86dFYti/Hh0RCpEJ2T9zfNdsjfRB3SoXCLkSuX3UhAkJIZr\n7QAFBMk06/XBI0imWa8P3qpBOu3L7grXsvrHMRVBGlivD96KQWry4P207y/NI0gD6/XBWzFI\nlcuO/dmmS529ny9n/AzoiFSsGKQsOGl7fn9JEeMHiVkxSKPrnt5fBEWQjHBS334Aq/nxV6Tf\n+UPHQf956x4j1f1nwNg5RrI+EKhPJI2rvz9WBC/5+dtbHVt++KuiPoFU3o/0uVPVnUfKyj3n\nkeahPoFU3iEbEUEaWK/PsmTu2aDsNyd7/oGOWO5Xg7RyE9g6gkSQkmH6D8UxEsdIqTDdf1ud\ntfvgnDdBGlCfyDbPIx0iBEl8CYvUwr6Y/fgib1/Ken3rWXPX7pzNva/x/CD9910ECb1Vj5HO\nc282RJD844u8fWhZd7LhMPP2dwTJP77I24eWxGftth4kpIIgESSBX3/8DwSJIAn8+uN/IEim\ng2R9oH67PjunLwgSQRKgPkFLBIkgedQnaIkgzQ+SwrUTRnZd0kSQ5i5oPUjG69s6gjR3QeMD\n1Xp90EKQog5U6/VJEVSPIEUdqNbrkyJIHkGKOlCt1ydlPUgcI81d0PhAtV6fFEEStESQ5g9U\n6/VJESRBSwRp/kC1Xt/WEaS5CxofqNbr2zqCNHdB4wPVen3QQpCiDlTr9UkRVI8gRR2o1uuT\nIkgeQYo6UK3XJ2U9SBwjzV3Q+EC1Xp8UQRK0RJDmD1Tr9UkRJEFLBGn+QLVe39YRpLkLGh+o\n1uvbOoI0d0HjA9V6fdBCkKIOVOv1SRFUjyBFHajW65MiSB5BijpQrdcnZT1IHCPNXdD4QLVe\nnxRBErREkOYPVOv1SREkQUsEaf5AtV7f1hGkuQsaH6jW69s6gjR3QeMD1Xp90EKQog5U6/VJ\nEVSPIEUdqNbrkyJIHkGKOlCt1ydlPUgcI81d0PhAtV6fFEEStESQ5g9U6/VJESRBSwRp/kC1\nXt/WEaS5CxofqNbr2zqCNHdB4wPVen3QQpCiDlTr9UkRVI8gRR2o1uuTIkgeQYo6UK3XJ2U9\nSBwjzV3Q+EC1Xp8UQRK0RJDmD1Tr9UkRJEFLBGn+QLVe39YRpLkLGh+o1uvbOoI0d0HjA9V6\nfdBCkKIOVOv1SRFUjyBFHajW65MiSB5BijpQrdcnZT1IHCPNXdD4QLVenxRBErREkOYPVOv1\nSREkQUsEaf5AtV7f1hGkuQsaH6jW69s6gjR3QeMD1Xp90EKQog5U6/VJEVSPIEUdqNbrkyJI\nHkGKOlCt1ydlPUgcI81d0PhAtV6fFEEStESQ5g9U6/VJESRBSwRp/kC1Xt/WEaS5CxofqNbr\n2zqCNHdB4wPVen3QQpCiDlTr9UkRVI8gRR2o1uuTIkgeQYo6UK3XJ2U9SBwjzV3Q+EC1Xp8U\nQRK0RJDmD1Tr9UkRJEFLBGn+QLVe39YRpLkLGh+o1uvbOoI0d0HjA9V6fdBCkKIOVOv1SRFU\njyBFHajW65MiSB5BijpQrdcnZT1IHCPNXdD4QLVenxRBErREkOYPVOv1SREkQUsEaf5AtV7f\n1hGkuQsaH6jW69s6gjR3QeMD1Xp90EKQog5U6/VJEVSPIEUdqNbrkyJI3vpBOuTOlbVSE9YH\nqvX6pKwHaZPHSK5bsXCdSqcJ6wPVen1SBEnQkihIlaua6/VSuYNKE9YHqvn6vm3uH3qh7QYp\nc037feNylSbMD1TqE9Untd0g+eegv89Fi56orA8E6pPVJ7XdIO18kDKVJqwPBOqT1ZeOVYNU\n7g+1O96+bar3sw0EifoSs2qQ7rttzmWNShPWBwL1yepLx5rnkc7nw6EsuymH6m2OCBL16djk\nMVKMJqwPBOqT1SdFkOYuaHwgUJ+sPimCNHdB4wOB+mT1SRGkuQsaHwjUJ6tPiiDNXdD4QKA+\nWX3pIEhRBwL1yepLB0GKOhBWvwj0WeL9lw6CFHUgfDtHvx4kjpHmLmh8IHw7RwRpLQQp6kCg\nPll9UgRp7oLGBwL1yeqTIkhzFzQ+EKhPVp8UQZq7oPGBQH2y+tJBkKIOBOoT1vdtcwciQYo8\nEKhv0/UFlc5eUrJKtCasdzT1bbu+oNLZS0pWidaE9Y6mvm3XF1Q6e0nJKtGasN7R1Lft+oJK\nZy8pWSVaE9Y7mvq2XV9Q6ewlJatEa8J6R1PftusLKp29pGSVaE1Y72jq23Z9QaWzl5SsEq0J\n6x1NfduuL6h09pKSVaI1Yb2jqW/b9QWVzl5Sskq0Jqx3NPVtu76g0tlLSlaJ1oT1jqa+bdcX\nVDp7Sckq0Zqw3tHUt+36gkpnLylZJVoT1jua+rZdX1Dp7CUlq0RrwnpHU9+26wsqnb2kZJVo\nTVjvaOrbdn1BpbOXlKwSrQnrHU19264vqHT2kpJVojVhvaOpb9v1BZXOXlKySrQmrHc09W27\nvqDS2UtKVonWxLffiZz6QKA+WX1BpbOXlKwSrYlv5yj1gUB9svqCSmcvKVklWhPfzlHqA4H6\nZPUFlc5eUrJKtCasdzT1bbu+oNLZS0pWidaE9Y6mvm3XF1Q6e0nJKtGasN7R1Lft+oJKZy8p\nWSVaE9Y7mvq2XV9Q6ewlJatEa8J6R1PftusLKp29pGSVaE1Y72jq23Z9QaWzl5SsEq0J6x1N\nfduuL6h09pKSVaI1Yb2jqW/b9QWVzl5Sskq0Jqx3NPVtu76g0tlLSlaJ1oT1jqa+bdcXVDp7\nSckq0Zqw3tHUt+36gkpnLylZJVoT1jua+rZdX1Dp7CUlq0RrwnpHU9+26wsqnb2kZJVoTVjv\naItdt6MAAAhLSURBVOrbdn1BpbOXlKwSrQnrHf3td3mkPlCt1xdUOntJySrRmkino5eJvn3j\n/We9vqDSBZ3/+SrRmkino5eJvv1v+1d9yfx9CdJKHW3Tt3NEkOIiSCv5do4IUlwEaSXW+896\nfUGlCzr/81WiNZFORy/DZIPt+oJKF3T+56tEayKdjl6GINmuL6h0Qed/vkq0JtLp6GUIku36\ngkoXdP7nq0RrIp2OXoYg2a4vqHRB53++SrQm0unoZQiS7fqCShd0/uerRGsinY62yXr/Wa8v\nqHRB53++SrQm0ulom6z337dPcxEkIwPBOuv99+0cESQjA0Hq54+Rvm1+Ty7o/M9XidaE9YEg\n9fNBMl5fUOmCzv98lWhNpNPRyxAk2/UFlS7o/M9XidZEOh29DEGyXV9Q6YLO/3yVaE2k09HL\nECTb9QWVLuj8z1eJ1kQ6HW2T9f6zXl9Q6YLO/3yVaE2k09E2We8/6/UFlS7o/M9XidZEOh1t\nk/X+s15fUOmCzv98lWhNpNPRy3CMZLu+oNIFnf/5KtGaSKejlyFItusLKl3Q+Z+vEq2JdDp6\nGYJku76g0gWd//kq0ZpIp6OXIUi26wsqXdD5n68SrYlvX4pFkOIOVOv1BZUu6PzPVzHYBFrf\nfh4iSHERJHQIksxqQSKxMr++6xlUuuDBfb6KwSZWbmihX6+PIMkQpMGv10eQZAjS4NfrI0gy\nBGnw6/URJBnr40dN7Onj1BEkma2PD8ykcKZqpScqggQoWDVIp33ZxbysTrGa+BCJlbHef+vV\nt2KQmjx4ySyiNPEx6wPBOuv9t8kgVS47nrvvLnXmqhhNfMz6QLDOev9tMkiZO9+/P7ssRhMf\nsz4QrLPef5sM0mgK5O98yLLJkn+2udaszTZZ7z879f34KxKgY91jpPrSfWfnGAnQseb0dxG8\npOZNlCaA71j3PFLVnUfKyr2V80iADq5sABQQJEABQQIUECRAAUECFBAkQAFBAhQQJEABQQIU\nECRAAUECFBAkQAFBAhQQJEABQQIUECRAAUECFBAkQIHRIAGJWTDK9YPzNdYfC/XJmK7PdHEf\nsv5YqE/GdH2mi/uQ9cdCfTKm6zNd3IesPxbqkzFdn+niPmT9sVCfjOn6TBf3IeuPhfpkTNdn\nurgPWX8s1Cdjuj7TxX3I+mOhPhnT9Zku7kPWHwv1yZiuz3RxH7L+WKhPxnR9pov7kPXHQn0y\npuszXdyHrD8W6pMxXZ/p4oBUECRAAUECFBAkQAFBAhQQJEABQQIUECRAAUECFBAkQAFBAhQQ\nJEABQQIUECRAAUECFBAkQMF2gpS57NslvNbfmr04fbuO1867zO3qb1fxSt9/edV8u5CXNhOk\n+tbRZsfB/RM5zt8u5JVqGKqXbxcyzfdfZrS+DQVp5yq3+3YRL/UfFFK54tuFvLB32e1ZqNlb\nHal9/10Ksx24nSDdduwyuw9m+MSdJR+8s4aLD9DO6JOR77jc7F6H0b/sx46uuj3hH79dxis+\nSEYP4yq3779pysN3K3nBB6k2GvTtBKlwp+vJ8Av/sGtnc5jees/swVvPB6lx+XcLeWkjQWq6\n5/rMWZ3VGY6Vq2/X8YLVXc67e4FmK7Va14eO3Ri1u283BKkw+sRvdnh6BGkluWtP0ZztvvB3\n/VxnzuaZJLPD0yNI67jcP47a5uyt//ufjR7FlfdjpNrmzrHPz8VoB24lSPt7kPbfLmWaHwhG\nn1D3vt9ORl/Tfb8d7R5mfrsAFfnwSnQxOg78QGiMzn/fzyMVRucVH+eRbO4bbyRIZ1cO31md\nx+0HQlNYfULddVc2XEqjQefKhnVU9/PdtdGRer9WzOYhSPsMxLV2IpsIUpZNfWvJMAwMX718\nLJ0rrJ498FfPGz0Cbm0iSMC3ESRAAUECFBAkQAFBAhQQJEABQQIUECRAAUECFBAkQAFBAhQQ\nJEABQQIUECRAAUECFBAkQAFBAhQQJEABQQIUECRAAUECFBAkQAFBAhQQJEABQQIUECRAAUEC\nFBAkQAFBAhQQJEABQQIUECRAAUECFBAkc7L5H+TqAjFLwj/R/9bUt1DU/16sQ5DMoP+t2bnK\n7T5ZgQxZwB/BmtuOXfbRX4UgWcAfwZijq66V6z5e3Lkmd+Xtm0PuskP3v3XZfjb6eI0+SI3L\nr/7r7TfVfbHHyoiIIBlTuNP15Ir2W+dusbnFoeyOgdpf7fvDoXGShlek8rbitQ3i/vabvV8j\nWBkxESRbmm7KLnPNtU1I0X6p2y9N0c5AuPal6vi0Lzf8WPdHVjt3uf0mO1/PWbtwsDJiIki2\nHLuXm37fzvWvMWWXqqbbyetMB+ma9+nL29+0uanbNf6sjDgIki15F55zd8AzJGQ0wX2p98WL\nIB1uO3W3vcJ9sOLTyoiHDjblch/4l8kgFROp8D92e4X70YoEaTV0sCn7+8Afv7AMdi4/1JcX\nQbrtENbXfPRSxtz4auhnU/L2BeXavjLlwXTcfaag+83LIJ1dcW4DOBxcdbMPJdMM6yBIlpzv\nkwKFO/uEHNs5uNsRUNkn5PzqGKmNYdYF0c/a1aOVERNBsqS6v37UrronpD8wyi7t//dO4UqP\nINWuPyvrXLdKOV4ZMREkS7Is/PaekEPu3K6Lwu6WkVM9fn15BKnpDq2635S3o6nnlRERQdqQ\n2vWHWMwwrI8u35DC9a9CBGl9dPlmPC6pI0jro8s3I7sfOhGk9dHlgAKCBCggSIACggQoIEiA\nAoIEKCBIgAKCBCggSIACggQoIEiAAoIEKCBIgAKCBCggSIACggQoIEiAAoIEKCBIgAKCBCgg\nSIACggQoIEiAAoIEKCBIgAKCBCggSICC/wH1oCVOr9D0CwAAAABJRU5ErkJggg==",
      "text/plain": [
       "Plot with title \"Fig. 2: Smart Scale Score Distribution by Area Type\""
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "boxplot(SMART.SCALE.Score~Area.Type,data=data, main=\"Fig. 2: Smart Scale Score Distribution by Area Type\",\n",
    "        xlab=\"Area Type\", ylab=\"Smart Scale Score\", ylim=c(0, 30), col=\"yellow\") "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In Fig. 2 above, area type A has the greatest amount of variance while area type D has the least. Area types C and D have similar distribution profiles, with area type D having the lowest overall median.\n",
    "\n",
    "Below we shift our perspective to look at SSS with respect to district. We begin with a basic summarization of scores by district:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>District</th><th scope=col>Num_Projects</th><th scope=col>Mean_Smart_Scale_Score</th><th scope=col>Median_Smart_Scale_Score</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>Bristol       </td><td>42            </td><td> 1.160738     </td><td>0.4695        </td></tr>\n",
       "\t<tr><td>Culpeper      </td><td>35            </td><td> 3.597029     </td><td>2.1690        </td></tr>\n",
       "\t<tr><td>Fredericksburg</td><td>25            </td><td> 8.902320     </td><td>1.8640        </td></tr>\n",
       "\t<tr><td>Hampton Roads </td><td>52            </td><td>23.536212     </td><td>2.0735        </td></tr>\n",
       "\t<tr><td>Lynchburg     </td><td>28            </td><td> 9.429893     </td><td>1.7990        </td></tr>\n",
       "\t<tr><td>NOVA          </td><td>58            </td><td>11.844241     </td><td>3.3630        </td></tr>\n",
       "\t<tr><td>Richmond      </td><td>72            </td><td> 7.878069     </td><td>3.2880        </td></tr>\n",
       "\t<tr><td>Salem         </td><td>50            </td><td>10.277980     </td><td>2.5410        </td></tr>\n",
       "\t<tr><td>Staunton      </td><td>42            </td><td>10.463905     </td><td>2.7075        </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|llll}\n",
       " District & Num\\_Projects & Mean\\_Smart\\_Scale\\_Score & Median\\_Smart\\_Scale\\_Score\\\\\n",
       "\\hline\n",
       "\t Bristol        & 42             &  1.160738      & 0.4695        \\\\\n",
       "\t Culpeper       & 35             &  3.597029      & 2.1690        \\\\\n",
       "\t Fredericksburg & 25             &  8.902320      & 1.8640        \\\\\n",
       "\t Hampton Roads  & 52             & 23.536212      & 2.0735        \\\\\n",
       "\t Lynchburg      & 28             &  9.429893      & 1.7990        \\\\\n",
       "\t NOVA           & 58             & 11.844241      & 3.3630        \\\\\n",
       "\t Richmond       & 72             &  7.878069      & 3.2880        \\\\\n",
       "\t Salem          & 50             & 10.277980      & 2.5410        \\\\\n",
       "\t Staunton       & 42             & 10.463905      & 2.7075        \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "District | Num_Projects | Mean_Smart_Scale_Score | Median_Smart_Scale_Score | \n",
       "|---|---|---|---|---|---|---|---|---|\n",
       "| Bristol        | 42             |  1.160738      | 0.4695         | \n",
       "| Culpeper       | 35             |  3.597029      | 2.1690         | \n",
       "| Fredericksburg | 25             |  8.902320      | 1.8640         | \n",
       "| Hampton Roads  | 52             | 23.536212      | 2.0735         | \n",
       "| Lynchburg      | 28             |  9.429893      | 1.7990         | \n",
       "| NOVA           | 58             | 11.844241      | 3.3630         | \n",
       "| Richmond       | 72             |  7.878069      | 3.2880         | \n",
       "| Salem          | 50             | 10.277980      | 2.5410         | \n",
       "| Staunton       | 42             | 10.463905      | 2.7075         | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "  District       Num_Projects Mean_Smart_Scale_Score Median_Smart_Scale_Score\n",
       "1 Bristol        42            1.160738              0.4695                  \n",
       "2 Culpeper       35            3.597029              2.1690                  \n",
       "3 Fredericksburg 25            8.902320              1.8640                  \n",
       "4 Hampton Roads  52           23.536212              2.0735                  \n",
       "5 Lynchburg      28            9.429893              1.7990                  \n",
       "6 NOVA           58           11.844241              3.3630                  \n",
       "7 Richmond       72            7.878069              3.2880                  \n",
       "8 Salem          50           10.277980              2.5410                  \n",
       "9 Staunton       42           10.463905              2.7075                  "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "dt2 <- data.frame(table(data$District))\n",
    "m2 <- data.frame(aggregate(data$SMART.SCALE.Score ~ data$District, FUN=mean))\n",
    "mm2 <- data.frame(aggregate(data$SMART.SCALE.Score ~ data$District, FUN=median))\n",
    "dt2 <- cbind(dt2, m2$data.SMART.SCALE.Score, mm2$data.SMART.SCALE.Score)\n",
    "colnames(dt2) <- c('District', 'Num_Projects', 'Mean_Smart_Scale_Score', 'Median_Smart_Scale_Score')\n",
    "dt2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In a manner similar to above, we see the extreme high-scoring projects inflating the means. In Fig. 3 below the distribution of scores by district is shown."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0gAAANICAMAAADKOT/pAAAAM1BMVEUAAABNTU1oaGh8fHyM\njIyampqnp6eysrK9vb3Hx8fQ0NDZ2dnh4eHp6enw8PD//wD////3XPWHAAAACXBIWXMAABJ0\nAAASdAHeZh94AAAgAElEQVR4nO2di4KivBIGA17G0dHl/Z92uXhDAwhp4EusOudfd1bSHTop\nIYiOKwAgGLd2BwBSAJEADEAkAAMQCcAARAIwAJEADEAkAAMQCcAARAIwAJEADEAkAAMQCcAA\nRAIwAJEADEAkAAMQCcAARAIwAJEADEAkAAMQCcAARAIwAJEADEAkAAMQCcAARAIwAJEADEAk\nAAMQCcAARAIwAJEADEAkAAMQCcAARAIwAJEADEAkAAMQCcAARAIwAJEADEAkAAMQCcAARAIw\nIBqR3BPNj8Ntjhvn8kPHc7usfHL/N6Eru2mRersc0J2XHDWb/bkv6cseXDcaqumuJ2Jndz7Y\nZHqP+zdalJRF2l+HqOc5537GduSUtTN/HKmvy9O785bjxrEz6esefCbSrdVMIk3qcf9Gi5Kw\nSMfb5vv35w6PYKfxHZkWqafLAd15y9GO5Ev69m8fiXR7dj6Rxvd4bLIZiUmkkS12zm0vxa9z\n2ftzeanXpSjO5bnfNqwjn0fq2YOA7nhzXI65d7f7O/KZSBO688Em03o8fqPZiFWk24/nUpfN\nwVvE/Ta7FB0vt7efLk/PHnKXl6+Jh8xtroeF362rFi7n6xbn3N1OwbojFZd9OSU2R2+E6waZ\ny24rgk+ClAuo8t93x/u2VUf6A9VuHh//cPkpHXXb3+J+CHiO1C5Ctju3YtWPnnOB117VS9Jd\na5XXivgQJXsqYFCPn+r0NjALE7lIp6Z8m+4SllN/996+LP/2+STqGqMcp/3jJGPzfMpRCXHb\n6DnUa6Rz1myx90V43uC5UW+Qe5Ttrat5vfDzBXp07djsd/MPt02rZq1pmd/+qVWE9jlWh0hv\nvbq+yPy1u/OIuLutgo7PZ9tBPX6qEyJ9iF+kW7m7SlgdEOrj0kv7elGS7X5vL+f3KNljdpTb\nbC61idv7Fr9vyV4j3Xt09ER43iDr6c5zkHofnubstSP+QI+ulQe3/P4Pu7rBpZzUh/a0bCK5\nu0jPIV2vSN5eVeza3XlEPN2k2z7LH9Tjpzq9DczCxCSSeypW81AtgMq59pt1lXDrbi+Dnicq\n8vulonLKV/M5/6sfivow8XRK1mxRvBn9Eqlsm/3VEyD3RChucl2aifJJkOqayaFs8ePupz51\nR7yBnjr3nNQ1HWlN1fddck01j1kTsiXS64+eXlVtd+3qtCNey3FuemHQ41adWCN9hlek7VWT\n364iPl4yX6kWtzX763an1kM7RvNnx+XZVqRrjy75z9kTodngOhW2nwXZVS/KFfv7qU/dEW+g\njmlZvdDsjq8bve7S7edj6wjqF6mjV5e30j1F/Gn2bn9rGtzjdrER6SO8It1XrV1FvFzPDnyc\nf3f1qcGheJ0s92jn3/3G3cf30pXqKdLrs+0Irf3IPgtyT3xud8QbqGNa/jQb7tovBa+75NrV\n7BXJ0ytPddoRL01Xs1vT4B53JFuFmER6//Fl6L1cWmcSL5y3rdOC9sPv7SjxSaprpPaznggv\nLwjDQR4/vE3Z90CPH9onRberjdm58Ez6aSJ19eppf14i1keQ49syanKPEWkCXpF6j0j7bdb1\n3NOLYmtIWg/V+WK+O/z1ivQaqfWsL0Lm6WlvkOfX/uwpvy/QU8vWNbBylv42l882Ree07Hq9\n7xLJ26u3uf0c8Vjl37RXrSE9RqQJeEXqXSNdR+zseaNvd39VbJ+ZtB6u72v0H/xeI22eT9t9\nEbaeqx+9Qbbvq5H6R1+gp85tWu/K1Bx3HUePx27X1y0ea6TKgpNXpK5evc3t54iV/MeX4Qjp\n8YY10ni8It2u2jlfEX/qazrV3QJv9whV15zqNwmr60lP86D1cP3Bd0S6dEZqXUjyRai6fKof\nNp8F8Vwfqxv5At2fPG1b17Dz+7rqdvR4eqv6ebereV91orrlL6sLd8qeivJo1dWrN5GeI17P\n11r3E4b0+P2q3aVYi7hFenkf6WWb+5OXt+ce74I0l1l9IjUGHrOX2NXP++5I96SHjgj3DU6f\nBXm8q/vyzrIv0FOc54NhdU51vl82u+2BT6SG+q23Xau2r606euXeRHqKWB9wW5caAnvcqtPL\nwCxM5CJdb0zdOJ9It/e9j57293mQPT/bejjdtqin6r19Nb823ZFOTzcleCPc7qVtjXlPkMeT\nLx9j8AV6nZW3rW9L983zHrxOy/Jh+9T2fA1/ffa1lb9X7lWk54j1gebl3YigHrfq9DowyxK5\nSMXfrr7XyvlEut6M5r9qXX8AyG1/fBeArw9V7Gz3d67H/tF+277s9BKpTuq2x6IzQn1/2PZl\nfdMT5Prk811tRVeg25zcXPf6vnW92tgcWnvgEak4bK5nmbfa/j4vylqtvL16E6kVsT4Jb+95\nUI/bdXodmEWJRqReLj13DoMQh2QHKm6RXHNV6C/44wewCH/ZmsuYWYlbpMeK139HHSjRDFTo\np+lFiVuk+832qb7OJUU9UIGfpZclbpGKy0/9QYkdx6MIyOtrF4kSuUgAGiASgAGIBGAAIgEY\ngEgABiASgAGIBGAAIgEYgEgABiASgAGIBGAAIgEYgEgABiASgAGIBGAAIgEYgEgABiASgAGI\nBGAAIgEYgEgABiASgAGIBGAAIgEYgEgABiASgAGIBGAAIgEYgEgABiASgAGIBGAAIgEYgEgA\nBiASgAELiOQAImPCLLcXZ4UUkZFeRdzLY+QgUhykVxFEQiSwwD39GT+IBCtRLyumrC0kQSRY\njUlLdFEQKQ6oiDiIFAdURBxEigMqIg4ixQEVEQeR4oCKiINIAAYgEoABiARgACLFARURB5Hi\ngIqIg0hxQEXEQaQ4oCLiIFIcUBFx1hDpsHEu/5k1BcCyrCGSOxfFpTcKIkFkrCLSBZEgMVY5\ntcs5tRsLFRFnDZF+sh+XHd/Chn23UeJQkW4k5ss6a6Ty5C6bNUVyUJEuRL76YRWRivJwhEij\noCJdiHwZ0SqnduV/76d2pimSg4p0oPL1eLwhC1GDSCungDRApJVTQCJ88RpJIkVkUJEuvviq\nnUSKyKAi3Xzt+0gSKSKDioiDSHFARcRBpDigIrMTdoaISABF+DULRAIowq+iIxKshsTltobg\n93URKQ4SrIjIG0ANiPQlJFgRkVsSGhDpS0ivIio3yV1hjfQdpFcRNZG4avcVpFcRMZF4Hwki\nRWqNFAwiwUpIXbULBpFgNYTeRwoGkeKAioiDSHFARcRBpDigIuIgUhxQEXEQKQ6oiDiIBGAA\nIgEYgEgABiBSHFARcRApDqiIOIgUB1REHESKAyrSjcQte4gUB1SkC5GbyBEJ4kbkY02IBFGj\n8kFbRILxSKxKGhBp5RSRoVQRkVVJAyKtnCIylCoisiq5ItIbRIoDoYqoHAOuiBwfESkOhCoi\nJpLIig2R4kCoInIiSYBIMBaRVYkWiARjEVmVaIFIMB6JVYkWiBQHVEQcRIoDKiIOIsUBFREH\nkeKAioiDSHGQZEVSumaBSLASaV1FR6R5SelF15i03tdFpDlJ60XXlsTuNEKkObF70U2lIg8Q\nCZE+xXCuJFKRJxAJkT4FkfpgjWTfixVSLAEi9ZHWAhKR5oQ1Ui8pXdJEpDlJ60UXekCkeUnp\nRRd6QCQAAxApDqiIOIgUB1REHESKAyoiDiLFARURB5HigIqIg0gABiASgAGIBGAAIsUBFREH\nkeKAioiDSHFARcRBpDigIuIgUhxQEXEQCcAARAIwAJEADECkOKAi4iBSHFARcRApDqiIOIgU\nB1REHESKAyoiDiIBGIBIAAYgEoABiBQHVEQcRIoDKiIOIsUBFREHkeKAioizhkiuOGbZcdYU\nyUFFxFlHpKyo/j9jCoBlWUWk3xyRIC3WEOn8kxec2kFSrCHS0W09UdwzoSmSg4qIs4ZI2fm8\n74/CtHmFioizjkhFsUWkUVARcVY5tcua/2ZMkRxURBzekI0DKiIOIgEYgEgABiASgAGIFAdU\nRBxEigMqIg4ixQEVEQeR4oCKiINIcUBFxEEkAAMQCcAARAIwAJHigIqIg0hxQEXEQaQ4oCLi\nIFIcUBFxECkOqIg4iARgACIBGIBIAAYgUhxQEXEQKQ6oiDiIFAdURBxEigMqIg4ixQEVEQeR\nAAxAJAADEAnAAESKAyoiDiLFARURB5HigIqIg0hxQEXEQaQ4oCLiIBKAAYgEYAAiARiASHFA\nRcRBpDigIuIgUhxQEXEQKQ6oiDiIFAdURBxEAjAAkQBqnAuZdogEUFQa3f6Y2H6RJoIpIoOK\nzIx7+nN6+7mbCKaIDCoyL+7lcXKAeZsIpogMKjIviPQlUJF5QaQvgYrMDGskAAO4agdgAu8j\nAawOIsUBFREHkeKAioiDSHFARcRBpDigIuIgUhxQEXEQCcAARAIwAJEADECkOKAi4iBSHFCR\nFq6PdXq0SBPBFJFBRXpQKA4ixQEV6UGhOIgUB1SkB4XiIBJEj8J0QSQAAxAJwABEigMqIg4i\nxQEV6UGhOIgUB1SkB4XiIFIcUJEeFIqDSHGQZEWsbudRKA4iwUqEfpPcUySDGKEgEqxE6Heb\naoFIsA7uemqXyFgjUhykVxF3PbVLZM8QKQ7Sq8jVIAuRFIqDSHGQXkUQCZFWIL2KGJ7aKRQH\nkeIgvYoYXmxQKA4iwUrYXf5WmC6IBCth94asAogEq5HMte8CkWKBioiDSHFARXpQKA4ixQEV\n6UGhOIgUB1SkB4XirCHSYeNc/jNriuTQqojYVQKFzqwhkjsXxaU3ikJloAu569YKXVlFpAsi\nxUxaHyQyYpVTu5xTOx+9v2JB5/cvuJdHqFhDpKPbeqKs/5s51sb9mwgirc8aImXn874/yneO\nESJNRKEn64hUFFtEeiMSkfTWSApdWeXULmv+mzFFlEQjElft3uENWR1iEYn3kTwgkg7xiCSG\nwv4jkg6IFDGIpAMiRQwi6YBIEYNIOiDSRBT2H5F0QKSJKOw/IumASBNR2H9E0gGRJqKw/8Ei\nHbfVm3Pbs1F/fCm+BUSaiML+h4q0ae7VdpmpSQqVWR5EiphAkQ5uc6lEOridWZcKREKk6AgU\nKXOXOe5h/M6ZgUgREyhSfVqHSDYg0kQU9j9QpPx6RPpzuVmXCo3KLA8iTURh/23WSMfMHcy6\nVGhUZnkQaSIK+x961W57/ZKFjVWH3lN8DYg0EYX9N3kfyW1/jbrjTREzYz4Ch0gTUdh/7myY\nk3HXYRApYgJF2u7NetKVImbGfU0IIkWMweXvGUhkZriXx6HNESleDC5/z0AiMwORlkFh/wNF\numw3J7O++FNEDCItg8L+B5/azfItwwqVsYA10iIo7D8izQlX7RZBYf+5/D0vvI+0AAr7j0g6\nIFLEBIv0u+HOBiMQKWJMPiHLvXYmIFLEBN/9nR2Lgru/TUCkiSjsf/Absn/1I59HMgCRJqKw\n/1a3CHH5OxxEmojC/psdkXp/cVhIiu8BkSaisP+skXRApIko7D9X7XRApIgJfx+JT8ha8XUi\nif0GzSC4s0GHLxNJ7nc6B4FIOnybSE9/WkRal1CRLvvqcl22t/18n0Jllue7RBr5Ya1PQq1J\noEjn7HqA5kv0w0GkwFBrEijSxu2qY9Fl77ZWPXpN8T0gUmCoNeHOBh2+SyTWSM9NsuuXn1wQ\nKZxvE4mrdg/2rv7yk9PGmX7BXSrVHceXicT7SM8/cGeDHV8nUkrY3NmwMb3TDpEQaRQK+88b\nsjog0kQU9h+RdECkiSjsf4hIl33911PuMk7tDECkiSjsf4hIzV0NRy42GIFIE1HY/wCRql97\nWT5k2V9x2TjTD1IoVGZ5ECliAkTauOr+upP7qf80PSR958xApIgJEKl5N23vTo8frPjOmYFI\nERMsUu6efrDiO2cGIk1EYf8DRMqrU7uz21V/v/AtQuEg0kQU9j9ApH11sWHnqm8RKg6NT1Yo\nVGZ5EGkiCvsfINIlu1/3Prjr99sZoVCZ5UGkiSjsf9AbsjvX3PTtnO3N3xKVWR5EmojC/pvc\nIuS2xr9IVqEyy4NIEcO9djogUsQgkg6IFDGIpAMiTURh/xFJB0SaiML+I5IOiDSRyfvv+pi9\nD4g0E4g0EaP9DwqDSDog0kQQacUUgiDSuqwr0nFbnU1uTb/6G5EQaQVWFWnTLMv4En0DECli\nAkWqPm5eicTd3wYg0kQU9j9QpOq7v2f4DmeFyiwPIk1EYf8DRapP6xDJBkSaSAJX7fLrEenP\n5SG96EvxPSDSRBIQ6bpGOmbO9Csiv3NmINJEEhCp2PLbKMxApHVZ/30ktzX9ekhEQqQ14M6G\nRECkiEEkHRBpIgr7HyBS0F3n5r1KAESaiML+I5I5kz/kgkgTSeGq3TwkNDNG7AoiTQSRVkyx\nFIg0PymJdNqG9OKjFHGCSNGwqkh71kj9IFI0rCnSw6NjSC/6UnwPiBQxwR+j+K1+c99540y/\ntPg7ZwYiTURh/w0+RvFTHo3++NWX4SDSRBT230CkY3Xn95g1UnW3eNZ7KqhQGSNYI81PAlft\ntuWp3dnlxWmcSFlR9P6Gv4RmBiLNTwIiHSuB6i9AGfGdDe43RyTfpog0jQREKhdI5R+7cb9o\n7PyTF5zaeTYVFqn3vqdEhkviDdlwEhyZikREeu3q4hkXIDqR8uHbN0NTRAkiRUyoSJd9tdjJ\n9pcxEQbPA79znL5ApHnOEBWmS6BI5+z6ZVx802o4XyCSdJggAkXauF11LLrsHTet+klzjWSE\nlkhrrpHuh2NuWu0AkXpQMMAoTPC9ds3i6IJIHSBSDwoGGIUJvvt7U92tetqMeiNpVIrISVOk\nhAbowaqXvzfXCy58QWQHiBQN676P9Ft9QeTG9AuLEx2nQRApYnTubFg4hSDfJ5JWmCAQSQdE\nWjdMECEiXfb1X0+5yzi16yLNNZIRWiKttkZq7mo4crGhD0TqQcEAozABIlW/G6l8yLK/4rJx\npr+PItqZ8Q4i9aBggFGYAJE2rrq/7uR+6j/5zgY/aYqU0AA9WEuk5maGffP9QdzZ0AEiRcO6\nIuXu6QcrkhynQRApYgJEyqtTu3PzZQ0X1/sdDFNTfBffJ5JWmCACRNpXFxt2zVesHsZ8+cnn\nKb4LRFo3TBABIl2y+3Xvg3N/hp2SqIwRaa6RjNASabX3kS63bw9y475FaESK+EGkHhQMMApj\ncouQ25p+8zciIVJ0YbjXbmbSFCmhAXqASMogUjQgUiIgUsQgkg7fJ5JWmCAQSQdEWjdMEIg0\nM2mukYzQEok1kjKI1IOCAUZhEGlmEKkHBQOMwgSKdL/pO+OmVT9pipTQAD1QEOnMxyg6QKRo\nWEukY+tXcuQhvTDoVQIgUsSEHJHyZ49M77b7znH6PpG0wgRhtUayRaEyy4NI64YJIlCkre3H\nJ3wpIifNNZIRWiIpXGwwJtqZ8Q4i9aBggFGYQJFyN+qXx05JETmI1IOCAUZhAkW6bDfGn+l7\nSxE5aYqU0AA9WPfULvhXUg+liBxEigZESgREihjutdPh+0TSChMEIumASOuGCcJKpNM2tCeD\nKeIkzTWSEVoirfoxij1rpH4QqQcFA4zCBIr08OgY0ou+FJGDSD0oGGAUJlCkzP1WvyfpvHHc\ntOonTZESGqAHa98i9FMejf74RWMdIFI0rC3S0R2sb7pLcpwGQaSICb37uzy1O7u8OCFSON8n\nklaYIAJFOlYCbaqLDfx+pGAQad0wQYRe/v6pftrxa106SXONZISWSKu+jzQP0c6MdxCpBwUD\njMIg0swgUg8KBhiFQaSZSVOkhAbowbq3CGXcItQLIkWDxi1CiBQMIkVM8BuyB7OudKT4Hr5P\nJK0wQRjc2TADCpVZHkRaN0wQwad2fItQP2mukYzQEmnViw2bzTkk/Scp4gaRelAwwChMqEhH\nLjb0g0g9KBhgFCZQpB+u2g2QpkgJDdCDNUXKuGo3ACJFw9qfR5qBJMdpEESKmOBTO67amfF9\nImmFCSL4YxR897cZiLRumCCCT+242NBPmmskI7REWneNhEj9IFIPCgYYheFjFDODSD0oGGAU\nBpFmJk2REhqgB6uJdNnXfz3lLjN+NymhcUKkaFhNpKxeGDU3CZl+P2Sa4zQIIkVMgEgHt6ne\nRMqyv+Kycb8r9yoBvk8krTBBBIi0cdWN3yf3U//JVxYHg0jrhgkiQKTmive++fp8Ln93kOYa\nyQgtkdZaIzXu5O7pByuinRnvIFIPCgYYhQkQKa9O7c7NdxVfXBbSi64UCYBIPSgYYBQmQKR9\ndbFh1/yGsQPf/d1BmiIlNEAP1hLpkt2vex+c+wvpRVeKBECkaFjvDdnbl+c7vkTfAkSKGJNb\nhNzW+LMU3zlO3yeSVpgguNdOB0RaN0wQiDQzaa6RjNASibu/lUGkHhQMMAqDSDODSD0oGGAU\nBpFmJk2REhqgB4ikDCJFAyIlAiJFDCI9Ze1jifxfJ5JWmCAQqYsVZicirRomiDVEOmycy39m\nTWHA8oMckUhGfN7x3rMFo1OJ6NZI1ccvLr1RFGYGIs3PCJEWqE18Il0QybspInVviUjvHPIY\nTu2MSFOkxTMi0udhF75GthCIZJMQkd7JB11JSKQRIFJPQvXarHJEGvwUoIJIXP7u6+rSYeRr\no3Nqt3CKQRCpr6tLh5GvDSJ1wVW7+dESKbo1kkSKQRBpfhBpZhRmBiLNDyLNTLQz4500RUry\nqh0iKYNINgkRaQoJiTQCROpJqF4bROqCy999XV06jHxtEKkLROrr6tJh5GuDSF1w1W5+tERi\njTQLiDQ/iDQzCjMDkeYHkWYm2pnxTpoiJXnVDpGUQSSbhIg0hYREGgEi9SRUrw0idcHl776u\nLh1GvjaI1AUi9XV16TDytUGkLrhqNz9aIrFGmgVEmh9EmhmFmYFI84NIMxPtzHgnTZGSvGqH\nSMqsLZLVF2NP3isbjGoz31eIpydSvF8vucCrrlVpohVpthKnJlJtkYlKSV7+FhNp8TUSIo1r\niUhdKay6unQYRJpE8Lm7QR+TvGondtaLSDODSLajrAoizQwi2Y6yKog0MwprJCPSFImrdm8d\nnLBP45ssl8Luqp0RWiKJXWwYkdCmNog0oq2SRqNApJ6ENrVBpOXh8vf8cVgjzQwi2Y7yI4VV\nV5cOg0iTSEikNNdIRiDSzCiMd5IiiYFIM6MwMxBpfhBpZqKdGe+kKRJX7d46OGGfxjcRTLEU\nWiItf7HB6BNAiDSJhEQaQZoi2eyUVpigioQ0EUwxCJe/TeJoGYBIy4NIJnG0DECk5Unyql20\nrw5aYXz7OaE045sIphgkSZGWR8sARFoeRDJBywBEipY0RWKNNL0iIU0EUyyFlkhcbDAJE1SR\nkCaCKQRBpO6d0goTVJGQJoIpBon2AldvCquufr6lzU5phQmqSEgTwRSDIJJJHC0DEGl5krxq\nN6Y3Kd4kh0jLk6RIYzpu05skw/jKNaHC45sIphgEkUx6k2QYX7kmVHh8E8EUS4FICYbxlWtC\nhcc3EUxhidZyorenI3bKpjdJhgmrbEATwRSWyI/yI0WsO6UVJqyyAU0EUwwS7ZxLcqe0woRV\nNqCJYIpBop1zSe6UVpiwygY0EUwxSLRzLsmd0goTVtmAJoIpBol2zhmhtVNaYXzlmlDh8U0E\nUwyCSEo7pRXGV64JFR7fRDCFJfKjHP9OaYXxlWtChcc3EUxhifwoP1LEulNaYcIqG9BEMIUl\n8qP8SBHrTmmFCatsQBPBFINEO+eS3CmtMGGVDWgimGKQaOdckjulFSassgFNBFMMEu2cS3Kn\ntMKEVTagiWCKQaKdc0Zo7ZRWGF+5JlR4fBPBFIMgktJOaYXxlWtChcc3EUzR5JH6/AMixRLG\nV64JFR7fRDBFk0dqeFgjxRImrLIBTQRTNHmkhgeRYgkTVtmAJoIpmjxSw4NIsYQJq2xAE8EU\nTR6p4UGkWMKEVTagiWCKJo/U8CBSLGHCKhvQRDBFk0dqeLhqF0sYX7kmVHh8E8EUTR6p4UGk\nWML4yjWhwuObCKZo8kgNDyLFEsZXrgkVHt9EMEWTR2p4WCPFEiassgFNBFM0eaSGB5FiCRNW\n2YAmgimaPFLDg0ixhAmrbECTJVO83gjXu63U8CBSLGHCKhvQZLkUtUUfq6Q1PIgUS5iwygY0\nWS6FG9Vea3i4ahdLGF+5JlR4fJPFUriXx6HNpYYHkWIJ4yvXhAqPb7JYCkQyQ2untML4yjWh\nwuObLJYCkT6rz7fWBpHGtWSN1Fue760NIn3ckqt2vTv88ZZiO6UVJqyyAU2WTMH7SH07/PGW\nYjulFSassgFNBFM0eaSGB5FiCRNW2YAmgimaPFLDw1W7WML4yjWhwuObvAY4Ztlx1hSf9kRq\neBApljC+ck2o8PgmrwGyovr/jCk+7YnU8CBSLGF85ZpQ4fFNXgL85h6Rur9/cT60hoc1Uixh\nwiob0KTN+ScvOLWbLUzvDn95bdISSSJFk0dqeBApljBhlQ1oIpiiySM1PIgUS5iwygY0EUzR\n5JEaHkSKJUxYZQOaCKZo8kgNj1WYFH/FhlYY31Qa2sCkiWCKJo/U8BAmljC+qTS0gUkTwRRN\nHqnhIUwsYXxTaWgDkyaCKZo8UsNDmFjC+KbS0AYmTQRTNHmkhocwsYTxTaWhDUyaLJni6z9G\nQZjZw/im0tAGJk2WS8EH+wgzfxjfVBrawKTJcincqPZaw0OYWML4ptLQBiZNFkvhXh6HNpca\nHsLEEsY3lYY2MGmyWApEIswCYXxTaWgDkyaLpUAkwiwQxjeVhjYwabJcCtZIhJk/jG8qDW1g\n0mS5FOM+GKg1PISJJYxvKg1tYNJkuRSIRJj5w/im0tAGJk2WS8GpHWHmD+ObSkMbmDRZLAUX\nGwizQBjfVBrawKTJYilc0dwihEiEmTGMbyoNbWDSZLEU7nqLECIRZsYwvqk0tIFJk+VSsEYi\nzPxhfFNpaAOTJoul4IhEmAXC+KbS0AYmTRZLwRqJMAuE8U2loQ1MmiyWgqt2hFkgjG8qDW1g\n0tnI9YoAAAsWSURBVGS5FKyRCDN/GN9UGtrApMlyKfhgH2HmD+ObSkMbmDRZMgUfNSfM3GF8\nU2loA5MmgimaPFLDIxaG75nsDOObSkMbmDQRTNHkkRoesTBGIqXoo28qDW1g0kQwRZNHanjE\nwiBSZxjfVBrawKSJYIomj9TwpBkGkQybLJnidQh6t5UaHsLEEsY3lYY2MGmyXAoufxNm/jC+\nqTS0gUmT5VLwhixh5g/jm0pDG5g0WSyFe3kc2lxqeAgTSxjfVBrawKTJYikQiTALhPFNpaEN\nTJoslgKRCLNAGN9UGtrApMlyKVgjEWb+ML6pNLSBSZPlUvjeg+jZWmp4CBNLGN9UGtrApMly\nKRCJMH1hbN4e9k2loQ1MmiyXglM7wvSFQaRxDbnYQBhvGEQa1xCRCDNjGN9UGtrApMliKRAp\n1TA2hxJEGteSNVJyYRBpCiFX7W5/fLS10SgTZvYwiDSFkBSvtevd1miUCfNdYXxTaWgDkyaC\nKZo8UsNDmFjC+KbS0AYmTQRTNHmkhocwsYTxTaWhDUyaCKZo8kgND2FiCeObSkMbmDQRTNHk\nkVrCEiaWML6pNLSBSRPBFE0eRCLMhDC+qTS0gUkTwRRNHi2RtHpDmM4wvqk0tIFJE8EUTR6p\n4UGkWML4ptLQBiZNBFM0eaSGB5FiCeObSkMbmDQRTNHkkRoewsQSxjeVhjYwaTJzis9frV8b\nGg0Ph5LuMCnWxjeVhjYwabJoihGN05wsUr0xGipEmgQiJdMbo6FCpEnEJpJWGESauzZBuxLS\nZNEUXy+SVBgjtHrj6+AiTQRTNHlshocw3WG0Rqo/R1AHF2kimKLJYzM8hOkOozVS/TmCOrhI\nE8EUTR6b4SFMd5j++i89Uka9sWmsLhJrJKUwRkOl1RubxoikPXW1whgNlVZvbBojkvbU1Qpj\nNFRaKzZfBxdpsmgKRFIKYzRUiDQJREomjBFavfF1cJEmgimaPDbDQ5juMFoj1Z8jqIOLNBFM\n0eSxGR7CdIfRGqn+HEEdXKSJYIomj83wEKY7TH/9lx4po97YNFYXiTWSTZjZ7u+cMlSINAlE\nWj+MVY1twsjvVKQiGb1aak1drTBhA2QdRmvF5uvgIk3MU0ytKyJ9HGYJEGlmFhPJ5simZYD8\nnJuC0UjN2MFFmpinMBIpqA/hvUEkKb5yjTT/XBmzaYoHNqvixAMiKc+VaA9sVnu1QBgjEGmW\nuWLECiItsJzQEkkhDCLpsMChxGr/FaauVphYRZr/RXd5EGndMEEg0uQ+mDN5p5YXyQhEmpmE\nREoyjBgKexWrSLOfBDF1v46vXCMh0qokuVeIZCSS0TnhyF2JkiRfZb5SpPnXSNM7L+wjV+3m\nChOpSLM1VkVhrtjHWTuM4WveGiIdNs7lP7OlQKQZwsxznF97rwxZQyR3LopL/2olKHxIY1US\nmXLSPgaxikgXn0ghh9U5RmdtjObcPFN3bfR2apVTu3zWUzuA5UnwYgPA8iASgAGIBGAAIgEY\ngEgABiASgAGIBGAAIgEYgEgABiASgAGIBGAAIgEYgEgABiASgAGIBGAAIgEYgEgABiASgAGI\nBGAAIgEYgEgABoiKBBAZE2a5vTimaH2XYpJhxLqjFUY231i0yppkGLHuaIWRzTcWrbImGUas\nO1phZPONRausSYYR645WGNl8Y9Eqa5JhxLqjFUY231i0yppkGLHuaIWRzTcWrbImGUasO1ph\nZPONRausSYYR645WGNl8Y9Eqa5JhxLqjFUY231i0yppkGLHuaIWRzTcWrbImGUasO1phZPON\nRausSYYR645WGNl8AEmCSAAGIBKAAYgEYAAiARiASAAGIBKAAYgEYAAiARiASAAGIBKAAYgE\nYAAiARiASAAGIBKAAYgEYAAiARiASAAGIBKAAYgEYAAiARiASAAGIBKAAYgEYEAMIv2t3QE7\nEtoVTf56fpoVTZGc2zaPb39+ytbl52JCu2v650au9Rf3cVTn2aoV17l913Md8YbCj8H9Ve23\nVaWzS3HJiqM7jgtQ9797cMb07pi7/G90K0+IdvMFZ7eoSEU5vG5TuGM5zuWAZdlpXE93h+KU\nV2Nd/bfPjtUkuT1WgYfSV39sqo2P2d7dGrm/PDu6v23dqTJYE35TbdIRxRU/zu3c9mfTpL/t\nyu/9+VP5c7VNk2fb36vNuThvysZVP6ogH3SiJ1xeRizrdNgV+9/id1/strtxAZr/3ztT7+Xp\nun8j++POxV9W1aGoq93EGkcdwlXsmpr8Fq1yz4yoSO6n/ONSTpDqNbt8pfnLx/X0dtyoRSpb\n55cywPWxCjyUvmp3qTYuU7tbo7If5+plr+qUq7vUbFZu0hGlMdKd3cFdW+S3dtdX9PrnZheb\nne3bp+O+2JcyFlU/qiAfdKInXH50dZ3KuLtidyyPS9m4AM0R6d6Zei9v+zeyP4ftoT6BKMPU\nQ1XHGkcTotmhZtyLVrlnRlSk28vduTpBq0syXqRbS3ed07fH4UjXc7jXRvWAuOxYdcrdx6vz\nHKs54mzzOoprtbgmOBzrf222ue5sb6+ySzNn6x/dJ53oCVeeC13nXZGV//srd3PUmqLu//1F\noWjv5ej+/B12VR1uQzVyvG8hmmI+jfuj3DOjKlLWTMTsJ5t0RCpPWY6bpyNS9Vp7e/xYpGrj\n6pX11qg5Iv3uq049DgadL77NIJ5Ot4nxdkQqo9ev4NU2dZ6f3kNCuf3Pdl94j0hjjwB1uPx2\naldsD5vipzy/+xkXoMifj0j1/yYekbK/4pxVdSjqak8RqQ5Rpq6DcERqaE7tqvPcanmyGb9G\nujQXG37rgu6y07E8C7k9fixStfGpXrs0jVy5ADiWZ1G/zYLr1ISvN+mIUs5U93ObZGX6267c\nRNofqzVFtU0V5Ld/tV++VF/KhYCrF1b12uiDTvTt5N/tYkMZ41Dkl6I6hR0ToOx/U5SHSKfr\n/o3sT/lqkh+rOhR1taeIVIcoC1wF+b2K9FzumdEUyRT38mjP3yeLC5P0PYk+6gTMBSIFsy2X\nTSO6EUTndbUPOwFz8QUiAcwPIgEYgEgABiASgAGIBGAAIgEYgEgABiASgAGIBGAAIgEYgEgA\nBiASgAGIBGAAIgEYgEgABiASgAGIBGAAIgEYgEgABiASgAGIBGAAIgEYgEgABiASgAGIBGAA\nIgEYgEgABiASgAGIBGAAIgEYgEgABiCSNPVv6Xb5/tL80Hru2Nqu8ylYBESSxl3JzsWrLblr\nbdf5FCwCJZemEeS8cZuOp3pawZJQcmluSuTvv6cZkaSg5NLclDi63fWH48a5zfF60lc9XHK3\nvW63z9zmfH8KloSCS3MT4uLy5odDs2Y6PETaOrdvttvUq6kLIq0BBZfmLkQjTVFk7q8ofm9a\nVX9uLtenf6u/7kqrOLVbAUouzZtI7r5Yuol0uv2wrf56cRkirQEll+ZNpH15Lvf393jqusHt\n6ZdWsBiUXJqbEufq+nf9w092fVsJkaSg5NLclPh9Wvoc9/nzGum+HSKtCSWX5vE+0unVlDeR\nNqyRVoSSS9O6s6H6IXe/96t2T7cNNVfGN5dyDbW/PQVLgkjSvN1r99v8fKqUejr4tN5Huj4F\nS4JI0jTabH6uPxTXOxuqS96n/FWk+pLe+f4ULAkiARiASAAGIBKAAYgEYAAiARiASAAGIBKA\nAYgEYAAiARiASAAGIBKAAYgEYAAiARiASAAGIBKAAYgEYAAiARiASAAGIBKAAYgEYAAiARiA\nSAAGIBKAAYgEYAAiARiASAAGIBKAAYgEYAAiARiASAAGIBKAAf8Bd6yQEvzmPvEAAAAASUVO\nRK5CYII=",
      "text/plain": [
       "Plot with title \"Fig. 3: Smart Scale Score Distribution by District\""
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "boxplot(SMART.SCALE.Score~District,data=data, main=\"Fig. 3: Smart Scale Score Distribution by District\",\n",
    "        xlab=\"District\", ylab=\"Smart Scale Score\", ylim=c(0, 30), col=\"yellow\", cex.axis=0.4) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In Fig. 3 the Hampton Roads district has significantly more variance than all other districts, while NOVA and Richmond have the highest medians scores.\n",
    "\n",
    "### 3.2 Descriptive Statistics for Minority High-Scoring Projects (SSS > 30)\n",
    "\n",
    "Up to now the exploration has been focused on the majority of projects which tend toward the lower range of SSS. We turn our attention to the scores falling on the far right of Fig. 1 ranging from 30 and above. We first tabulate the number of projects falling in this range by area type:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Area_Type</th><th scope=col>Num_Projects</th><th scope=col>Mean_Smart_Scale_Score</th><th scope=col>Median_Smart_Scale_Score</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>A        </td><td>13       </td><td>105.97615</td><td>61.1800  </td></tr>\n",
       "\t<tr><td>B        </td><td> 4       </td><td> 58.23500</td><td>47.1480  </td></tr>\n",
       "\t<tr><td>C        </td><td> 6       </td><td> 47.45583</td><td>41.3365  </td></tr>\n",
       "\t<tr><td>D        </td><td>10       </td><td> 59.80840</td><td>62.2070  </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|llll}\n",
       " Area\\_Type & Num\\_Projects & Mean\\_Smart\\_Scale\\_Score & Median\\_Smart\\_Scale\\_Score\\\\\n",
       "\\hline\n",
       "\t A         & 13        & 105.97615 & 61.1800  \\\\\n",
       "\t B         &  4        &  58.23500 & 47.1480  \\\\\n",
       "\t C         &  6        &  47.45583 & 41.3365  \\\\\n",
       "\t D         & 10        &  59.80840 & 62.2070  \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "Area_Type | Num_Projects | Mean_Smart_Scale_Score | Median_Smart_Scale_Score | \n",
       "|---|---|---|---|\n",
       "| A         | 13        | 105.97615 | 61.1800   | \n",
       "| B         |  4        |  58.23500 | 47.1480   | \n",
       "| C         |  6        |  47.45583 | 41.3365   | \n",
       "| D         | 10        |  59.80840 | 62.2070   | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "  Area_Type Num_Projects Mean_Smart_Scale_Score Median_Smart_Scale_Score\n",
       "1 A         13           105.97615              61.1800                 \n",
       "2 B          4            58.23500              47.1480                 \n",
       "3 C          6            47.45583              41.3365                 \n",
       "4 D         10            59.80840              62.2070                 "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "high.ss <- subset(data, SMART.SCALE.Score >= 30)\n",
    "dt3 <- data.frame(table(high.ss$Area.Type))\n",
    "m3 <- data.frame(aggregate(high.ss$SMART.SCALE.Score ~ high.ss$Area.Type, FUN=mean))\n",
    "mm3 <- data.frame(aggregate(high.ss$SMART.SCALE.Score ~ high.ss$Area.Type, FUN=median))\n",
    "dt3 <- cbind(dt3, m3$high.ss.SMART.SCALE.Score, mm3$high.ss.SMART.SCALE.Score)\n",
    "colnames(dt3) <- c('Area_Type', 'Num_Projects', 'Mean_Smart_Scale_Score', 'Median_Smart_Scale_Score')\n",
    "dt3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see above that area types A and D have the majority of these high-scoring projects. A similar breakdown by district is shown below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>District</th><th scope=col>Num_Projects</th><th scope=col>Mean_Smart_Scale_Score</th><th scope=col>Median_Smart_Scale_Score</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>Bristol       </td><td> 0            </td><td>       NA     </td><td>     NA       </td></tr>\n",
       "\t<tr><td>Culpeper      </td><td> 0            </td><td>       NA     </td><td>     NA       </td></tr>\n",
       "\t<tr><td>Fredericksburg</td><td> 2            </td><td> 57.72750     </td><td>57.7275       </td></tr>\n",
       "\t<tr><td>Hampton Roads </td><td>11            </td><td> 95.42873     </td><td>59.8350       </td></tr>\n",
       "\t<tr><td>Lynchburg     </td><td> 3            </td><td> 67.47600     </td><td>79.6260       </td></tr>\n",
       "\t<tr><td>NOVA          </td><td> 4            </td><td>108.68375     </td><td>86.2600       </td></tr>\n",
       "\t<tr><td>Richmond      </td><td> 4            </td><td> 50.86100     </td><td>40.4150       </td></tr>\n",
       "\t<tr><td>Salem         </td><td> 4            </td><td> 58.85550     </td><td>48.1505       </td></tr>\n",
       "\t<tr><td>Staunton      </td><td> 5            </td><td> 50.44980     </td><td>42.0710       </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|llll}\n",
       " District & Num\\_Projects & Mean\\_Smart\\_Scale\\_Score & Median\\_Smart\\_Scale\\_Score\\\\\n",
       "\\hline\n",
       "\t Bristol        &  0             &        NA      &      NA       \\\\\n",
       "\t Culpeper       &  0             &        NA      &      NA       \\\\\n",
       "\t Fredericksburg &  2             &  57.72750      & 57.7275       \\\\\n",
       "\t Hampton Roads  & 11             &  95.42873      & 59.8350       \\\\\n",
       "\t Lynchburg      &  3             &  67.47600      & 79.6260       \\\\\n",
       "\t NOVA           &  4             & 108.68375      & 86.2600       \\\\\n",
       "\t Richmond       &  4             &  50.86100      & 40.4150       \\\\\n",
       "\t Salem          &  4             &  58.85550      & 48.1505       \\\\\n",
       "\t Staunton       &  5             &  50.44980      & 42.0710       \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "District | Num_Projects | Mean_Smart_Scale_Score | Median_Smart_Scale_Score | \n",
       "|---|---|---|---|---|---|---|---|---|\n",
       "| Bristol        |  0             |        NA      |      NA        | \n",
       "| Culpeper       |  0             |        NA      |      NA        | \n",
       "| Fredericksburg |  2             |  57.72750      | 57.7275        | \n",
       "| Hampton Roads  | 11             |  95.42873      | 59.8350        | \n",
       "| Lynchburg      |  3             |  67.47600      | 79.6260        | \n",
       "| NOVA           |  4             | 108.68375      | 86.2600        | \n",
       "| Richmond       |  4             |  50.86100      | 40.4150        | \n",
       "| Salem          |  4             |  58.85550      | 48.1505        | \n",
       "| Staunton       |  5             |  50.44980      | 42.0710        | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "  District       Num_Projects Mean_Smart_Scale_Score Median_Smart_Scale_Score\n",
       "1 Bristol         0                  NA                   NA                 \n",
       "2 Culpeper        0                  NA                   NA                 \n",
       "3 Fredericksburg  2            57.72750              57.7275                 \n",
       "4 Hampton Roads  11            95.42873              59.8350                 \n",
       "5 Lynchburg       3            67.47600              79.6260                 \n",
       "6 NOVA            4           108.68375              86.2600                 \n",
       "7 Richmond        4            50.86100              40.4150                 \n",
       "8 Salem           4            58.85550              48.1505                 \n",
       "9 Staunton        5            50.44980              42.0710                 "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "dt4 <- data.frame(District=unique(data$District))\n",
    "dt4 <- merge(dt4, data.frame(table(high.ss$District)), by.x='District', by.y='Var1', all.x=TRUE)\n",
    "m4 <- data.frame(aggregate(high.ss$SMART.SCALE.Score ~ high.ss$District, FUN=mean))\n",
    "dt4 <- merge(dt4, m4, by.x='District', by.y='high.ss.District', all.x=TRUE)\n",
    "mm4 <- data.frame(aggregate(high.ss$SMART.SCALE.Score ~ high.ss$District, FUN=median))\n",
    "dt4 <- merge(dt4, mm4, by.x='District', by.y='high.ss.District', all.x=TRUE)\n",
    "colnames(dt4) <- c('District', 'Num_Projects', 'Mean_Smart_Scale_Score', 'Median_Smart_Scale_Score')\n",
    "dt4"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Because there are so few projects at this score range, we will omit further drill-down to look at actual score distributions.\n",
    "\n",
    "## IV. Examination of Scoring Consistency\n",
    "\n",
    "An important aspect of this analysis is to see if the scoring methodology as specified in the technical guide was applied consistently to the projects. As previously discussed, in the technical guide (Table 4.2, p. 36) there are six major factors used to compute the Smart Scale score which are themselves calculated from more basic measurements:\n",
    "\n",
    "* Congestion Mitigation\n",
    "* Economic Development\n",
    "* Accessibility\n",
    "* Safety\n",
    "* Environmental Quality\n",
    "* Land Use\n",
    "\n",
    "Of these six factors only Land Use is directly given in the dataset. All other major metrics are calculated from their basic constituent measurements as specified on pages 27-31. Within the data there are four area types: A, B, C, and D. Each of these area types has a unique factor weighting that results in its Smart Scale Score. These weightings are specified for each area in Table 4.2 on page 36 within the technical guide. It is noted that for areas C and D, the Land Use factor is not used (or equivalently, carries 0% weight).\n",
    "\n",
    "The examination of scoring consistency will be accomplished by independently calculating the SSS score of each project according to the specifications in the technical guide. These independently-calculated scores will then be compared to the scores provided in the dataset by fitting a linear regression, and the degree of fit will be used to determine the degree of consistency.\n",
    "\n",
    "### 4.1 Comparing Assigned Smart Scale Scores to Scores Calculated According to the Technical Guide\n",
    "\n",
    "Pages 37-40 of the technical guide specify how to calculate the final smart scale score, and in this section we implement these calculations in an attempt to replicate the assigned Smart Scale Scores. In the dataset all basic measurements needed for each factor weight calculation are present except for those pertaining to the Land Use factor. Instead, a Land Use factor score is provided. All steps in the calculation process are summarized in Table 4.6 on p. 40 of the technical guide.\n",
    "\n",
    "#### 4.1.1 Replicating the Technical Guide Calculations on the Data\n",
    "\n",
    "The Smart Scale Score calculations specified technical guide are replicated on the dataset below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "#------ Step 1: Normalization of Measure Weights\n",
    "measure.cols <- c('Throughput.Score', 'Delay.Score',\n",
    "                  'Econ.Dev.Support.Score', 'Intermodal.Access.Score', 'Travel.Time.Reliability.Score',\n",
    "                  'Access.to.Jobs', 'Disadvantaged.Access.to.Jobs', 'Multimodal.Access.Score',\n",
    "                  'Crash.Frequency.Score', 'Crash.Rate.Score',\n",
    "                  'Air.Quality.Score', 'Enviro.Impact.Score')\n",
    "for(i in measure.cols) {data[[i]] <- 100 * data[[i]] / max(data[[i]])}\n",
    "\n",
    "#------ Step 2: Apply Measure Weights\n",
    "# For a given data frame, a set of weights, and a set of selected columns, this function returns a new\n",
    "# vector corresponding to the weighted averages of the columns.\n",
    "weighted.column.average <- function(data.frame, weight.vec, selected.columns) {\n",
    "  v <- 0 \n",
    "  for(i in 1:length(weight.vec)) {\n",
    "    v <- v + (weight.vec[i] * data.frame[[selected.columns[i]]])\n",
    "  }\n",
    "  v\n",
    "}\n",
    "\n",
    "# Calculate the composite factor scores (except Land Use, which is already given in the data) from \n",
    "# component scores according to the Smart Scale November 2017 technical guide, p. 40.\n",
    "data$Congestion.Score <- weighted.column.average(data, c(0.5, 0.5), c('Throughput.Score', 'Delay.Score'))\n",
    "data$Economic.Score <- weighted.column.average(data, c(0.6, 0.2, 0.2), \n",
    "                                               c('Econ.Dev.Support.Score', 'Intermodal.Access.Score', \n",
    "                                                 'Travel.Time.Reliability.Score'))\n",
    "data$Accessibility.Score <- weighted.column.average(data, c(0.6, 0.2, 0.2), \n",
    "                                                    c('Access.to.Jobs', 'Disadvantaged.Access.to.Jobs', \n",
    "                                                      'Multimodal.Access.Score'))\n",
    "data$Safety.Score <- weighted.column.average(data, c(0.5, 0.5), c('Crash.Frequency.Score', 'Crash.Rate.Score'))\n",
    "data$Environmental.Score <- weighted.column.average(data, c(0.5, 0.5), c('Air.Quality.Score', 'Enviro.Impact.Score'))\n",
    "# Land use score already in data - no component scores listed.\n",
    "\n",
    "#------ Step 3: Apply Factor Weights\n",
    "# Given a dataframe with calculated factor scores, this function computes the project value according to \n",
    "# Table 4.2, p.36 of the technical guide.\n",
    "calc.proj.value <- function(dataset, weight.a=c(0.45, 0.05, 0.15, 0.05, 0.1, 0.2),\n",
    "                                     weight.b=c(0.15, 0.2, 0.25, 0.2, 0.1, 0.1),\n",
    "                                     weight.c=c(0.15, 0.25, 0.25, 0.25, 0.1),\n",
    "                                     weight.d=c(0.10, 0.35, 0.15, 0.3, 0.1)) {\n",
    "  factor.cols.ab <- c('Congestion.Score', 'Economic.Score', 'Accessibility.Score', 'Safety.Score',\n",
    "                      'Environmental.Score', 'Land.Use.Score')\n",
    "  factor.cols.cd <- c('Congestion.Score', 'Economic.Score', 'Accessibility.Score', 'Safety.Score',\n",
    "                      'Environmental.Score')\n",
    "  sa <- weighted.column.average(dataset, weight.a, factor.cols.ab)\n",
    "  sb <- weighted.column.average(dataset, weight.b, factor.cols.ab)\n",
    "  sc <- weighted.column.average(dataset, weight.c, factor.cols.cd)\n",
    "  sd <- weighted.column.average(dataset, weight.d, factor.cols.cd)\n",
    "  scores <- c()\n",
    "  for(i in 1:nrow(dataset)) {\n",
    "    if(dataset[i, 'Area.Type'] == 'A') {scores[i] <- sa[i]}\n",
    "    else if(dataset[i, 'Area.Type'] == 'B') {scores[i] <- sb[i]}\n",
    "    else if(dataset[i, 'Area.Type'] == 'C') {scores[i] <- sc[i]}\n",
    "    else {scores[i] <- sd[i]}\n",
    "  }\n",
    "  scores\n",
    "}\n",
    "\n",
    "# Add a new column to data with our calculated project benefit scores.\n",
    "data['Calc.Project.Benefit.Score'] <- calc.proj.value(data)\n",
    "\n",
    "#------ Step 4: Divide Project Benefit by Smart Scale Funded Cost. \n",
    "data['Calc.SMART.SCALE.Score'] <- 1e7 * data[['Calc.Project.Benefit.Score']] / data[['SMART.SCALE.Request']] \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 4.1.2 Comparing the Calculated versus Assigned Smart Scale Scores\n",
    "\n",
    "In order to assess the whether the calculation methodology as specified in the technical guide was applied consistently across the projects, we will use a linear regression model that fits the assigned SSS to the calculated SSS. If the two scores are perfectly consistent across the projects then they will fall perfectly along the line of best fit. Since the factor weightings differ across area type, we assess the fit for each of the four area types below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0gAAANICAMAAADKOT/pAAAAPFBMVEUAAAAzMzMzZv9NTU1o\naGh8fHyMjIyampqnp6eysrK9vb3Hx8fKysrQ0NDZ2dnh4eHp6enr6+vw8PD///8C09X6AAAA\nCXBIWXMAABJ0AAASdAHeZh94AAAgAElEQVR4nO2djXarrBKG+WpMbHvSZMf7v9cj/oKCMjDA\noDNr7TY7HWCc1wdREUXLxsYWbCJ3AGxsZzAGiY0NwRgkNjYEY5DY2BCMQWJjQzAGiY0NwRgk\nNjYEY5DY2BCMQWJjQzAGiY0NwVBAEpPJj3a3Wiz+xkpAbblGvvZ8PuY67n/L15+feyXqH8dK\nDa3s/LGRf3w8fetGMZBGo+vjb/3nKCKZWkut0l+tNfSqRdV8zJ/MQiYE6Wf627Pz3ARj/HKv\nLcfI1/X+VWodc+pe1fBF9XGrdhuW/W/N8Mdqs1umNJBGs/NKkTgiGVpLrdJnqPc1tTO08jZ/\nMgqJBJKD02tO60PcxWP9Z+OXIe1Z6q1+5vKfRtymr2/i0YnzrkXjXrVjSI9xw38qv6pxDKTR\n6Py95GewOCIZWkut0o+s8SG+x/8+xM/83faTUUhskPqP3Zbenuugb9X0jRCf8aNW8DOJ+Krq\nLn/d/vfoO56/e9cRaJkbHD9DguWv7ou7qN/9f+dy68Z6ayq12eUvk2r97+4gfm+Xvrv/OUex\nie89eY8Vrfriqvob/1N5yo9iII0WpdYlY4iktdZbcpXu8mD06jzUdkRt/mQUMgZIn8pwWP8W\nv+M3z65TewyH8cVn+VKIWvZ6fR23dhhOdKbGPhZrev/frh/psjUd8Ody63p7+4hmKa/2dXdl\nXNHH3+gSLVGs4+u973aQGiVezxEJhoE0mo8RTZtAJK01aelVqrT+Q6nR9MkoJOo50tjQd4ft\np9ZFei39gZTiuR4hLF8OcvRZbeTR9CZ++zGH1t5Yp+we7t0Qukvbp+0P+Eu5db29fY+7wWb0\n/a7Erfntu8uu+BT/krwlinV8jWFrN/npfz3noUMGA2k0OacRadNaepX0A3EX5c94amv6ZBQy\nBkg30W3qWw+66noidTy33qTlSyFLd3X03/bH2vfzuzZpNByQRX9e+pIt3vRyhsbufeWzcq+l\nzs/3TXZef+0Y/2vVC81RrOMzbe0mP/2v9xJWegNpNDnf13W0MUTatJZepRVIf7Ldvh7TJ6OQ\nMYZ2q6ik9QOC4Zvx8Ktf/VG+1KVsZVewHoRM/5FdaN9DqgdgzXvdmBLas+8qW8X71Txq2aUZ\nD+dzFOv4VlsrVhGof9ylLbKBNBq//9VP6mOJtGktvUrrdDxv4j6ciZk+mYRMBJIS+GPqaNQa\nlC/XKXiI28/zbdRI9jONLv8qQevG1NCaXqOV/kPXud2SJYrTgiTW2zRkY7FYIm1aS6/S+Gdt\nc98Dw6ZPiUAyHEaVwKv+TO2jB618OZS7CbXG9mPW6Cma/jRRDC3Warl1vUrJsfxNGeqK8fxR\nj793XD59FkmWdtyHdpRA2tfIHHAskTatpVfppl+16+1vPkJuPiUCqe9EDCd2/Td/49HhIZTb\nWuqXU0fUyIN9f9nxb32aKJT8jB1Wfyr5rZZb19vbNPru//Na7sDJmP/6K0T3/kT8M8RfdWOI\n6YT2Tzu1XdpZvI/yk/kcSfu4q9H4c8jGbLFE2raWXqVGbttjvvxx646gn/7kzvQp2TmScml1\nm9pm3KefYrkurH05pnWo49VfKLENlfph+2//hRwcV1q5db29TdeDpv8t18lvQyvyznV/0bQe\ndrfOvudPSm+9tLO6sLqTn8xX7bSPuxotB6dXG18kvTVp6VUapi30xzvp89P/T+4zpk8xr9rp\nH+XNvl+LSNU0yqqq5c/ql1P33Y3I6165/oNFo2G8IH/W4vHWy60bk6bdodCHDT+1mGZVveWN\nwwH7qr8Bokaxia/zvrcuIFG6j7Sr0bhrV8P1stgi6a1Jy6CSnGs3OPY+P7fubGto0PAp2n0k\no4kUM2L+ht7K/eSjcfbEP6HxndgSz5JoBBcpq0oujcaa2bCuU44BGudJWSFWr2dIHFrlOnUY\nXaK8c+1WllAjD5EyquRg8ebarWwcqb4jVK2bEOMJKyCdf667M7pEeWd/ryyZRl4iZVTJweLN\n/l5bN5ocB8NxrZqunkDS+XTshbElyvw80tpSaeQnUjaVHCzi80hsbFc3BomNDcEYJDY2BGOQ\n2NgQjEFiY0MwBomNDcEYJDY2BGOQ2NgQjEFiY0MwDJD+rWzzxbF5FElVJlszCMpYNXINENsv\nW8ORAlwSzCCRbQZBGatGheyn6f0YpGhlGCREPwYJIBLhvZVwaAwSCT8GKVoZBgnRj0ECiER4\nbyUcGoNEwo9BilaGQUL0Y5AAIhHeWwmHxiCR8GOQopVhkBD9GCSASIT3VsKhMUgk/BikaGUY\nJEQ/BgkgEuG9lXBoDBIJPwYpWhkGCdGvdJCEEJrfkmAGiUwzi0jjN5jmGSCDpFm/FqDqtyTY\nDaSqGhZ+Xv82ikR4byUcWquINH6DaZ4BMkiqjUv0K35Lgp1AqsYf699mkQjvrYRDE19fK5IQ\nMWKQGKTsZRI18/XVk6RWgmmeATJIqgWD1BuDFLPMyBGDlNkv8jmSNDNI/3WGKfhFbeIoqJJB\nEafzWJr7KQG/2FftBlH4iBSpjOQo+KpdtchyqBHR/TS/X4L7SDy0i1ZGcrQtAuWoZZDC/Rik\naGXiN9Mfj4JBahkkBL8Dx/+t/IAg8VW7mGV6jmKCxOexWNYNHGzGIOVuZuCIj0g0/HYdR6kU\nPyBIPLMhXplJHAaJhN+eo8pRwDmS3byDDipyDpBmcRgkEn47jhpHDBKtZhZxGCQSfnZHnSMG\niVQzijgMEgk/q+OKIwaJUjPbk1etEqjxzIZwP5vjmiMGiVAzhpNXrRJM89woBmmwDUcMEp1m\nTCevWiWY5rlRDFJvW44YJDLNGE9etUowzXOjGCRpBo4YJCrNmE9etW8wzXOjGKRBqQ1HDBKN\nZrbiMEgk/LaOZo4YJBLNGMRhkEj4bRwtHDFIFJoxicMgkfBbO1owYpDyNyPEzphb+wbTPDfq\n6iBZOWKQcjdj4YhBouGnO9o5YpAyNzNwJLZ/YJBI+GmOOxwxSJmbGThikKj6qY57HDFImZvZ\nLl9nLYKgjFUjCvspST/FcZcjBilvM1aOGCQafovjPkcMUtZm7BwxSDT8ZscDjhiknM0cn7xq\n32Ca50ZdFSTbbdi135JgBildMw4nr9o3mOa5URcF6ZgjBilfMy4nr9o3mOa5UdcEyYEjBilb\nM04nr9o3mOa5UZcEyQEjBilbM24nr9o3mOa5UVcEyYkjBilPM64nr9o3mOa5URcEyY0jBilL\nM84nr9o3mOa5UdcDyZEjBilHM+4nr9o3mOa5UZcDyZUjBilDM4CTV+0bTPPcqKuB5MzRMUg/\ndyHa+hUgEoOkl4GcvGrfWC1cIwbJaO4cHYH0ufUvYBTij0FCagZ08qp9YzEMjRgkk0mpYBVa\nQXqIplOo/RU1g4TTDOzkVfvGYhgaMUhbG85kkUCSrwOe/jFICM0AT161byyGoRGDtLHxihCD\nFKtMUDPQk1ftGwYpod90ZRV3aNeIB4OE0Az45FX7xmIYGjFIK5uVQgLpU4neqjeDFN4MgCPI\nxQYEjRgk3RalkEBq2++bELfmA9CIQbKUgXAEuvztoRHbnvWX68KMb8hGawbEEQikcI34iKSa\nqhTSEamGjLstIjFIfRkYRwCQMDRikBTTlEICqfI5QvluXFgR4iABOQKAhKERg7SYrhQSSK+6\ngZzCGkVikP65PtSy24wt3RgaMUizrZRCAknM5i8Sg+T8UMtuM7Z0Y2jEIE22VopBilXGpxk4\nRwxSHr+NUmiXvz3Md+PCihAGyYMjvmqXxW+rVE6Q2HRDuCuBbrj7n7cfMZAMPR4WSJ+Gb8gG\nFvE5HoGOSAgaMUjSTEohgfTmKUKhRSAPtew2Y0s3hkYM0v5Lll0rtIL0EHUnz7vmSau+RUAP\ntew2Y0s3hkYMkm3kgHbVTv/NIAGLwB5q2W3Glm4MjRgk2wicQYpVBlQE+FDLbjO2dDNIGH62\nM1ke2sUqAykyqcNDu1QV+vpZrwjxxYZYZQBFwE+H7TZjSzdfbAj3s19Z5cvfscq4F4E/Hbbb\njDXffPk71G/nDgUWSD7mG0tYEXIgKerEBQlBo2uDtHenj0GKVca1iM/TYbvNIChj1ejSIO3e\nMccC6d5/IW58jgQs4vV02G4z1nwjaHRlkPZnniCB1AzXVAVftQMW8Xs6bLcZW7oxNLowSAcz\nuJBAqoZ1cF98HwlWxPPpsN0itnRjaHRZkA5fDIIEEt+Q9Sri+3TYbhFbuvmGrL/fIUdYIN3F\n4yOvr/La35Ai3k+H7RaxpRtDo4uCdMwRFkjzzT7IO0N8YwkrQgck/6fDdovY0o2h0TVBcnnA\nBeuq3XizD7S6hm8sYUXIgBTwdNhuEWu+ETS6JEhOD4phgeRjvrGEFaECUsjTYbtFEJSxanRF\nkNweuGSQYpXZLxL2dNhuEQRlrBpdECTHB5dRQOpOYLufv5W489DOrYjl7DUiSDgaXQ8k1wUA\nUECq5BXVv/5EFjIj0jeWsCIUQLJdBYoIEo5GlwPJeSENDJB+RN1pc6vlzfOGQTouEvxQy24R\nY66RNLoaSO4L0mCAVAv5wJicefIRFYN0WCT8oZbdIsZcI2l0MZAACzthgNTfKv/tOzqe2XBc\nBOGhlt0ixlwfalT1Nn1gkKRBFkjDAKl/zUEjXgySSxGMh1p2ixhz7aZRNf6zanQpkEALpGGA\ndBdy+H1r5cksTxE6KILyUMtuEWOunTSq5h8MEnSBNJyLDY/2Kb674Xctfhik3SI4D7XsFjHm\n2kmj7QHpwiABF0jDAKl/y6+8qCrEDcDRFUFCeqhlt4gx1y4aDQek5RTpv84gep7JvlItxK6N\ns1+34TYf6MLqFUHCeqhlt4g52Q4aVfoPk0ZXOSLNN/qSHpF8zTeWsCLZQBIC7aGW3SLeelSG\nT54BFg7ScsOcQcpeZl3EgaO8IFWmj54Blg2Sz7pODFKsMqsiA0cidjP/QkHioZ3fuk4MUqwy\na5B6jsoASTk2eQZYMkh+6zoxSLHKrIoMHNEHqdUmNlwQJM91nRikWGX0Ik4cZb7Y4KDR+UHy\nXdeJQYpVRivixhGDlKxCm5/3uk5oIP3chWhryLoaFwIJ++mw3SL2hIdrdHaQ/Nd1QgLpc5OP\njLViWIOQQdKKOCzmhNHM/I3FMDQ6OUgB6zohgfQQjZxV/MuTVrdF3DmKCxKGRqcGySRUapDk\n1PzpH4OkFQFwFBckDI3ODJJRKAYpe5mxCOThMAYpWYVbv8B1nXCHdg2/jWJVBMRRkqFdkEbn\nBSl0ORqsiw38DlljERhHkS828Dtk7X42odJf/v7md8huiwA5in35O1yjs4IUvq4TGkge5h10\nUJGEIEE5ig1SuEYnBQlhXae4IE3ztta/jSKdDiQwRwxSsgo1P4x1nTBAEqppHI0/1r+vARKc\no4ggWTVikHDWdWKQYpXx4IhBSlah4oezrlP8c6SLguTDEQ/tklW4+CGt65QLpLOvUJNqEZrY\nFrS74PlFa/ho3klykBrbsGF+4vJSRySv41HsI5JVo+uCdDh/KzVIjXX8fUmQQIvd+jdjLnLI\nEYM0+R3Pg0wNUiVetXh/6s0UfW2ttKuABFvs1rsZSxEbE1aNLguSw7ghNUhdL/ctnu1nPUVf\nX5nmIiABF7v1bcZWxMaETaPLguQy/s4A0lOuKb0aNlTKr8uANOpDESSTRpcFyek8NjVId/H7\nFjf5akWNI/WdO63y+8wggRe79WvGXsTGhFmjy4Lkdj0oNUhSnVqex179MQqPxW59mtkrYks3\nhkbnAcnxumryy9/Pm3zgBbaKvnfQQUWillEuA5EDCUOj04CEvhwNFkg+5h10UJGYZdTLqfRA\nQtDoLCDhL0fDICGW8VvsFtzMQREEZawanQSkCMvRYIF0778Qtys/Ieu52C20maMi1nwjaHQO\nkGIsR4MEUjNcCrr0xQbfxW7DyriDhKHRKUCCzDtJDVI13C5/XXgVIe/FbsPKuIOEodEZQALN\nO0kN0iTOdUHyX+w2rIw7SBganQAk2LyT1CDdxePTtp/msiutBix2G1bGHSQvjc5mX9Seb1mB\n9J6WeoKs0B5lN8pTJmSx27Ay7iBhaFT6EWm+P0H1iNR1dP1ST5ALQicCKWix27AygKt2CBoV\nDhJ83klykHzMO+igIhHKhC12G1YGABKCRmWD5HG7nEFKWMZ8W4JBQvRDqVDViSZI3Qls9/O3\nEvdLDu0CV40OK+MKEo5GJYPkNe8kLUiVvKL615/IQtbD9Q46qAh6mdBVo8PKuIKEo1HBIPnN\nO0kK0o+oO21utbx5Dpla7B10UBHsMsHLr4eVcQQJSaNyQfKcd5IUpFp0o4W3nHnyEdVGipOD\nFL78elgZR5CQNCoWJN95J0lB6m+V//Yd3eVmNiAsvx5WxhEkJI1KBcl73klSkCr5n0a8LggS\nxvLrYWUcQULSqFCQ/OedJAXpLuTw+9bKk9lLTRHaXRaNFkhIGpUJUsC8k8QXGx7tU3x3w+9a\nrlJzGZD2lxekBRKSRkWCFDLvJClI/TsV5UVVIW4AjUoH6WCZTlogIWlUIkhB806SgtS+bsNt\nPti6GoWDhLZqdFgZR5CQNCoQpLB5J2lB8jTvoIOKIJXBWzU6rIwrSDgalQeSRScGKbgITpnj\nx/4ZJEQ/7wpDb5czSHHLYC6/HlaGQdrzC75dziBFLYO6/HpYGQZpxy/8djmDFLMM7vLrYWUY\nJLsfwu1yBiliGeTl18PKMEhWP4zb5QxSvDLYy6+HlWGQbH4ot8tTgzTN36rOP/sbffn1sDLu\nIGFoVA5ISLfLk4JUCcXODhL+8uthZRxBQtKoGJAOOModoBmkH0Wjs8+1i7D8elgZR5CQNCoF\npCOOcgdoBqkFTs23iFQCSDGWXw8rAx/ahWhUCEh4t8tTg1RDFma3iEQfpMN+DqeZwCK2dGNo\nVAZIiLfLU4NUXeGIBOKIHkgYGhUBEubt8tQgvWrYAp4mkaiDBOOIHkgYGpUAEurt8vTnSKe/\nagfDiCBIF7lqh3u7nEHCLgPliEHC9HOvEPl2eWqQvMw76KAifmXAHNEDCUMj8iBh3y5nkHDL\nwDlikDD9XB3Rb5cnB6k59dDOgyOCIO1rVEkbf5cKEv7t8tQgNT7j72JM6pM7hnA70KhSfi0k\nRd7/cPfTCLfL099HetXi/amH9/06mnfQQUXgZXyORwSPSAcalQ+S1EkgN5z+ql37LZ7t54wL\nRPbHo/jNeJaBXLXb06hSfxcJUs+RcCSJMEhPORnyhOdI/fGIZmjmIjsg7Wg0nyItIP3XGUDO\nvPY1clTaycXmrea/b3GT798B1BFlN8IuM4zrSIZmKWJL975GE0CFHpEWjtwOSVSPSFKdWm4F\nZGKkd9BBRWBlxvMjiqHZitjS7aBRsSApHJUNUvu8te1DwJbx9A46qAiozHSdgWBo1iLWfB9r\nVCpIo0zuHNEFyce8gw4qAikzX6+jF5q9iJ8YJQ/tZpmcOWKQwosAyizXvcmFtlPET4xK+VcY\nSOrtCZIBbv2WvF8BJEUgaqHtFfFUY5rRUNzMBu02H8UADX5LglcgfR7nmyLk09F5NBNYxh0k\nDI0I7qf67XKCAZr8rCDdT/cYxZdXRwduJrgM5PL3GUFaTTuhF6DRzwqSEL8Adcwikdpbv/w6\nOmgz4WXcQcLQiNx+up6+RS5As5+iip7v28nWbFg/Vk4otMMitnRjaERtP91Mg6QWoMVvSfBK\nlPet+YSKRGhv9dYnrEhckDA0orWfGlbRoBWg1W9J8Lp3+z3TOZK/PmFFIl+1Q9CI1H5qWo2G\nVIB2PytIp7rYEKBPWBG+2ABx/DI9NkEpwB0/K0hnuthgevyISGhORWzpPtnFhml2XdyG0x+R\n4BoRBcn4GB+N0NyK2NKNoRGd/XSepRq34eTnSPfHSRaIND8OSyI0xyLWfCNoRGY/tT01QSbA\nfb8lwWdd187yWDmF0FyL2NJ9pnXteplMk72pBHjgd3KQhLAtz5A9NECRC4C0PDYRu+HkQzsf\n8w46qIi9jJ2j7KFBiiAoY9WIxn66sxoNjQAP/ZYEnxCkgSPzAy0MEizAqPvp3qpOJAI89lsS\nrIH0aeR/fytxB53NegcdVMQO0nACi9UOMZBwNKKwn+6ujkYhQAc/C0iVHHb/ydF3BZmE4h10\nUBFrmb3lnMoHCUej/PupaTpDxIaTgvQj6k6bWy0X84Qs2uAddFARW5ndZdGKBwlJo+z76QFH\n+QN08zODVItutPCWi9N8RMCr5/PurfvLCxYPEpJGuffTI46yB+joZwapv57623d0xV7+Ruro\nwopEBAlJo8z76fHi0UWD1L+ctBGvgkHC6ujCikQECUmjvPupwyLsRYN0F3L4fWvlyWyRa38f\nDhjKBwlJo6z7qQNHZYP00w29n+K7G37Xcm3p4kA65qh8kJA0yrmfunBUNkifaryoKsQNoBEV\nkBw4Kh8kJI0y7qdOHJUNUvu6Dbf5YCsWEwEJtaMLKxLzPhKORvn2UzeOCgfJ07yDDiqyKoPb\n0YUViQkSjkbZ9lNHjhgkUNBBRfQyyB1dWBEGyWauHDFIoKCDimhlsDu6sCIMksWcOWKQQEEH\nFVHLoHd0YUUYJLO5c8QggYIOKqKUwe/owoowSEaTMuUChEFyKBOhowsrQhAkAtZzdFo7AUgA\njq4LkmeAeH7jXb5LHJGEav4iJd5bIRyVDxKSRsn302kVSAaJKkggjhgkYIBYfl/T0y2XAEna\nvZbPu9R3gEZZQYJxVD5ISBol3k+XVSCvAlI/uVh+DVHJO+igIn0ZIEenAAlDo7T7qbIK5FVA\nGocLn0KGdlCOTgEShkZJ91N1FcirgFSLYdiw6e2Gx5qdXvSbbm8Fc3QKkKwaEQVJWwXyKiC9\nq+EQXK3WeqqWt81XR6+eT7a3wjk6BUg2jYiCpKt0FZDaT3MT4va9WumpagmC5MHRKUCyaEQU\npJVKlwHJZvRA8uHoHCD5mGeAwX5rlRgk9YcC0n+dYQrubOeecIJukfc/m9+mt7sOSD93Idr6\ntf6a2hHJ63h0liOSRSN6IG1VugpIn1t/x1yIP9ogQSYSh7TjVyQuSFaNyIFk6O2uAtJDNPI+\nxe9mqSdaIPUKXRUkq0bUQDKNGq4CkrzJN/2jCxJsIrF/O75F4oJk1YgYSMbRN4NECKRRIQaJ\nNEjms9irgDQOGxq5SvsWJBozGyaFrgqSVSNSIFmuBl0FpA/CXfPIe+us0FVBwtAo+n5qu6p6\nFZDa9lveNW9C7prH3VsXha4KEoZGsfdT692J64DkYd5BexRRFLouSOEaRd5P7Xf5GCR3kWLu\nrapCDJK/RnH30/GxcrwK0/mhgdRUlB8113q6y4KEoFHM/VRMj5VjVZjSDwukhvSaDX4z8uHt\nhBaJCxKGRhH304mja4NUgV66YxYp2t7qOSMf3E5wkbggYWgUbz+dObo2SKBeziJSrL3Vd0Z+\nYBlyIGFoFG8/3efoMiBNC2uEiBRpb/WekR9YhhxIGBpF208POLoMSO+qhtzmM4oUZ2/1n5Ef\nWIYcSBgaRdpPv3bPjyI2jOaHN7QjerEhYEZ+YBlyINFdIPJrWi8Iq8L0fqcGSYivkBn5gWUY\nJFc/o0opGkb0Q7uP5GPeQbsVsXB0WZAwNIqxn+7dPoraMKbfiUGyccQgBWgUYT89PD2K1TCq\nH/7Qrno4n9F6B+1SRIw3JiBlfNrBLJJqaOevEf5+enS5LlrDuH74IHXmqpJ30A5FhH3CCYNk\n1Eh9Zsz68CX6furI0WVAek7L4f4CHhzzDtqhyM7ErauCtK/R/PRypX3tGaCr3zD4PuboMiDd\npjcd3AB30L2DdiiycwZ7VZD2NcoCkvviaFcBaRIGtCaAd9DHRTAmEgeWIQeSg0bVmqO4IAEW\nGbwKSPObDur2T/Z4mUHavaJ6VZAcNKq0U6TYq+Hyorft3tsohPMk4yi7kbT9nu6qIB1rNI/s\nUlxs+AIt1nmVI1L7mdcDEI0jR9FAOhgxXBWkY42qzYd4IH3BFuu8DEg+5h30fpGjkfdlQTqy\nyvDJM8BDvy/gYp0M0o55B71b5PAMlkEyW7X8jA4SeI3BS4DUr8xOZdLq8ZWgS4J0rFG1/LIv\n4om0/8HXGGSQEoOEOJE4sExhIFXT5brd1XBx9j+PNQYvAZKveQdtLeIyIf+aICFphLH/CZ81\nBhkkd5HCdyMnjhikAI0Q9j+VIwZpA9JPNxz4E9V3gEiBu5H1sQmEZs4BEoJG4fufPimfQVqB\n9NMNvPsbfhCVvIM2FRHCDaPrgoShUfD+t5pMzCBtJq3+df9+XmI1VSsZSONjE6AyPu1ELRIX\nJAyNQve/9SqQDNJ20uoTNPMbHSTbY3w7ZXzaiVokLkgYGgXuf5vVghikzUqr74d4yRF4JpCc\nHvz3buYUIGFoFLb/bVfdYpBWIH3LyZCys3OeaIcLEoCjy4KEoVHQ/mdYvY5B2ryNQlTPrtOD\naIQIEoSjy4KEoVHI/md6qpxBInUfCfCAmF8z5wAJQaOA/c/4VDmDRAkkGEcMUoBG/vtf4NJo\nlwEp4w1ZIEfXBSnjDVnbTT4Gic4NWShHlwUp4w1Z681yBonMDVkwR5cFKd8NWftNPgaJyg1Z\nOEeXBSnfDVmENQavAlKuG7IeHF0WpGw3ZHdW92aQMG7Ihhsv6ASwXDdk91bJZ5AI3JB1n+8d\n1IxvGXJHpEw3ZHuRghfrvAxIPuYddG++HF0XJASN4Psf0hqDDJK7SLDdaHiMz3FWkH8z/mUY\nJGlYawxeDKS/JsXFBjlOGDlikMAgBWkE3P8On1lmkAwgPR+VSHHVrp/4+GW/pOq2KdHLkAQp\nVCNIX+eyhgaDtAbp+ZC79+MJ0MgTpJmjfmEp1+0ENxNchh5ICBpB+rp/DoMGBkkDaVBIjO/f\nSQDSfEWVQXIFCUcjQF+3fqrcMfpAx6JBGvs50B3zIJAWiRgkR5CQNPIYNDBIBj8bSPdPC5x6\n4g/SOANyUIhBcnDvGkoAACAASURBVAYJRSP4oGH/LJZByndE0q6oMkgUj0jzPb7Dq0EMUq5z\npNVtWAaJ5jkS+hqD1wBp1unxFxmktUQMEvyqXYhGTgEC5pwwSHnuI22WgWSQACAhaOQ+aEBe\n0elKILXRZzaIzRNiDBIMJLBGcJtEitrIeSzHXDuxWjjab8+7OkhBGjkMGkBzTviIlAEkA0cM\nEjWQgHO3GKT0IJk4YpCIgQSdA8kgZQDJJBGDRAok8BxIBik1SMbjEYNECyT4HEgGKTFIFo4Y\nJEIgDfOCYHMgGaS0INk4YpDogKRylHs/Te9XBkjb+0fQ8DOUuRhIQuvsGKR9vyXBaUGyXgpi\nkIiApHOUez9N71cGSAgrdaYvcymQVhzl3k/T+xUAkvX8CBJ+hjJXAmnNUe79NL0ffZB2OWKQ\nSIC0XYuGQdr3WxKcCqSDqcQMEgGQDGs6MUj7fkuC04Ak7NfrYOFnKHMZkExrozFI+35LgpOA\ntBl6+29n+jIXAcm8XBCDtO+3JDgFSMccMUi5QbJoxCDt+y0JTgCSA0cMUmaQbBoxSPt+S4Kj\ngyRcOGKQ8oJk1YhB2vdbEhwbJDeOGKSsINk1YpD2/ZYERwbJkSMGKSdIOxoxSPt+S4Ljg+T0\nqCWDlA+kPY0YpH2/JcFxQTq8fwQNP0OZs4O029cxSPt+S4KjguTMEYOUC6SDsTeDtO+3JDgm\nSK7jOsh2pi9zapCOzmEZpH2/JcERQQJwxCDlAenwWhCDtO+3JDgeSBCOGKQsIB1fU2WQ9v2W\nBEcDCfZ6WAYpA0hi1ih4QxikWCABX7PMIKUHyYGj3Ptpej9qIEFfV84gJQfJhaPc+2l6P2Ig\nub3B1yP8DGVOChLu1C0GKQpIYI4YpNQgIU/dYpD8QKo6s4pkfkQMJ/wMZUoFyVEjnA1hkLxA\nquYfBpFAl72h4WcoUyhIOBoxSPt+S77RQfLiiEFKCpK7RgzSvt+S7wggeXDEICUGCXvqFoMU\nCtJ/nWl/63s7wa8ezW2sUVqLMbQDXWaA9QMZypzwiPTPWSM+Iu37LfmOcrEBzBGDlBakf9gz\nThgkfJCcRfIJP0OZU4JUyH6a3o8SSIT3VsKhMUgk/BikaGUYJEQ/Bmml0t5dc8J7K+HQsEHa\n16iQ/TS9X2KQdPMOOqgIgxSgUSH7aXo/BilaGQYJ0Y9BAohEeG8lHBqDRMKPQYpWhkFC9GOQ\nACIR3lsJh8YgkfBjkKKVYZAQ/RgkgEiE91bCoTFIJPwYpGhlGCREPwYJYv8du3AzmZqJ1Bp6\n9MUFyCBdq5lIrTFIDNK1monUGoPEIF2rmUitMUgM0rWaidQag8TP7bOxIRiDxMaGYAwSGxuC\nMUhsbAjGILGxIRg+SPojzhGs0pqJ1Ny6+rKbsbSKWVeFVelUEVJ9lVofQoW2vQ8dpNWiG/hW\nLS1U0ZpbV192M5ZWUevCI3NVMVadOHVZ977iQKpa26bgttKm2MMZJEuNkYJEqOk0ILXWTYnR\n0Ima0RrEbaXC5wgdJMyDG4MEbuhEzWgN4oOEforUtlhBzmczGNUxSF7tRG9G32vKBAl1v4+Q\nikr/gVgZg+TaDh+RXOtDrTMCSKtPYZUxSPSaKR8k1P10rggvyMr4Mag2BgnYSuxm1tWXCJKe\nKJzqIoDEQ7udFiLuepXeFoO0W5kyDsaoTx1U49SIF2AykM4xs0G7Hx6vmfJnNmwShVKj+ju8\nQr1elNoSzGxgY7uiMUhsbAjGILGxIRiDxMaGYAwSGxuCMUhsbAjGILGxIRiDxMaGYAwSGxuC\nMUhsbAiWE6SfWoj6V/9OmAPSv/4xzfVYXD4/90rUPw4BLPXslzFF5dxK0cYaOVo+kN6V6K3W\nvnUSyeg0f/kaK64+RyEs9RyUMTTo3krBxho5Wz6QKvF4t+2zElqHgSDSTTy6xL1r0RyFsNRz\nUMbQoHsrBRtr5GzZQPoV9/73U8hj99+96zfk5vbpeN+V//Q/+0+Tk+xkuv9+HqJPVJ+q+5LI\n8cNnLHfv/va+iftHb+ZV1WM96zJz86uolvbsJd7S5z01sC5TmrFG7pYNpLv4Gz68WqlUb82w\n7Z/+iHxfizQ7jcntvW7t6H+fE34Xz7mV/nvxe+t+PPRmavFQRFLLzM3rUS3t2UsMH+RIom9g\nXaY0Y43cLRtI2pH4Jrrz2dekRtNt3t/cx80i6U7tt8xeIwcdTTeG/9Rzhd3A/tb8vsdWHl2/\n2nn+rmtotBjUMnPzeoNLe/YSMpK2HmRt2k2Z0ow1AuQqsLx/w3rL7+d3vajxUX0WuVSnzqv/\n411+6nL1Vq4Ifcve7fbXl33LH59tDW89BqXM3Pw6qqk9e4kxktvUwLpMYcYaAXIVWN6/Ya3l\nejhCq4q0G5E0p1YI7b96ha/mUcu+avRsjTWsYtDKmKKa/mcvsd6rDGVKMtbI3bJpPI+/2+73\nQ9x+nu8DkXSnfZH6L6q1SKsatkWWMtLWURkTrpU4GUiskbtl03i6IvRXPcZNmi7hrIYNunaz\n03hIbtvNsEGMxRU5W2MNikhqGb15JSotfmMJddjQtpsypRlr5G75lJ7vUbzkJv2Np6LyX9Od\n/vXnj1V3TF6+1p16r07qWp4u1p92OZHtzia7DvTTzNeUFJG0GhSR1DJz89uohvbsJdQT2bbd\nlCnOWCNnywfS+zZfuew2Rh0CvKcLkv3X30uWFqdquowpJV5dWm3Hiqv3WiS9mbGeTRm9+dl9\nac9eQr20Kn3WZUoz1sjZco49no9qmsf1EF3nMWfvVY+3zJpKfC8H/Nnpp0/uu/+/LP2+qzf7\nhgliVSMP5bpIq2bGetZlluZ196U9ewnlZl/vsy5TnLFGjlb2IJ6NjYgxSGxsCMYgsbEhGIPE\nxoZgDBIbG4IxSGxsCMYgsbEhGIPExoZgDBIbG4IxSGxsCMYgsbEhGIPExoZgDBIbG4IxSGxs\nCMYgsbEhGIPExoZgDBIbG4IxSGxsCMYgsbEhGIPExoZgQSBNa+vNK1WYrZ7+Znpr1VANqDXX\nmLfNPR9zLXdluYvA91HZA/qrtYZetbISx+rT42mpI9BAGo2L7jRvUzWg1gDeq+ZSa9TKpYja\nYAESgPQzrdZifGvV8EIBl40AgmRo7q9Sa5lVCn0flTWgj9BWenrNa0OZPlVxFhsCaTQ7r/WI\no5GpudQatf2SXvJXmACBIDk4vaasmt9aJRdUusu3a+C12Nqaq37mGj7N8iaP0PdRWYP6kTU+\nxPf434cMZfxu+8n4rshwA2k0/Hw/NvtrHI1MzaXWqF/Qa1j4K0gALJCGVcVqcXuuQ75Vwzf6\nW6u0op9JxtWLn+aXSK1a/AwJlr+EfEtV3Y8N9BdG6c311mgNLxFMuk3rtd2XZdqGn/q7r/Ze\nn6V3x3d5MHrNLzoYa6zNn9rKcx/ZN5BG0/cD/NE10pvrLblGbVX9TT19iACYII1rVupVfovf\n4RvtrVVK0WfX1T36g/v6xU/LS6TWLTa9/2+/MuFjOuCvXzKlNifto75vR+3tNu+janSRVu++\nsr8+ay1Spe2dSo2mT91GxXi1H0ij6ftXT3Z0jfTmpKXXqG2UiAMEQDlHaodYvvWXSfX2mnsD\n28FVCvTsxw3rFz8tL5FSWhxrlZmXrIhhTelm88KobXPf456wGX/r76OqN0tPr999ZX191jY/\nWiT1UESYP3VJ+LbUE2J+Gq02KbZG84f0GinNBwmACdL6ZVLSqq4j2gdpLt5uX/w0vURKd27H\nQdPwto7X+G6B1Qujts3dxxdXDfZQFnvW3kf1Vl85Nw3jtXdf2V+fZdq0JZI/2W5fj+lTV0+M\nF5L5abTapNgaKYPh1BopzQcJgDm0MwDTjwd2QRqPyarfrPz8EqlNi7IT7btIdZCkeW+bU8J4\nDmMJxV9/H9Vqazbvvlr/f27BFIMyQpLvGx4CMHyydjVB5qeRHkp0jUwZS6SR+scQASKDpISt\nvbVqscfU+2xTsLxEatti1880+g6wStC2OTW8ZjwJ0BK6vI9K35rtu69gIGlXg97z+cDmUyaQ\nxHp75HFS65tjabRtLr1G6h+pgGQ4iCpha2+tWqzqz/A+S4KWFz/1//+YRXqKpj+VF0Ob9eaF\nUdvmtIzelPGw/j6qaRt6x+WT8u4r6+uz1nbTr9oN8cwn5ptPKUDa12j6/q7fooil0ba59Bqp\ncVABaXyF06bK8Rv1rVWz/Y23Jx79San8tLz4aXmJlKHFmxg7rP5U8nv7wqhNc3f15b4v5Y04\n6vuolhdiqW/Q0t99ZX191trkm7SHO0VjyM/20585mD5FO0fSPu5rNPyUN3bUv0XTaNtceo2U\niMmcIymXVk2ZVd9aNTs04wis673G75YXPzX2g3A/bv/tv5CD40or126bkzZdEZr+N19aVd9H\n1V82rYddToxv0Fq/+8r++qyVDdMWxPQSoJ85HtOnaFfttI/7Gs0Hpz+1aDSN9OakpddIiTjn\nVTv9o7zZ92sFSXtr1fRdNXVG3YcJuPnFT/NLpMwtDp1X12j/2qjtC6OU5qRp9yj0gYPyPqq3\nvHfYO81v0Fq/+2rn9Vkrk3Pt6mWn/Ll1I/mhQcOnJPeR9jUadsZb89GKRtNIb05aBo2W9vLd\nRzLXGGeqy9r+ht4KMqxtnH2jnK0cNhplZoO5KdbIVGO2mQ3ruuQhunGekxVm9XCrG5TOynXy\ncA6QIs21WxlrZLWMc+1WNo5Tt3Pw8U1Ms7pB6fxzTVUOkCLN/l4Za2S1jLO/19aN98excGyr\npgsssHQ+HXviDCDFeh5pbayRxXI+j8TGxjYYg8TGhmAMEhsbgjFIbGwIxiCxsSEYg8TGhmAM\nEhsbgjFIbGwIxiCxsSFYCEj/Vrb5YscAvpGqJRDCviuKvom2hEAys4SwJJpBIlttMEQJNSKQ\nTAYJ0bWwEBik0kNYEs0gka02GKKEGhFIJoOE6FpYCAxS6SEsiWaQyFYbDFFCjQgkk0FCdC0s\nBAap9BCWRDNIZKsNhiihRgSSySAhuhYWAoNUeghLohkkstUGQ5RQIwLJZJAQXQsLgUEqPYQl\n0QwS2WqDIUqoEYFkMkiIroWFwCAVFoIQYvXH2RgkstUGQ5RQIwLJTBBCvySg/sfZGCQq1e70\ndvQ1IpDM+CGM6/Rrf5yNQSJS7V5vR18jAslMANLXF4NEMgR11LAnEn2NCCQzfghfX71I2h9n\nY5BoVLsrEn2NCCQzeggGiRgkGiFsOjsGiW4IJo4YJBohrDs7PkeiG4KUyPDH2RgkAtX2IvFV\nO8ohmDlikGiEsNvZMUiEQrBIxCDRCGG3s2OQ6IRgk4hBohHCbmfHIJEJwSoRg0QjhN3OjkGi\nEoJdIgaJRgi7nR2DRCOErx2JGCQaIex2dgwSiRB2OWKQaITQOopEXyMCyYwTwj5HDBKNEFxF\nYstlUiE3TwYpX7X7HPERKX8IPUe7Fc3Gr3XJZV/uvV2ABextuffi/CH0PV0CkAK25JyJh/ge\nDL7/8REpewiDQgwS6RC+DkcNDFLmEMaejkGiHIJDZ8cg5Q1hGjEwSIRDcOnsGKSsIcwjbwaJ\nbghOnR2DlDOE5QyWQSIbgltnxyBlDEG5EsQgUQ3BsbNjkPKFoF5RZZCIhuDa2TFI2ULQ7kww\nSDRDcO7sGKRcIeh3+BgkkiG4d3YMUqYQVnfKGSSKIQA6OwYpTwjrGScMEsEQIJ0dg5QlhM3M\nLQaJXgigzo5BSh6CEGI7A5JBIhcCrLNjkFKHMGC0nknMIFELAdjZMUiJQzBzxCBRCwHa2TFI\niUMYOBJOvvMfZ2OQElUL7uwYpMQhGFfJZ5CIhWB8jI9BIhSChSMGiVIIlsdhGSQ6Idg4YpAI\nhWB7rJxBIhPCzoK3uxXNxiDFr9a6PAODRCME0+0jt2qXRG9A+rkL0dYvQiLRSzzM17Oz2wGJ\nnkZF67nLkSdIn5t8cVwrxB8dkcglHubr29lZQaKoUcl6DhwZzo4cqrWC9BBNp1D7K2o6IlFL\nPMzXu7OzgkRRo5L1HDhCBqkTaP5HRSRqiQf57q/uvV+RxShqVLKe1ut1DtUySIlCOFjde78i\nBilFCAcchQ3tGvGgIxKxxEN8j1b33q/IYhQ1KlfPI458LzZUorfqTUckWomH+B6tSux5sYGg\nRsXqebhwtPfl7++bELfm46ARg3Tge7y69/5frUZPo1L1PObIGySABWxJqYmH+AZ2dnxDNn4I\nDhx5glS7jLsZJBff0M7OChJFjcrU04UjT5AqyBEqYEvKTDzEN7izs4JEUaMi9XTiyBOkV924\nnMIySEe+4Z2dFSSKGhWop3UGJKRaK0hiNjoiEUk8xNdRJD+QKGpUnp6uHDFIGUNA6eyAIFW9\nTR8YpANXZ45Qr9rlFYlE4iG+OJ2dz1W7avy3WNiWRPGlEIIzRhEuf2cTiULiIdW6i4QNUjX/\nYJB2DcCRL0ifxnKzL59IBBIPqRars7ODtKuRLhGDZDIIR54gva3TT/KJlD/xWTo7K0hWjYa+\nbhl9/9eZrZIrm5QIvdLNpNW6k+ddbyZEskhuFkUk3Q40mn/0FtAlEOiVoviCjkf+V+303xRE\nyp14iCtMJN+rdnsarT55b0k838wh9F0dVrWeIOURqSDtgSLhglSZPvpuSUTfvCH0XV18kGzD\nhpwilaM9VCQ/kPY14qHdrg1Dhvgg2U5kc4pUjPZgkXAvNiwaKd2e35ZE9c0Zwjj0jg+S7dJq\nTpFK0R4uEu7l71EZ7Z45g7SY8vajBCBZLKdIhWjvIZInSBDz2ZLIvplCUN8ilg+knCKVob2P\nSAxSOl/tbXwJQLr3X4gbofUAitDeSyRPkAhqRF9PbRnI+CA1wzVVQWmFmhK0V+4fRQeJokb0\n9dSWgYwPUjWsg/vixyhArup92OggUdSIvp7acqrxQbLekGWQdkybzxAdJIoakddTX5Y4Pkh3\n8fjI66uU1pUmr70+Lyg6SBQ1Iq7n+gmx+CDNN/tc3hkCatLfl7j2ASL5gURRI9p6bp60THDV\nbrzZ57S6BqhJf1/a2oeI5HnVjqBGpPXcPrGcACSIgZr09yWtfZBIniAR1IiynoYZ+QxSJt/d\nzi5CtcNfMQzWZBLfxCGYnmyJDVJ3Atv9/K3EndKwgbD2gSL5gERTI7p6Gp8Qiw1SJa+o/vUn\nsi4rtEcLL0m1CCGEiuQDEk2NyOppftIyMkg/ou60udXy5nlDRySy2geL5AESUY2o6ml5Yjky\nSLWQD4zJmScfUW0kYZBWFi6SB0hENSKqp+3J/8gg9bfKf/uOjtJdc6LaI4jkARJRjWjqaV1B\nIzJI/WsOGvEiJhJN7TFE8gCJqEYk9bSvRBMZpLuQw+9bK09mCU0/Iak9ikgeIBHViKKeOys6\nRb/Y8Gif4rsbftfih45IFLXHEcnrYgNJjQjqubcyWmSQ+rf8youqQtwcNLowSEgieYBEVCN6\neu6uMBgZpPZ1G27zOV1YvTBIWCJ5gERUI3J67q/UGRskoEULL0m1/iGgieQDEk2NqOl5sOIt\ng5TJd60RkkgMUiTfo5WjGaRMviuNsERikKL4HkrEIOXyXVxRRWKQ8H2H1YIO3mTAIGXy1Tq7\nCNWa/4phAdEVqeeAEaJGDFKEEFw0YpAyVdv7DhwJF1eUEJZEM0jurk4cMUiZqu19v/TVgqKH\nsCR6A9LPXYi2dllX42ogOWGUBCR6GhHRc8CIAkifm4yjFcMahDREIqK9I0fxQaKoEQ09XTlK\nANJDNHJW8S+lCZE0tHflKD5IFDUioWcvkQtHCUCSU/Onf1REIqG9M0fxQaKoEQU93SVikHL5\ntnFEYpAQfQESpRvaNZTedEBA+0giBQ3tSGmUX0+IRCkuNljeT3ptkGKJ5HmxAaLRVUxKlDWA\nzfDg2/h+0ktbdpHWBtAoAPOSOsZeogwhLInmG7KHBjoeJbmPRE+jzHr2EjFIeK6xRMqyZcEQ\nJdQor55DV0cHJKEaHZHyag/t7CKDRFSjrHqOQwYGCc81lkgMEqIvdrXT0JsOSFCLFl6SaiEi\n0QGJqEbZ9BweP8oWwpJoBmnHvjw6OwYpZbUKR8RAaugNG7Jp/+XT2SUAiaBGmfQcOBL5QrCC\n1BAcf+fS/surs4sPEkWNcoGkPn5ECqRKvGrx/tSUpuhn0l69fUQKJIoa5dFzfIyPIkhdL/ct\nnu2H0hT9PNprt2FJgURRoyx66hyRA+kp15SmNGzIor0+nYEaSOQ0yqHniiNaIN3F71vc5KsV\n6YiUQ/vVtCBSIFHUKIOem5lbpECS6tTyPJbQFP0M2q9FIgUSRY3S67mdAUkKpPZ5kw+8uK3Q\nHi28JNVCRCIFEkWNkutpmElMCySIRQsvSbUQkWiBRFCj1HqaZuQzSHiusURikBB9Eao1PtlC\nC6R7/4W4XfcJ2cDOLgFIBDVKq6f5CTFSIDXDpSBSJ7JptQ/t7OKDRFGjpHpanrQkBVI13C5/\nUbq0mlT74M4uPkgUNUqpp+2JZVIgTeJQEiml9uGdXXyQKGqUUE/rk/+kQLqLx6dtPw2l6Sfp\ntLeukk8KJIoapdPTvoIGKZDe01JPLiu0RwsvSbVbXytHtECiqFEyPXdWoiEFUtfR9Us9OS2Z\nFi28JNVufO0c0QKJokap9Nxb0YkWSBCLFl6Sate+SBrxfaR41e50dalCWP9xNgZpNCyNGKRo\n1e5zRAmk7gS2+/lbiTulYUMS7dE0ig0STY1S6HnAESWQKnlF9a8/kXVZDzdaeEmq3YiUPgQf\nkGhqlEDPwwVv6YD0I+pOm1stb567TC2OFl6SaiEi0QGJqEbx9TxeOJoOSLXoRgtvOfPkIyo6\nIsXXHrOziwwSUY2i6+mwADsdkPpb5b99R0fprnl07VE7u8ggEdUotp4OHBECqZL/acSLmEix\ntcft7CKDRFSjyHq6cEQIpLuQw+9bK09mCU0/iaw9cmcXGSSiGsXV04kjQiD9dEPvp/juht+1\nXKWGikhxtcfu7KJfbCCpUVQ93TgiBFL/TkV5UVWIm4NGpwAJvbOLDBJRjWLq6cgRIZDa1224\nzee2rsYZQMLv7GLfR6KpUUQ9XTmiBBLQooWXpFqISIRAslolbfx9IpCcOSoCpKwixdM+RmeX\nDyTl1yJSoi2Jpqc7R2WAlFOkaNpH6ewYJMxqARwxSIiusUQqAKRK/X0WkCAclQBSXpHiVHs0\nlzhBCLggTaPvRaP/OoNVQsykRLljcDU3kM4n0ldJIrnYpM2Jjkig41E5R6RsIsWoFnQ8KuKI\nNAl1HpD6ri5vCIe+VpCm+VvqxbncIkWotsco++jSEySKGuEn8wuqESGQKqEYHZHwqx0OR0WC\nZNXoZEO7L7BGhED6UTTS5nGdbGg3DuuKBGlXo9V5bLEgDS8rB2pECKTWNjU/r0jY1U6nR0WC\nZNVovlle/syGgaOsIbj6WkGqzQuzZxUJuVofkUiBZNHIbIm2BDOZA0ciZwjOvlaQKsjcu2jh\nxax2uVxXKkgUNUIFqeeocJBetdsCnoWCpAwaygWJokaYyRw4Khyknat2JwBJ46hYkChqhJhM\nlSPyu9Q1QVIH35lCgLheEySNI/K7lBUkkEULL1K16uA7UwgQ1wBlkmuEttWrGSfUd6kl0VcC\nSe/sGCTELcHa6vXMLeq71JLoNUgNvWEDVrUrjsoFiaBGSFu9mQFJfJfafRkzOZFwqv1ac1Qs\nSBQ1wtnq7Uxi0rvUvx2QKvGqxftTD+/7pSESSrWG6d6lgkRRI5StNszIp7xL9X+0gdT1ct/i\n2X4oLT6IUa3psYlSQaKoEcZWm55sIbxLDX/cAekpJ0NSGjYgVGt8+qhgkMhphLDVoRrRAuku\nft/iJt+/Q0ek8GrNT/GVChJFjcK3OlgjWiBJdWp5HusyMTJaeMjVWp6GLRUkihoFb3W4RrRA\nap+3tn0It2U8o4WHW63tqfJSQaKoUehWI2hEDCSIRQsPtVrr6gzFgkRQo7Cttq+gQXOXUv44\n2+lBsq9ywiDhbUnQVu+sRENyl1L/ONu5QdKne2cJwd81GKKEGoVs9d6KTvR2qdUfZ1uB9HnQ\nu2vuX+0uR8WCRFGjgK3eXRiN3C61/qMNpDvB6Sfe1a4em8gRQoirLe0UNfLf6v0FBqntUps/\n2kAS4tdBnVJA0h+byBFCiKst7RQ18t7qg4U6qe1Smz/aQLoRXA/Au9rNNNX0IYS42tJOUSPf\nrT5a8JbaLrX5ow2k9635UBPJr9pxXGfnqFiQKGrkudWHC0eT2qVMf7SB1P7SG397VXvMUbEg\nUdTIb6uPF2CntEsZ/2gDieKJrE+1B9cZUoQQ7GpLO0WNvLba4UUGhHYp8x9tIFE8kfUCaf86\nQ4oQgl1taaeokc9WO3BEaZcy/9EG0p3giaxPtYfjuvghBLva0k5RI4+tduGI0i5l/qMNpPb+\nILf4ILza7WPlyUNAcLXmnaBG8K124ojOLmX7ow0kimumgat146hYkChqBN5qN47I7FLWP54Z\npL25W4lCQHE9M0iOElHZpex/tIEEsmjhBVXrqlGxIFHUCLjVMTRikBBdO19njRgkxC2BbXUU\njSiB9Gnkf38rcXc6m40WXkC17hoVChJNjUBbHUcjSiBVctj9J0fflcsklGjh+VcL0KhQkGhq\nBPGNpBEhkH5E3Wlzq+Vini4LAkQLz9sVolGZIEE1Imbj22FzhxHBVJBq0Y0W3nJxmo/Yvnq+\nADupRqpBNUrUJTj6Lm9ZzhYCbrVLolWQ+uupv31HR+nSqrMrTKMyj0hENXLzBXOUfZc68jWD\n1L+ctBEvYiK5uvbHo7whoFZrzDlRjRxBcpgBGTkE5GrNIN2FHH7fWnkyS2hdaUfXvq+jnniI\nqzHnRDVy83WbcRI1BORqzSD9dEPvp/juht+1XFuaikhursOYgXriIa7GnBPVyMl35IjAhcPI\nIH2q8aKqEDcHjSiBNI29qSce4mrMOVGNXHyn49H5QWpft+E2n+OF1WjhgV3nc1jqiYe4mpNO\nUyMH3/ky2Bf5kAAAEFlJREFUwwVAAlq08ICuX1+xRSIEEk2NDn2Vy3UMUiaRABqRTzzENRii\nhBod+SbQiEEKdNVuTVBPPMQ1GKKEGh34anePGKRMIgE0Ip94iGswRAk12vdNohGDFOS6ulNO\nPfEQ12CIEmq047tZgJ1ByiTSjutmkXzqiYe4BkOUUCO77/ZFBgxSJpF2OrvNzC3qiYe4BkOU\nUCOrr2GBwfODJFSjI5K9s9tO3aKeeIirMedENbKDlEwjBsnb1TB1i3riIa6nACmdRoRAknav\n5fMu9d1Bo9wgjRpdCySiGpl9zQuwXwKkfnKx/NpFpWjhObjOGiURiRRIFDUy+lpeZHAJkMbh\nwofSsMGu0b9UIpECiaJGJt/hKb7tQxOXAKkWw7CBUm9n0sj2sgnqiYe42tJOUSOT75flMb5L\ngPSuhtPYymWtp2jhHbpaH7WknniIqy3tFDUy+FoXjr4ESO2nuQlx+3Z6JVy08A5drY9aUk88\nxNWad4IabX2Np7BpQ0hQrR0kiEUL78jV/sgy9cRDXAOUSa7Rxtd8Cps0hBTVLokuEaSdR/+p\nJx7iGgxRQo103+20oOQhJKp2SfQGpJ+7EG39IiSSsbMLr7ZgkOhppPnuc3QNkD63/o65EH90\nRFq55hCJFEgUNVJ9j97fewmQHqKR9yl+KS31pLkerDBIPfEQV1vaKWqkgXSwet0lQJI3+aZ/\nVERSXY9W6qSeeIirLe0UNVp8LdMZUoaQsNpSQTpc8ZZ64iGuRYJ0zNE1QBqHDY1cpZ2KSIvr\n8cLR1BMPcbWlnaJGk68DR9cA6UPwrvns6rAAO/XEQ1xtaaeo0QySw+relwCpbb/lXfOG0l3z\nydWBI/KJh7jaE09Po8nXZXXvi4AEsGjhGV1dOCKfeIhrgDLJNRp9XThikHKJtGiEWi2DhLgl\ng68TRxcBqanIPcY8a4RbbbkgEdSo982rES2QGoLrAQA0Ip94iKuVI4IaSd/MGtECqTK/dKfq\nbPq9vLk0WngbV1eNyCce4mpjg6JG4mB+XYIQaIFk7uWq6Yf++t9o4a1dnTUin3iIqw0kghqB\nXg57CZCmhTUIifSvddeIfOIhrjaQ6GkEe8nyJUB6V7XtNl+11igVSACNyCce4moDiZ5GtmVO\nAquN55tiaGc9ka204fd/ndmUxjWpUZqWCjF6Go2XvZO0RdWcQZpHDWlPZCFjhkgh5KrWKhlF\njY5vH0UNAeib4D6S1arNhxQiATkin3iIq6syFDQCcHRxkCrDp2jhzdZzdKbEQ1z9OUquEYVk\n0gJpGTZUD/WMtlp+JhQJrBH5xENcbbxQ1IhAMqmC1NmiUrX8Uq4KRQtvNLhG5BMPcT0GiY5G\nBJJJC6TntBzur/rgWDVdCtJumscVablVfqbEQ1xtINHRSLk1kT+ZtEC6TW86uDk8yhwtvH/a\nkk5nSjzE1ZZ2Mhqpt/jyJ5MWSJMwTmsCRAtPX9LpTImHuNrSTkUj7VZ5/mTSAml+00Hd/ske\nLxtIyiPLZ0o8xNWWdioaaVNO8ieTFkjKmw6EeZJxGpHUexNnSjzE1ZZ2Shr51HsmPa0gtZ95\nPQDRHGgUUSTtHt+ZEg9xteadgEabxybyJ5MYSBCLFp5+r/xMiYe4BigTW6Pt40f5k8kgrSxE\nI/KJh7gGQxRNI8NjfPmTSQikfmX23I8xB2lEPvEQV2POCWnk+zr5M+lJF6QwjcgnHuJKGKTt\nMpD5k0kIJKjFCC9MI/KJh7gGKBNXoy/Dcqr5k8kgaWaYk3+mxENcgyGKpJGJIwLJJAbST9W2\nf6L6zijS+tmWMyUe4mrNe16NzI+I5U8mLZB+uoF3f8PPRSX88IzPiJ0p8RBXK0dZNbI8aZk/\nmbRAuom/7t/PS6zW0EgikuVx2DMlHuJqS3t2jeCbghoCgm98kLrO7uk0qziCSLbHys+UeIir\nLe25NfLYFMwQMHzjg1SJ90O85Ag8uUjW5RnOlHiIqy3tmTXy2RTEEFB844P0LSdDys7ucBIX\ntkgoGpFPPMTVlnaKGhFIJi2Q2kZUz67Tc9EIVSQcjcgnHuJqzTtBjQgkkxhIEEMMD0kj8omH\nuAYoE0Gjg1Xy8yeTQfp3sCrxmRIPcQ2GCFOjo7dN5E8mMZCy3OxD04h84iGu1rxn0MgwBdK7\n3jPpaQUpy82+g9VUz5R4iKuVowwaHb6tPH8yaYGU42bf0arEZ0o8xNWW9lwa7S5LnD+ZtEBK\nfrPP4U1vZ0o8xNWW9gw3ZI+X986fTFogpb7Z5/LGxDMlHuJqS3v6G7IOy+TnTyYtkBLf7Ds8\nh/WrFt2XFEjJb8ge93UUkkkLpMQ3+w7PYf2qRfclBVLqG7IuHBFIJjGQIBYensOY4VyJh7gG\nKIOtUeimBIaA7nsykL6cODpV4iGuwRChaRS8KUEhRPBNBtJfk+BE1pGjUyUe4nqQ/XQahW9K\nSAgxfNOA9HxUIsEVoRgakU88xHUv9dQ0IpBMYiA9H3Kdp8fTQaMwkZwlOlXiIa47FJHTiEAy\nKYE0KCTG9+9EFSmORuQTD3Hdo4iaRgSSSQiksZ9zumPuIZJ6OhRJI/KJh7gacx5ZI9UgGhFI\nJimQ7p/WceoJ3NS1QaVGURo5vcXVSDXWCGTJjkjjKrv9Z1Bfd6oeDOJq1ivJEcmwSn7Qpvi6\nktfTDFLU8bcCUt/XIW1JGl86ICU5Rxo5KiuZlECadXr8RQSp7+sumniIqzXv0TQabToelZVM\nYiC10e5RTK9QGMYMF008xHUv9VHvI01TictKJj2Q2kh3zVWOrpp4iOtB9mPNbBg1EgySg+8R\nSG4GDk85h71o4iGuAcqEayQH4WUl80IgCeXxo4smHuIaDFGIRgK2JQSSeR2QVI2umniIazBE\nARpBt4RAMi8DksbRVRMPcQ2GKKFGBJJ5FZB0ja6aeIhrMEQJNSKQzIuAtNLoqomHuAZDlFAj\nAsm8Bkhrja6aeIhrMEQJNSKQzCuANEqkLnNy0cRDXIMhSqgRgWReACSDRldNPMQ1GKKEGhFI\n5vlBMml01cRDXIMhSqgRgWReASTDKicXTTzENRgid41EqEYEknl2kMwaXTXxENdgiOAaeW8J\ngWSeHCQLR1dNPMQ1GCJXjf6ZOSosmecGycbRVRMPcQ2GyFEjG0eFJfPUIFk5umriIa7BECXU\niEAyzwzSfClou5rqRRMPcQ2GKKFGBJJ5YpB2NLpq4iGuwRAl1IhAMs8L0jxkMK3ufdHEQ1yD\nIXLQ6B+SRgSSeVqQrENv6KacKfEQ12CIABoZRSormWcFaZ+jqyYe4hoMkbtGwZ0dgWSeFKQD\njq6aeIhrMETOGoV3dgSSeU6Qjji6auIhrsEQOWqEsSUEknlGkPYuBcE35UyJh7gGQ7Sr0XFf\nV1gyTwiSA0dXTTzENRiifY0OOSosmecDyYWjqyYe4hoM0a5GxxwVlszTgXR0CgvelDMlHuIa\nDNGeRi4ilZXMs4HkxtFVEw9xDYboQCPMzo5AMk8GkiNHV008xDUYon2NUDs7Ask8F0iuHF01\n8RDXYIh2NTo4PwJuCYFkngokZ46umniIazBECTUikMwzgeSu0VUTD3ENhiihRgSSWRRIVWd2\nkQAaXTXxEFdPdLJoRCCZJYFUzT9MIkE0umriIa6eHGXRiEAyTwMSSKOrJh7iGgGkaBoRSOZZ\nQIJpdNXEQ1zxQYqnEYFkFgrSf51pf+tvlkuN/PRnQzLWKI8hHpEGw+8TztSDQVz9BM2jEYFk\nFnpE2oj0D6TRVRMPccUHKZ5GBJJ5GpD+ASS6bOIhrhFA+gfhqLBkngck8ltNLoT0IJ04mQwS\nomthITBIpYfgCdLBXXPqW00uhBgg5dGIQDKLAkm3aOElqZZACFFAyqIRgWQySIiuhYXAIJUe\nwpJoBolstcEQJdSIQDIZJETXwkJgkEoPYUk0g0S22mCIEmpEIJkMEqJrYSEwSKWHsCSaQSJb\nbTBECTUikEwGCdG1sBAYpNJDWBKNOaP+v2MXQtWWFi6OcTIjVcsgFVItknEyI1XLIBVSLZJx\nMiNVyyAVUi2ScTIjVcsgFVItknEyI1XLj++zsSEYg8TGhmAMEhsbgjFIbGwIxiCxsSEYHkj6\nI86oVVbIdU/14VZbqdVipwLHYgQWSaTSNEIDabXoBmqV6ICu6keumiZEbbwNjrHVpWnEIOHX\nzCBh1VuQRpRBmuqNxFEskchyFA9zfJGK06gIkGKdIrUtcshRzujwLCpIkU6R2rYMjaiDFCOZ\nU32xREKvFssidnbYdRenUQkg6R/wKo4k0uoTIYsJkv4Brd5iNCIOUrw9M4ZIlfEjHYvd2ZUA\nUiyNaINURak62rAhTrh4FrOzK2VoF0sj0iApl0CxcxnlRHYRiSRHMTs77K0uTiPKMxu029DY\nFau/0erVq6dnkWafxBCpNI14rh0bG4IxSGxsCMYgsbEhGIPExoZgDBIbG4IxSGxsCMYgsbEh\nGIPExoZgDBIbG4IxSGxsCJYLpJ9aiPpX/06Yg9G//jHN7FhcPj/3StQ/DgEs9eyXMUXl3ErR\nxhoBLA9I70r0VuuxuIhkdJq/fI0VV5+jEJZ6DsoYGnRvpWBjjUCWB6RKPN5t+6yE1mEgiHQT\njy5x71o0RyEs9RyUMTTo3krBxhqBLAtIv+Le/34Keez+u3f9htzcPh3vu/Kf/mf/aXKSnUz3\n389D9InqU3VfEjl++Izl7t3f3jdx/+jNvKp6rGddZm5+FdXSnr3EW/q8pwbWZUoz1ghmWUC6\ni7/hw6uVSvXWDNv+6Y/I97VIs9OY3N7r1o7+9znhd/GcW+m/F7+37sdDb6YWD0UktczcvB7V\n0p69xPBBjiT6BtZlSjPWCGZZQNKOxDfRnc++JjWabvP+5j5uFkl3ar9l9ho56Gi6Mfynnivs\nBva35vc9tvLo+tXO83ddQ6PFoJaZm9cbXNqzl5CRtPUga9NuypRmrBEwXwh1wBvVW30/v+tF\njY/qs8ilOnVe/R/v8lOXq7dyRehb9m63v77sW/74bGt46zEoZebm11FN7dlLjJHcpgbWZQoz\n1giYL4Q64I1qrdbDEVpVpN2IpDm1Qmj/1St8NY9a9lWjZ2usYRWDVsYU1fQ/e4n1XmUoU5Kx\nRjDLovM8/m673w9x+3m+D0TSnfZF6r+o1iKtatgWWcpIW0dlTLhW4mQgsUYwy6LzdEXor3qM\nmzRdwlkNG3TtZqfxkNy2m2GDGIsrcrbGGhSR1DJ680pUWvzGEuqwoW03ZUoz1ghmedSe71G8\n5Cb9jaei8l/Tnf71549Vd0xevtadeq9O6lqeLtafdjmR7c4muw7008zXlBSRtBoUkdQyc/Pb\nqIb27CXUE9m23ZQpzlgjkOUB6X2br1x2G6MOAd7TBcn+6+8lS4tTNV3GlBKvLq22Y8XVey2S\n3sxYz6aM3vzsvrRnL6FeWpU+6zKlGWsEslzjj+ejmuZxPUTXeczZe9XjLbOmEt/LAX92+umT\n++7/L0u/7+rNvmGCWNXIQ7ku0qqZsZ51maV53X1pz15CudnX+6zLFGesEcDKHsizsRExBomN\nDcEYJDY2BGOQ2NgQjEFiY0MwBomNDcEYJDY2BGOQ2NgQjEFiY0Ow/wPTHCiQtlsOkAAAAABJ\nRU5ErkJggg==",
      "text/plain": [
       "plot without title"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Given a dataset this function builds a scatterplot with the linear regression trend line.\n",
    "ss.scatter.plot <- function(dataset, title='Calculated vs. Assigned Smart Scale Scores', title.size=22) {\n",
    "  ggplot(data=dataset,\n",
    "         aes(x=Calc.SMART.SCALE.Score, y=SMART.SCALE.Score)) +\n",
    "    geom_point() + \n",
    "    geom_smooth(method=lm) +\n",
    "    labs(x='Calculated Smart Scale Score', y='Assigned Smart Scale Score', \n",
    "         title=title) +\n",
    "    theme(plot.title = element_text(size=title.size))\n",
    "}\n",
    "\n",
    "# Build the linear models for each area type.\n",
    "a.lm <- lm(SMART.SCALE.Score ~ Calc.SMART.SCALE.Score, data=subset(data, Area.Type == 'A'))\n",
    "b.lm <- lm(SMART.SCALE.Score ~ Calc.SMART.SCALE.Score, data=subset(data, Area.Type == 'B'))\n",
    "c.lm <- lm(SMART.SCALE.Score ~ Calc.SMART.SCALE.Score, data=subset(data, Area.Type == 'C'))\n",
    "d.lm <- lm(SMART.SCALE.Score ~ Calc.SMART.SCALE.Score, data=subset(data, Area.Type == 'D'))\n",
    "\n",
    "# Build plots for each of the four area types.\n",
    "sp.a <- ss.scatter.plot(subset(data, Area.Type == 'A'), \n",
    "                        paste('Fig. 4A: Area Type A (R-Squared = ', round(summary(a.lm)$r.squared, 3), ')', sep=''), 11)\n",
    "sp.b <- ss.scatter.plot(subset(data, Area.Type == 'B'), \n",
    "                        paste('Fig. 4B: Area Type B (R-Squared = ', round(summary(b.lm)$r.squared, 3), ')', sep=''), 11)\n",
    "sp.c <- ss.scatter.plot(subset(data, Area.Type == 'C'), \n",
    "                        paste('Fig. 4C: Area Type C (R-Squared = ', round(summary(c.lm)$r.squared, 3), ')', sep=''), 11)\n",
    "sp.d <- ss.scatter.plot(subset(data, Area.Type == 'D'), \n",
    "                        paste('Fig. 4D: Area Type D (R-Squared = ', round(summary(d.lm)$r.squared, 3), ')', sep=''), 11)\n",
    "\n",
    "# Display the four plots within a single grid.\n",
    "grid.arrange(sp.a, sp.b, sp.c, sp.d, nrow=2) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The R-squared values are nearly 1.0 for all area types, indicating a near perfect correspondance between calculated and assigned scores. This indicates that the Smart Scale methodology was consistenly applied accross projects."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4.2 Synopsis and Conclusions\n",
    "\n",
    "In this section we examined whether the Smart Scale scoring methodology was consistently applied to the projects in the dataset by independently implementing the calculation method specified in the technical guide and comparing the results to the assigned Smart Scale scores. The degree of agreement was determined using linear regression fitting the assigned score to the calculated score. The results suggest that the methodology was applied consistently and in conformance to the November 2017 technical guide.\n",
    "\n",
    "## V. The Effect of Prioritizing Accessibility on Fredericksburg Projects\n",
    "\n",
    "The Fredericksburg area is unique in that a large percentage of its workforce commutes outside the region. Accordingly, in this section we will examine the effect of two hypothetical scenarios which give priority to the Accessibility factor in calculating the Smart Scale Scores of projects in the Fredericksburg district. Specifically, three scenarios will be compared: 1) the current SSS calculation method as specified in the November 2017 technical guide, 2) a SSS calculation method that give priority to the Accessibility factor, and 3) a SSS calculation method that gives priority to the Accessibility factor while simultaneously increasing the Access to Jobs metric by a factor of 1.5 for projects in the Fredericksburg region (explained further below). In each case we will look at the number of Fredericksburg projects appearing in different percentile ranges by Smart Scale Score. These scenarios will be referred to as 1, 2, and 3, respectively.\n",
    "\n",
    "The Accessibility factor is composed of the weighted average of three measures: 1) Access to Jobs (60%), 2) Access to Jobs for Disadvantaged Populations (20%), and 3) Access to Multimodal Choices (20%) (Nov. 2017 Technical Guide, p. 29). In Appendix C of the technical guide, the \"Access to Jobs\" measure is defined as the average increase in accessibility to jobs within 45 minutes. According to a recent study, approximately 37% of the workforce in Fredericksburg commutes outside the region, with more than 50% of these commuters commuting for an hour or more [1]. The vast majority of these commute to Northern VA, DC, Richmond, and Maryland - nearly all of which are beyond 45 minutes. A more appropriate measure for Fredericksburg would thus be to use a 75-minute reach or greater. Since there is no way to know how this reach would change the Access to Jobs without the GIS tool, an estimated increase of 1.5 times the current numbers will be used in scenario 3.\n",
    "\n",
    "To begin, we inspect the number of Fredericksburg projects in different areas:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Area.Type</th><th scope=col>Count</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>A </td><td>14</td></tr>\n",
       "\t<tr><td>B </td><td> 0</td></tr>\n",
       "\t<tr><td>C </td><td> 0</td></tr>\n",
       "\t<tr><td>D </td><td>11</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|ll}\n",
       " Area.Type & Count\\\\\n",
       "\\hline\n",
       "\t A  & 14\\\\\n",
       "\t B  &  0\\\\\n",
       "\t C  &  0\\\\\n",
       "\t D  & 11\\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "Area.Type | Count | \n",
       "|---|---|---|---|\n",
       "| A  | 14 | \n",
       "| B  |  0 | \n",
       "| C  |  0 | \n",
       "| D  | 11 | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "  Area.Type Count\n",
       "1 A         14   \n",
       "2 B          0   \n",
       "3 C          0   \n",
       "4 D         11   "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fb <- subset(data, District=='Fredericksburg')\n",
    "dt5 <- data.frame(table(fb$Area.Type))\n",
    "colnames(dt5) <- c('Area.Type', 'Count')\n",
    "dt5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Thus, Fredericksburg has 14 projects in area type A and 11 in D for a total of 25 projects. In the technical guide, Table 4.2 on p. 36, the factor weights by area type are as follows:\n",
    "\n",
    "![alt text](./table-4.2.png)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since all Fredericksburg projects are within area types A and D, we will reassign priorities as follows. In type A, we will give Accessibility a weight of 45% and Congestion Mitigation a weight of 15%. In type D we will give Accessibility a weight of 35% and Economic Development a weight of 15%.\n",
    "\n",
    "Below we will recalculate the Smart Scale Scores and percentile ranges under each assumption and see how this changes the distribution of Fredericksburg projects. \n",
    "\n",
    "### 5.1 Calculation of Smart Scale Scores Under Scenarios 2 and 3\n",
    "\n",
    "We begin by calculating the new SSS for scenarios 2 and 3 below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "# ------ Recalculate for scenario 2\n",
    "# Recalculate Project Benefit with revised weightings and add the new column to the data.\n",
    "data['Calc.Project.Benefit.Score.Scenario.2'] <- calc.proj.value(data, weight.a=c(0.15, 0.05, 0.45, 0.05, 0.1, 0.2),\n",
    "                                                                       weight.d=c(0.10, 0.15, 0.35, 0.3, 0.1))\n",
    "# Add a new SSS column under the scenario 2 assumptions\n",
    "data['Calc.SMART.SCALE.Score.Scenario.2'] <- (1e7 * data[['Calc.Project.Benefit.Score.Scenario.2']] / \n",
    "                                              data[['SMART.SCALE.Request']]) \n",
    "\n",
    "# ------ Recalculate for scenario 3\n",
    "# First, adjust the Fredericksburg Access.to.Jobs values by the designated factor.\n",
    "a2j.factor <- 1.5\n",
    "for(i in 1:nrow(data)) {\n",
    "    if(data[i, 'District'] == 'Fredericksburg') {data[i, 'Access.to.Jobs'] = a2j.factor * data[i, 'Access.to.Jobs']}\n",
    "}    \n",
    "# Second, renormalize the Access.to.Jobs column.\n",
    "data$Access.to.Jobs <- 100 * data$Access.to.Jobs / max(data$Access.to.Jobs)\n",
    "    \n",
    "# Third, recalculate the Accessibility factor\n",
    "data$Accessibility.Score <- weighted.column.average(data, c(0.6, 0.2, 0.2), \n",
    "                                                    c('Access.to.Jobs', 'Disadvantaged.Access.to.Jobs', \n",
    "                                                      'Multimodal.Access.Score'))\n",
    "    \n",
    "# Fourth, add a new column to data with our calculated project benefit scores for scenario 3\n",
    "data['Calc.Project.Benefit.Score.Scenario.3'] <- calc.proj.value(data, weight.a=c(0.15, 0.05, 0.45, 0.05, 0.1, 0.2),\n",
    "                                                                       weight.d=c(0.10, 0.15, 0.35, 0.3, 0.1))\n",
    "\n",
    "# Fifth, add a new SSS column under the scenario 3 assumptions\n",
    "data['Calc.SMART.SCALE.Score.Scenario.3'] <- (1e7 * data[['Calc.Project.Benefit.Score.Scenario.3']] / \n",
    "                                              data[['SMART.SCALE.Request']]) \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5.2 Comparison of Scenarios\n",
    "\n",
    "Below we perform the comparisons of scenarios. We first look at the number of Fredericksburg projects falling into different Smart Scale Score percentile bands in each scenario. By ranking the discrete percentile bands from 1 being the lowest to 8 being the highest, we then compare the average percentile ranking of Smart Scale scores for the Fredericksburg projects under each scenario."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Percentile_Rank</th><th scope=col>Percentile_Band</th><th scope=col>Scenario_1_Fredericksburg_Project_Count</th><th scope=col>Scenario_2_Fredericksburg_Project_Count</th><th scope=col>Scenario_3_Fredericksburg_Project_Count</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>1      </td><td>0-5%   </td><td>0      </td><td>0      </td><td>0      </td></tr>\n",
       "\t<tr><td>2      </td><td>5-10%  </td><td>2      </td><td>2      </td><td>1      </td></tr>\n",
       "\t<tr><td>3      </td><td>10-25% </td><td>5      </td><td>5      </td><td>4      </td></tr>\n",
       "\t<tr><td>4      </td><td>25-50% </td><td>7      </td><td>8      </td><td>9      </td></tr>\n",
       "\t<tr><td>5      </td><td>50-75% </td><td>5      </td><td>4      </td><td>5      </td></tr>\n",
       "\t<tr><td>6      </td><td>75-90% </td><td>3      </td><td>3      </td><td>2      </td></tr>\n",
       "\t<tr><td>7      </td><td>90-95% </td><td>2      </td><td>2      </td><td>3      </td></tr>\n",
       "\t<tr><td>8      </td><td>95-100%</td><td>1      </td><td>1      </td><td>1      </td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllll}\n",
       " Percentile\\_Rank & Percentile\\_Band & Scenario\\_1\\_Fredericksburg\\_Project\\_Count & Scenario\\_2\\_Fredericksburg\\_Project\\_Count & Scenario\\_3\\_Fredericksburg\\_Project\\_Count\\\\\n",
       "\\hline\n",
       "\t 1         & 0-5\\%    & 0         & 0         & 0        \\\\\n",
       "\t 2         & 5-10\\%   & 2         & 2         & 1        \\\\\n",
       "\t 3         & 10-25\\%  & 5         & 5         & 4        \\\\\n",
       "\t 4         & 25-50\\%  & 7         & 8         & 9        \\\\\n",
       "\t 5         & 50-75\\%  & 5         & 4         & 5        \\\\\n",
       "\t 6         & 75-90\\%  & 3         & 3         & 2        \\\\\n",
       "\t 7         & 90-95\\%  & 2         & 2         & 3        \\\\\n",
       "\t 8         & 95-100\\% & 1         & 1         & 1        \\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "Percentile_Rank | Percentile_Band | Scenario_1_Fredericksburg_Project_Count | Scenario_2_Fredericksburg_Project_Count | Scenario_3_Fredericksburg_Project_Count | \n",
       "|---|---|---|---|---|---|---|---|\n",
       "| 1       | 0-5%    | 0       | 0       | 0       | \n",
       "| 2       | 5-10%   | 2       | 2       | 1       | \n",
       "| 3       | 10-25%  | 5       | 5       | 4       | \n",
       "| 4       | 25-50%  | 7       | 8       | 9       | \n",
       "| 5       | 50-75%  | 5       | 4       | 5       | \n",
       "| 6       | 75-90%  | 3       | 3       | 2       | \n",
       "| 7       | 90-95%  | 2       | 2       | 3       | \n",
       "| 8       | 95-100% | 1       | 1       | 1       | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "  Percentile_Rank Percentile_Band Scenario_1_Fredericksburg_Project_Count\n",
       "1 1               0-5%            0                                      \n",
       "2 2               5-10%           2                                      \n",
       "3 3               10-25%          5                                      \n",
       "4 4               25-50%          7                                      \n",
       "5 5               50-75%          5                                      \n",
       "6 6               75-90%          3                                      \n",
       "7 7               90-95%          2                                      \n",
       "8 8               95-100%         1                                      \n",
       "  Scenario_2_Fredericksburg_Project_Count\n",
       "1 0                                      \n",
       "2 2                                      \n",
       "3 5                                      \n",
       "4 8                                      \n",
       "5 4                                      \n",
       "6 3                                      \n",
       "7 2                                      \n",
       "8 1                                      \n",
       "  Scenario_3_Fredericksburg_Project_Count\n",
       "1 0                                      \n",
       "2 1                                      \n",
       "3 4                                      \n",
       "4 9                                      \n",
       "5 5                                      \n",
       "6 2                                      \n",
       "7 3                                      \n",
       "8 1                                      "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th scope=col>Scenario</th><th scope=col>Average_Fredericksburg_Percentile_Rank</th><th scope=col>Average_Fredericksburg_Smart_Scale_Score</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><td>1       </td><td>4.48    </td><td>8.902320</td></tr>\n",
       "\t<tr><td>2       </td><td>4.44    </td><td>8.645798</td></tr>\n",
       "\t<tr><td>3       </td><td>4.64    </td><td>8.861798</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lll}\n",
       " Scenario & Average\\_Fredericksburg\\_Percentile\\_Rank & Average\\_Fredericksburg\\_Smart\\_Scale\\_Score\\\\\n",
       "\\hline\n",
       "\t 1        & 4.48     & 8.902320\\\\\n",
       "\t 2        & 4.44     & 8.645798\\\\\n",
       "\t 3        & 4.64     & 8.861798\\\\\n",
       "\\end{tabular}\n"
      ],
      "text/markdown": [
       "\n",
       "Scenario | Average_Fredericksburg_Percentile_Rank | Average_Fredericksburg_Smart_Scale_Score | \n",
       "|---|---|---|\n",
       "| 1        | 4.48     | 8.902320 | \n",
       "| 2        | 4.44     | 8.645798 | \n",
       "| 3        | 4.64     | 8.861798 | \n",
       "\n",
       "\n"
      ],
      "text/plain": [
       "  Scenario Average_Fredericksburg_Percentile_Rank\n",
       "1 1        4.48                                  \n",
       "2 2        4.44                                  \n",
       "3 3        4.64                                  \n",
       "  Average_Fredericksburg_Smart_Scale_Score\n",
       "1 8.902320                                \n",
       "2 8.645798                                \n",
       "3 8.861798                                "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# For a list of ranges and a vector, this function tabulates the counts by range.\n",
    "tabulate.by.ranges <- function(v, ranges) {\n",
    "  counts <- c()\n",
    "  ranges <- sort(ranges)\n",
    "  i <- 1\n",
    "  for(i in 1:length(ranges)) {\n",
    "    if(i == 1) {\n",
    "      counts[1] <- sum(v <= ranges[1])\n",
    "    }\n",
    "    else {\n",
    "      counts[i] <- sum((v > ranges[i - 1]) & (v <= ranges[i]))\n",
    "    }\n",
    "  }\n",
    "  counts\n",
    "}\n",
    "\n",
    "# Compare the scenarios.\n",
    "sss1.ord <- sort(data$SMART.SCALE.Score)\n",
    "sss2.ord <- sort(data$Calc.SMART.SCALE.Score.Scenario.2)\n",
    "sss3.ord <- sort(data$Calc.SMART.SCALE.Score.Scenario.3)\n",
    "percentiles <- c(0.05, 0.1, 0.25, 0.5, 0.75, 0.9, 0.95, 1) \n",
    "cut.scores.1 <- sapply(percentiles, function(p){sss1.ord[ceiling(p * nrow(data))]})\n",
    "cut.scores.2 <- sapply(percentiles, function(p){sss2.ord[ceiling(p * nrow(data))]})\n",
    "cut.scores.3 <- sapply(percentiles, function(p){sss3.ord[ceiling(p * nrow(data))]})\n",
    "\n",
    "fb <- subset(data, District=='Fredericksburg')\n",
    "\n",
    "counts.1 <- tabulate.by.ranges(fb$SMART.SCALE.Score, cut.scores.1)\n",
    "counts.2 <- tabulate.by.ranges(fb$Calc.SMART.SCALE.Score.Scenario.2, cut.scores.1)\n",
    "counts.3 <- tabulate.by.ranges(fb$Calc.SMART.SCALE.Score.Scenario.3, cut.scores.3)\n",
    "\n",
    "perc.range <-  c('0-5%', '5-10%', '10-25%', '25-50%', \n",
    "                 '50-75%', '75-90%', '90-95%', '95-100%')\n",
    "\n",
    "scenarios <- data.frame(Percentile_Rank=1:length(perc.range),\n",
    "                        Percentile_Band=perc.range,\n",
    "                        Scenario_1_Fredericksburg_Project_Count=counts.1,\n",
    "                        Scenario_2_Fredericksburg_Project_Count=counts.2,\n",
    "                        Scenario_3_Fredericksburg_Project_Count=counts.3)\n",
    "scenarios\n",
    "\n",
    "avg.perc.rank <- data.frame(Scenario=1:3,\n",
    "                    Average_Fredericksburg_Percentile_Rank=c(sum(scenarios$Percentile_Rank * counts.1) / sum(counts.1),\n",
    "                                                             sum(scenarios$Percentile_Rank * counts.2) / sum(counts.2),\n",
    "                                                             sum(scenarios$Percentile_Rank * counts.3) / sum(counts.3)),\n",
    "                            Average_Fredericksburg_Smart_Scale_Score=c(mean(fb$SMART.SCALE.Score),\n",
    "                                                                       mean(fb$Calc.SMART.SCALE.Score.Scenario.2),\n",
    "                                                                       mean(fb$Calc.SMART.SCALE.Score.Scenario.3)))\n",
    "avg.perc.rank"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the first table above we see the counts of projects within the Fredericksburg district arranged by percentile under scenarios 1 (current), 2 (prioritization of Accessibility factor for area types A and D), and 3 (prioritization of Accessibility factor for area types A and D + 1.5x increased \"Access to Jobs\" measures for Fredericksburg district). In the second table we see the average Fredericksburg percentile rank and average Fredericksburg Smart Scale Scores for each scenario. Two observations are apparent:\n",
    "\n",
    "1. There is no major difference between scenarios 1 and 2. \n",
    "\n",
    "2. There is a modest upward movement in Fredericksburg project percentile rank in scenario 3 (which prioritizes the Accessibility factor and also increases the weighting factor by 1.5), with the 25-50% and 90-95% bands seeing increases in the number of Fredericksburg projects. Of the 25 Fredericksburg projects in the dataset, two projects moved up into the 25-50% band and one more moved up into 90-95% band. \n",
    "\n",
    "### 5.3 Synopsis and Conclusions\n",
    "\n",
    "The Smart Scale methodology (as specified in the Nov. 2017 technical guide) places a relatively low priority on the Accessibility factor in area types A and D and only considers a 45-minute radius. This does not entirely reflect the unique commuter conditions found in the Fredericksburg region. The value of some projects in this region may accordingly be underestimated and, depending on the cut score, these projects may be edged out of funding.\n",
    "\n",
    "## VI. Conclusions\n",
    "\n",
    "This study performed a statistical analysis on a set of 404 transportation projects scored using the Virginia Smart Scale methodology. Descriptive statistics were first provided which characterize the data over several project features. An analysis was then performed to see whether the Smart Scale methodology was applied consistently to these projects, which entailed replicating the score calculations as specified in the November 2017 Technical Guide and comparing them to the assigned scores. The results indicate that the scoring methodolgy was applied consistently to the projects. Finally, the study considered whether increasing the emphasis on Accessibility in the scoring calculation would increase the competitiveness of projects within the Fredericksburg area. The results indicate that increasing the emphasis on Accessibility and considering a larger radius would lead to small increases in average percentile ranking for some Fredericksburg projects. Accordingly, because Fredericksburg has such a large percentage of its workforce which commutes outside the area, the impact of these projects may be underestimated by using the current Smart Scale methodology.   \n",
    "\n",
    "## References\n",
    "\n",
    "[1] Hansen, B.A. (2016). \"Fredericksburg Region Commuter Workforce Skills Study\", University of Mary Washington study prepared for the Fredericksburg Regional Alliance and the Fredericksburg Regional Chamber of Commerce, Accessed at [http://www.umw.edu/news/wp-content/uploads/sites/7/2016/05/Report-on-Commuters-FINAL_mm.pdf](http://www.umw.edu/news/wp-content/uploads/sites/7/2016/05/Report-on-Commuters-FINAL_mm.pdf)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R",
   "language": "R",
   "name": "ir"
  },
  "language_info": {
   "codemirror_mode": "r",
   "file_extension": ".r",
   "mimetype": "text/x-r-source",
   "name": "R",
   "pygments_lexer": "r",
   "version": "3.4.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}