classification.Rmd

---
title: "Predicting 2016 Election Results in the USA"
author: "Ian Jeffries"
date: "5/16/2019"
output:
  pdf_document: default
  html_document: default
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

The following R Markdown outlines the steps used by Ian Jeffries to predict the 2016 election results in the USA 
by county, based on the demographics of that county. This is a supervised machine learning classification problem, 
and three models were tested to predict the election outcome: K Nearest Neighbor, Artifical Neural Networks, and
Decision Trees. 

## Install necessary packages

The following code will install the packages used in the project:

```{r packages, message=FALSE, warning=FALSE}
#list of packages used
packages <- c("dplyr", "tidyr", "ggplot2", "class", "rpart", "rpart.plot", "neuralnet", 
              "arules", "plyr", "mltools", "arulesViz", "plotly", "RCurl")

#check to see if package is already installed, if not, install
for(p in packages){
  if(!require(p, character.only = TRUE)) {
    install.packages(p)
    library(p, character.only = TRUE)
  } 
}
```
## Explore the Demographics Dataset

Import the data from my github page and explore.

```{r demographics}
#import our county features dataset
path <- "https://github.com/raw/ianjeffries/election-predictions/master/
data/county_facts.csv"

demographics <- read.csv(path, header=TRUE)

#See the names of the columns
names(demographics)
#See the structure of the dataset
str(demographics)

#find the state summary info and remove 
#(want only counties, this dataset includes state summaries)
demographics <- demographics[-which(demographics$state_abbreviation == ""), ]
```

## Add in Dictionary information

Add in the dictionary file to replace columm headers with readable values. 

```{r dictionary}
#pull in dictionary to get true names of the variables
path <- "https://github.com/raw/ianjeffries/election-predictions/master/
data/county_facts_dictionary.csv"

var_names <- read.csv(path, header = TRUE)

#view variable names
print(var_names)
```

## Pull in 2016 Election Results and Clean the Data

As of now, the US primarily operates in a two-party system. Because of this, only the 
democratic and republican candidates are of interest. 

```{r 2016_Results}
#add in the election results
path <- "https://github.com/raw/ianjeffries/election-predictions/master/
data/pres16results.csv"

election_results <- read.csv(path, header = TRUE)

#see the names of the columns
names(election_results)

#see structure of the dataset
str(election_results)

#drop columns we don't need
election_results <- election_results[-c(2, 5, 8)]

#filter by the republican and democratic candidates & remove null values
ER_clean <- election_results %>%
                filter(cand %in% c("Donald Trump", "Hillary Clinton") & county != "<NA>")

#change to wide dataset to merge the election results with the demographics dataframe
ER_clean <- spread(ER_clean, key = cand, value = votes)

#add in % won by
ER_clean$Percent_Rep <- round(ER_clean$`Donald Trump` / ER_clean$total_votes, 2)
ER_clean$Percent_Dem <- round(ER_clean$`Hillary Clinton` / ER_clean$total_votes, 2)
```

## Join the Demographic and Election Results Datasets

Now that cleanup is complete, the two datasets can be joined. 

```{r Joined_Data, warning=FALSE}
#join our election results with the demographics table
final_dataset <- left_join(demographics, select(ER_clean, Percent_Dem, 
                                                Percent_Rep, st, county), 
                           by = c("area_name" = "county", "state_abbreviation" = "st"))

#find number of counties with N/A election results
print(paste0("Number of counties with null election results: ", 
             nrow(final_dataset[which(complete.cases(final_dataset) == FALSE), ])))

#looks like these aren't counties (the majority) and 
#the N is small enough I feel comfortable removing them
final_dataset <- na.omit(final_dataset)
```

## Final Cleanup

The final cleanup steps are to remove irrelevant values and reassign column names for analysis. 
The data will also be normalized, with 0 or 1 values assigned to the target variable.

```{r cleanup}
#change any column in dataset based on % of the population to true percentage 
#(basically already normalized between 0 and 1)
final_dataset[,c(8:24, 28:29, 35, 41:46)] <- final_dataset[,c(8:24, 28:29, 35, 41:46)]/100

#drop clearly irrelevant values
final_dataset <- final_dataset[ , -c(5:8, 38)]

#assign headers to something more understandable, based on the dictionary dataframe
names(final_dataset) <- var_names$description[match(names(final_dataset), 
                                                    var_names$column_name)]

#assign labels that weren't in the dictionary
names(final_dataset[ , c(1, 2, 3, 50, 51)]) <- c("fips", "area_name", "state_abbreviation", 
                                                 "Percent_Dem", "Percent_Rep")

#normalize any data not described as a %
for (i in c(4, 21, 22, 23, 26, 27, 28, 29, 30, 32, 33, 
            34, 35, 42, 43, 44, 45, 46, 47, 48, 49)) {
  final_dataset[,i] <- ((final_dataset[,i] - min(final_dataset[,i])) / 
                          (max(final_dataset[,i]) - min(final_dataset[,i])) *
                          (1 - 0) + 0)
}

#add in 1 or 0 if they voted democratic or republican
final_dataset$Winner <- 0

for (i in 1:nrow(final_dataset)) {
    if (final_dataset[i, "percent_dem"] < .50) {
      final_dataset[i, "Winner"] <- 0 
      } else {
        final_dataset[i, "Winner"] <- 1 
      }}

#set winner column to factor for classification
final_dataset$Winner <- as.factor(final_dataset$Winner)
```

## Look for Variables to use in Classification

Rather than use all 45 variables in classification, each normalized variable was graphed against the percentage 
vote for the democratic candidate to see if there are any clear correlations. Any variables that didn't seem to 
hold a correlation were removed. 

```{r relevant_variables, fig.height=10, fig.width=9.5}
#create narrow dataset to understand relationship between attributes and voting preference
final_narrow <- gather(final_dataset, key = demographic_measure, value = stat, 4:49)

#create side by side scatterplots to look for possible correlations
ggplot(final_narrow, aes(x = percent_dem, y = stat)) +
  geom_point(alpha = .5) +
  geom_smooth() +
  facet_wrap(vars(demographic_measure), ncol = 5) +
  theme(strip.text = element_text(size = 5, face = "bold"))

#remove columns that don't seem to have a correlation to voting preference
final_dataset <- final_dataset[ , -c(1,5,7,10,12,13,16,28,30,34,35,
                                     37,39,42,43,44,45,47,48)]
```

## Classification Prep Work

Before the models can be created, training and test datasets must be set up. An accuracy matrix was 
also created to house the results from each model. 

```{r prep_work}
###### CREATE TRAINING AND TEST DATASETS ###### 

#take only the predictor values and class from our dataset
e_predictions <- final_dataset[ , c(3:30, 33)]

#set set to make our test data reproducable
set.seed(111692)

#get 80/20 sample of data
dist <- sample(2, size = nrow(e_predictions), replace = TRUE, prob = c(.8, .2))

#set test and training data
e_training <- e_predictions[dist == 1, ]
e_test <- e_predictions[dist == 2, ]

#show dimensions of training dataset
dim(e_training)

#show dimensions of test dataset
dim(e_test)

#create accuracy table for comparison
algorithm_results <- matrix(data = rep(0, 6), nrow = 2, ncol = 3, 
                            dimnames = list(c("Accuracy %:", "Precision %:"), 
                                            c("KNN", "DTREE", "ANN")))

#show accuracy table
print(algorithm_results)
```

## K-Nearest Neighbor Classification
```{r knn, fig.height=3.5, fig.width=5}
###### KNN CLASSIFICATION ######

#create KNN accuracy matrix to compare 15 K-Values
KNN_accuracy <- matrix(data= rep(0, 45), ncol = 3, dimnames = list(c(1:15), 
                                                                   c("K_Value", "Accuracy", 
                                                                              "Precision")))

#set each K value as a factor
KNN_accuracy[ , 1] <- as.factor(c(1:15))

#test multiple K values and store the results 
#(For KNN function, you need to remove the columns with actual classification)
for (i in 1:15) {
  KNN_results <- as.data.frame(knn(e_training[,1:28], e_test[,1:28], e_training$Winner, 
                                 k = i, prob = TRUE))
 
  #create a comparison table
  KNN_table <- table(KNN_results[,1], e_test$Winner)
  
  #find % accuracy
  KNN_accuracy[i , 2] <- round(((sum(diag(KNN_table)) / sum(KNN_table)) * 100), 2)
  
  #find % precision
  KNN_accuracy[i , 3] <- round((KNN_table[2, 2] / (KNN_table[2, 2] + 
                                                     KNN_table[2, 1]) * 100), 2)
}

#plot accuracy results
theme_update(plot.title = element_text(hjust = 0.5))
ggplot(as.data.frame(KNN_accuracy), aes(x = K_Value, y = Accuracy)) +
  geom_bar(stat = "identity", fill = "black") +
  ggtitle("K Value Accuracy Comparision") +
  scale_x_continuous(breaks = c(1:15)) +
  scale_y_continuous(limits = c(0, 100)) +
  geom_bar(data=subset(as.data.frame(KNN_accuracy), Accuracy==max(Accuracy)), 
           aes(K_Value, Accuracy),
           fill="green", stat="identity")

#plot precision results
ggplot(as.data.frame(KNN_accuracy), aes(x = K_Value, y = Precision)) +
  geom_bar(stat = "identity", fill = "black") +
  ggtitle("K Value Precision Comparision") +
  scale_x_continuous(breaks = c(1:15)) +
  scale_y_continuous(limits = c(0, 100)) +
  geom_bar(data=subset(as.data.frame(KNN_accuracy), Precision==max(Precision)), 
           aes(K_Value, Precision),
           fill="green", stat="identity")

#add to results table (highest k value at the time I ran)
algorithm_results[1:2, 1] <- KNN_accuracy[12, 2:3]
```

## Decision Tree
```{r dtree, fig.height=7, fig.width=9.5}
###### DECISION TREE CLASSIFICATION ######

#Train the decision tree
e_tree <- rpart(Winner ~ ., data = e_training)

#plot the branches of the tree 
rpart.plot(e_tree, type=1)

#returns a list of strings summarizing the branch path to each node
rpart.rules(e_tree, style = "tall") #useful to see what conditions predict a certain factor

#predict using that model on test data
tree_results <- as.data.frame(predict(e_tree, newdata = e_test, type = "class"))

#create a comparison table
tree_table <- table(tree_results[ , 1], e_test$Winner)

#create an accuracy table
tree_accuracy <- matrix(data= rep(0, 2), ncol = 2, dimnames = 
                          list("%",c("Accuracy", "Precision")))

#find accuracy
tree_accuracy[1, 1] <- round(((sum(diag(tree_table)) / sum(tree_table)) * 100), 2)

#find % precision
tree_accuracy[1 , 2] <- round((tree_table[2, 2] / (tree_table[2, 2] + 
                                                     tree_table[2, 1]) * 100), 2)

#add to results table
algorithm_results[1:2, 2] <- tree_accuracy
```

## Neural Networks

```{r neural-net, fig.height=3.5, fig.width=5, message=FALSE}
###### ANN CLASSIFICATION ######

#build the neural network 

#neural networks do not work for categorical variables (ie 1, 2 or 3 ranking) 
#To get around this we put each variable in it's own column, either a 1 or 0 
#for columns 1, 2 or 3 (use the model.matrix function to do this)

#first move categorical variables into their own columns
#(I believe this formula finds columns that are categorical classes)
e_training_m <- model.matrix(~0 + ., data = e_training)

#do the same for the test dataset
e_test_m <- model.matrix(~0 + ., data = e_test)

#create an accuracy table
ANN_accuracy <- matrix(data= rep(0, 39), ncol = 3, dimnames = list(c(3:15), 
                                                                   c("Nodes", "Accuracy", 
                                                                              "Precision")))

ANN_accuracy[ , 1] <- c(3:15)

#test multiple hidden nodes in our model
for (i in 3:15) {

  #build the artifical neural network using neuralnet function
  e_network <- neuralnet(Winner0 + Winner1 ~ Population_2014 + Percent_65_Years_and_Older + 
                           Percent_White_Alone + Percent_Black_or_African_American + 
                           Percent_Asian + Percent_Hispanic_or_Latino + 
                           Percent_White_Alone_not_Hispanic_or_Latino +
                           Percent_Foreign_born_persons + 
                           Percent_Language_other_than_English_spoken_at_home +
                           Percent_High_school_graduate_or_higher + 
                           Percent_Bachelors_degree_or_higher +
                           Percent_Veterans + Mean_travel_time_to_work_minutes + 
                           Housing_units + Percent_Homeownership_rate + 
                           Housing_units_in_multi_unit_structures + 
                           Median_value_of_owner_occupied_housing_units + 
                           Number_Households + Per_capita_money_income_in_past_12_months +
                           Percent_Persons_below_poverty_level + 
                           Private_nonfarm_establishments + Private_nonfarm_employment + 
                           Percent_Black_owned_firms + Percent_Asian_owned_firms + 
                           Percent_Hispanic_owned_firms + Percent_Women_owned_firms + 
                           Accommodation_and_food_services_sales_2007_1000 + 
                           Population_per_square_mile,
                         data = e_training_m,
                         hidden = i,
                         lifesign = "full",
                         linear.output = FALSE,
                         threshold = .1)
  
  #test the ANN on the test data
  ANN_results <- compute(e_network, e_test_m[ ,1:28]) 
                #^ computes result of each row, remove winner variables

  #extract the results
  prediction <- ANN_results$net.result
  
  #use max.col function to pull the original values from our test dataset into one column
  original_values <- max.col(e_test_m[ ,29:30]) - 1
  
  #use max.col function to get prediction values
  prediction_values <- max.col(prediction) - 1
  
  #create a dataframe with both values
  results <- data.frame(actual = original_values, prediction = prediction_values)
  
  #create a comparison table
  ANN_table <- table(results$prediction, e_test$Winner)
  
  #find accuracy
  ANN_accuracy[i-2, 2] <- round(((sum(diag(ANN_table)) / sum(ANN_table)) * 100), 2)
  
  #find % precision
  ANN_accuracy[i-2 , 3] <- round((ANN_table[2, 2] / (ANN_table[2, 2] + ANN_table[2, 1]) * 
                                    100), 2)
}

#plot the accuracy results
ggplot(as.data.frame(ANN_accuracy), aes(x = Nodes, y = Accuracy)) +
  geom_bar(stat = "identity", fill = "black") +
  ggtitle("Hidden Node Accuracy Comparision") +
  scale_x_continuous(breaks = c(3:15)) +
  scale_y_continuous(limits = c(0, 100)) +
  geom_bar(data=subset(as.data.frame(ANN_accuracy), Accuracy==max(Accuracy)), 
           aes(Nodes, Accuracy),
           fill="green", stat="identity")

#plot the precision results
ggplot(as.data.frame(ANN_accuracy), aes(x = Nodes, y = Precision)) +
  geom_bar(stat = "identity", fill = "black") +
  ggtitle("Hidden Node Precision Comparision") +
  scale_x_continuous(breaks = c(3:15)) +
  scale_y_continuous(limits = c(0, 100)) +
  geom_bar(data=subset(as.data.frame(ANN_accuracy), Precision==max(Precision)), 
           aes(Nodes, Precision),
           fill="green", stat="identity")

#add to results table
algorithm_results[1:2, 3] <- ANN_accuracy[1, 2:3]
print(algorithm_results)
```
## Model Comparison
```{r comparison, fig.height=5, fig.width=8}
###### COMPARISON OF CLASSIFICATION MODELS ######

#create final results (results kept changing as I used different work stations, 
#                       using values from last test for the write up)
final_results <- matrix(data= rep(0, 18), ncol = 3, dimnames = list(c(1:6), 
                                                                    c("Model", "Measure", 
                                                                      "Percentage")))
final_results <- as.data.frame(final_results)

#add results manually 
final_results$Model <- c("KNN", "DTREE", "ANN", "KNN", "DTREE", "ANN")

final_results$Measure <- c("Accuracy", "Accuracy", "Accuracy", 
                           "Precision", "Precision", "Precision")

final_results$Percentage <- c(94.97, 93.79, 94.23, 91.78, 79.12, 88.31)

#plot the accuracy results
ggplot(final_results, aes(x = Model, y = Percentage, fill = Measure)) +
  geom_bar(stat = "identity", position = position_dodge(width = .6), alpha = .8) +
  ggtitle("R Model Comparison") +
  scale_y_continuous(limits = c(0, 100)) + 
  geom_hline(yintercept = 94.97, linetype = "dashed", size = 1, color = "#E54445") + 
  geom_hline(yintercept = 91.78, linetype = "dashed", size = 1, color = "#5B94C2") +
  scale_fill_brewer(palette = "Set1")

#make narrow dataset even more narrow to plot by party
final_narrow2 <- gather(final_narrow, key = party, value = percent, 4:5)

#plot vote by population (to see if high populations voted democratic)
final_narrow2 %>%
  filter(demographic_measure == "Population_2014",
         party == "percent_dem") %>%
  ggplot(aes(x = percent, y = stat)) +
  geom_point(color = "#5B94C2") +
  geom_smooth() +
  xlab("Percent Voted Democratic") +
  ylab("Population Size (Normalized)") +
  ggtitle("Population by Percent Voted Democratic")

#plot top demographics
final_narrow2 %>%
  filter(demographic_measure == "Percent_White_Alone_not_Hispanic_or_Latino" |
         demographic_measure == "Percent_Black_or_African_American" |
         demographic_measure == "Housing_units_in_multi_unit_structures" |
         demographic_measure == "Per_capita_money_income_in_past_12_months" |
         demographic_measure == "Median_value_of_owner_occupied_housing_units") %>%
  ggplot(aes(x = percent, y = stat, color = party)) +
  geom_smooth() +
  facet_wrap(~demographic_measure, ncol = 2) +
  ggtitle("Voting Trends by Top Demographics") +
  xlab("Percent of Vote") +
  ylab("Normalized Demographic Value") +
  theme(strip.text = element_text(size=11), legend.position = c(0.75, 0.15)) + 
  scale_colour_manual(labels = c("Democrat", "Republican"), values = c("#5B94C2", 
                                                                       "#E54445")) +
  scale_x_continuous(limits = c(0, 1), breaks = c(0, 0.25, 0.5, 0.75, 1)) +
  scale_y_continuous(limits = c(0, 1), breaks = c(0, 0.25, 0.5, 0.75, 1))
```

## Conclusion

In conclusion, KNN was the best performing model to use when predicting election results by county. 
A full write-up of the results can be found on Ian Jeffries' github page.