IntrusionDetection.Rmd

---
title: "IntrusionDetection"
authors: Jacob Tarnow, Simon Kwong, Daniel Silva
output: html_document
---
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see <http://rmarkdown.rstudio.com>.

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

## Install Dependencies

```{r}
install.packages("multicore")
install.packages("data.table")
install.packages("foreach")
install.packages("doParallel")
install.packages("knitr")
```

## Call all libraries

```{r}
library(ISLR)
library(MASS)
library(e1071)
library(foreach)
library(parallel)
library(iterators)
library(doParallel)
library(data.table)
library(tree)
```

## Print the runtime.
The code snippet below is used to call at any point in our code base to calculate the amount of time it takes to run various blocks.

```{r}
print.time <- function(starttime) {
  starttime = proc.time() - starttime
  if (starttime[3] > 60) {
    cat(sprintf("%f minutes \n", starttime[3]/60))
  } else {
    cat(sprintf("%f seconds \n", starttime[3]))
  }
  rm(starttime)
}

```

## Pre-Processing: Multi-Core and Data Setup
Due to the fact that our dataset is so highly skewed and noisy from KDDCup, we need to create a new sample from pre-processing. In order to maximize our computation time and efficiency, we used multi-core processing. The below code sets everything up for processing.

```{r}
# Setup MultiProcessing ---------------------------------------------------
# This is setting up for multicore processing

packages <- c("foreach", "doParallel")

num_cores <- detectCores() - 2
# This Will require permissions (So let's not use)
# core_cluster <- makeCluster(num_cores)
registerDoParallel(cores = num_cores)
# getDoParWorkers()


# Data Setup --------------------------------------------------------------
# setwd("path to directory where the data files are stored")
getwd() # to check if the working directory is the correct path

data.load.time <- proc.time()

kddcup.data = fread("kddcup.data.csv")
kddcup.data.ten.percent = fread("kddcup.data_10_percent.csv")
kddcup.testdata = fread("kddcup.testdata.csv")
kddcup.testdata.ten.percent = fread("kddcup.testdata_10_percent.csv")

print.time(data.load.time)
rm(data.load.time)
```

## Creating the Target Variable
As our project determines intrustion detection and whether or not a TCP connection is going to be "Good" or "Bad" (binary), we created a list of the bad connections and used this to determine the representations of our target variable.

```{r}
bad_connections <- c(
  "back.",
  "buffer_overflow.",
  "ftp_write.",
  "guess_passwd.",
  "imap.",
  "ipsweep.",
  "land.",
  "loadmodule.",
  "multihop.",
  "neptune.",
  "nmap.",
  "perl.",
  "phf.",
  "pod.",
  "portsweep.",
  "rootkit.",
  "satan.",
  "smurf.",
  "spy.",
  "teardrop.",
  "warezclient.",
  "warezmaster."
)
#good_connections <- c()

# Good = 0 || Bad = 1
kddcup.data <- cbind(kddcup.data, rep(0, dim(kddcup.data)[1]))
kddcup.data.ten.percent <- cbind(kddcup.data.ten.percent, rep(0, dim(kddcup.data.ten.percent)[1]))

# Adding 2 Columns for the test data set
kddcup.testdata <- cbind(kddcup.testdata, rep(as.factor(NA), dim(kddcup.testdata)[1]))
kddcup.testdata <- cbind(kddcup.testdata, rep(as.factor(NA), dim(kddcup.testdata)[1]))

kddcup.testdata.ten.percent <- cbind(kddcup.testdata.ten.percent, rep(as.factor(NA), dim(kddcup.testdata.ten.percent)[1]))
kddcup.testdata.ten.percent <- cbind(kddcup.testdata.ten.percent, rep(as.factor(NA), dim(kddcup.testdata.ten.percent)[1]))

```

## Labeling our Data
Most of the data that was gathered from KDD Cup was un-labeled. In order to determine which predictors are significant and aid in plotting later on, we needed to re-label our data set. From our column names, then we can recall the target variable and tune it.

```{r}
column_names <- c(
  "duration",
  "protocol_type",
  "service",
  "flag",
  "src_bytes",
  "dst_bytes",
  "land",
  "wrong_fragment",
  "urgent",
  "hot",
  "num_failed_logins",
  "logged_in",
  "num_compromised",
  "root_shell",
  "su_attempted",
  "num_root",
  "num_file_creations",
  "num_shells",
  "num_access_files",
  "num_outbound_cmds",
  "is_host_login",
  "is_guest_login",
  "count",
  "srv_count",
  "serror_rate",
  "srv_serror_rate",
  "rerror_rate",
  "srv_rerror_rate",
  "same_srv_rate",
  "diff_srv_rate",
  "srv_diff_host_rate",
  "dst_host_count",
  "dst_host_srv_count",
  "dst_host_same_srv_rate",
  "dst_host_diff_srv_rate",
  "dst_host_same_src_port_rate",
  "dst_host_srv_diff_host_rate",
  "dst_host_serror_rate",
  "dst_host_srv_serror_rate",
  "dst_host_rerror_rate",
  "dst_host_srv_rerror_rate",
  "connection_type",
  "access_type"
)

colnames(kddcup.data) = column_names
colnames(kddcup.data.ten.percent) = column_names
colnames(kddcup.testdata) = column_names
colnames(kddcup.testdata.ten.percent) = column_names

# Sets all accesses as good, replace as bad if in bad_connections
# Good = 0 || Bad = 1
kddcup.data$access_type = 0
kddcup.data.ten.percent$access_type = 0

kddcup.data$access_type[kddcup.data$connection_type %in% bad_connections] = 1
kddcup.data.ten.percent$access_type[kddcup.data.ten.percent$connection_type %in% bad_connections] = 1

train = kddcup.data

##### Remove Random Period From The Dataset
# connection_types <- data.frame(kddcup.data.ten.percent$connection_type)
# kddcup.data.ten.percent$connection_type = substr(connection_types, 0, nchar(connection_types)-1)
```

## Feature Selection
Through our feature selection we were able to go through the data set and determine the unnecessary features that are not needed in our models. We determined these through histograms of each feature and gathering the density of each. Based off of the density of each feature, we can see which ones are more appropriate for our models. We figured that features with a good distribution of density would work well, instead of those that have more of a slow differential for density.

```{r}
# Feature Selection -------------------------------------------------------

names = names(train)
classes = sapply(train, class)
par(mfrow=c(1,1))
numeric.classes = c("integer", "numeric")
for (name in names[classes %in% numeric.classes]) {
  print(name)
  hist.predictor = hist(train[,name], xlab=name)
  hist.predictor
  print(summary(hist.predictor$density))
}

train = kddcup.data
train.src_bytes = train[train$src_bytes < 1100, ]


# based off of the density skews, we can drop the following predictors
# Here are just a handful of the unnecessary features

unnecessary.features = c(
  "connection_type", 
  "access_type", 
  "duration", 
  "dst_bytes", 
  "land",
  "wrong_fragment",
  "urgent",
  "hot",
  "num_failed_logins",
  "num_compromised",
  "root_shell",
  ""
)

# Statistically Significant After Sampling
used.feaures = c(
  "src_bytes",
  "logged_in",
  "flag",
  "num_root",
  "num_file_creations",
  "count",
  "srv_count",
  "serror_rate",
  "srv_serror_rate",
  "rerror_rate",
  "srv_rerror_rate",
  "same_srv_rate",
  "diff_srv_rate",
  "srv_diff_host_rate",
  "dst_host_count",
  "dst_host_srv_count",
  "dst_host_same_srv_rate",
  "dst_host_diff_srv_rate",
  "dst_host_same_src_port_rate",
  "dst_host_srv_diff_host_rate",
  "dst_host_serror_rate",
  "dst_host_srv_serror_rate",
  "dst_host_rerror_rate",
  "dst_host_srv_rerror_rate"
)
```


## Dimensionality Reduction
Our dimensionality reduction uses the unique function to gather unique data points from our most significant features: connection type, protocol type, service, and flag.

```{r}
# Dimensionality Reduction ------------------------------------------------
set.seed(666)
train = kddcup.data

connections = unique(train$connection_type)
protocols = unique(train$protocol_type)
services = unique(train$service)
flags = unique(train$flag)
```

## Quick and Dirty Sampling
Quick way to gain a decent sample. But later on you can grab the sample that we used. Takes roughly 3hrs to run on a decent CPU.

```{r}
# Quick And Dirty Sampling ------------------------------------------------
# Just Takes a 10 percent sample
# Still Very Skewed

new.train = train[0, ]
sample.size = 20000

for (connection in 1:length(connections)) {
  print(connections[connection])
  train.sample = kddcup.data[kddcup.data$connection_type == connections[connection], ]
  print(nrow(train.sample))
  if (nrow(train.sample) < sample.size) {
    train.sample = train.sample[sample(nrow(train.sample), size = nrow(train.sample), replace = nrow(train.sample) < sample.size), ]
  } else {
    train.sample = train.sample[sample(nrow(train.sample), size = 0.1*nrow(train.sample), replace = nrow(train.sample) < sample.size), ]
  }
  new.train = rbind(new.train, train.sample)
  rm(train.sample)
}

train = new.train
rm(new.train)

new.train = train[0, ]
sample.size = 20000

for (service in 1:length(services)) {
  print(services[service])
  train.sample = kddcup.data[kddcup.data$service == services[service], ]
  print(nrow(train.sample))
  if (nrow(train.sample) < sample.size) {
    train.sample = train.sample[sample(nrow(train.sample), size = nrow(train.sample), replace = nrow(train.sample) < sample.size), ]
  } else {
    train.sample = train.sample[sample(nrow(train.sample), size = 0.1*nrow(train.sample), replace = nrow(train.sample) < sample.size), ]
  }
  new.train = rbind(new.train, train.sample)
  rm(train.sample)
}

train = rbind(train, new.train)
rm(new.train)
```

## Very Normalized Sampling
Run this at your own risk!!
```{r}
# Very Normalized Sampling -------------------------------------------------
### Only Run This If You Have Too Much Time On Your Hands
# Creates a sample where All possible combinations of connections are trained
# In order for us not to have to run this everytime, we have exported the data in a .RData file that we can load
# when we restart R Studio's workspace

train = kddcup.data
new.train = train[0, ]
sample.size = 2000

print.sample <- function(c, p, s, f) {
  cat(sprintf("\"%s\" \"%s\" \"%s\" \"%s\"\n", connections[c], protocols[p], services[s], flags[f]))
}

# Undersampling | Oversampling
sample.train <- function(c, p, s, f) {
  train.sample = train[train$connection_type == connections[c] & 
                                              train$protocol_type == protocols[p] &
                                              train$service == services[s] &
                                              train$flag == flags[f], ]
  if (nrow(train.sample) < sample.size) {
    train.sample = train.sample[sample(nrow(train.sample), size = nrow(train.sample), replace = TRUE), ]
  } else {
    train.sample = train.sample[sample(nrow(train.sample), size = 0.1*nrow(train.sample), FALSE), ]
  }
  
  return(train.sample)
}

nestedfor.time <- proc.time()

# Sequential version
for (connection in 1:length(connections)) {
  for (protocol in 1:length(protocols)) {
    for (service in 1:length(services)) {
      for (flag in 1:length(flags)) {
        print.sample(connection, protocol, service, flag)
        train.sample = sample.train(connection, protocol, service, flag)
        new.train = rbind(new.train, train.sample)
        rm(train.sample)
      }
    }
  }
}

# Multicore Version
# Only works on certain OS systems that support multicore processing

# foreach(connection=1:length(connections)) %:%
#   foreach(protocol=1:length(protocols)) %:%
#   foreach(service=1:length(services)) %:%
#   foreach(flag=1:length(flags)) %dopar% {
#       print.sample(connection, protocol, service, flag)
#       train.sample = sample.train(connection, protocol, service, flag)
#       new.train = rbind(new.train, train.sample)
#       rm(train.sample)
#   }

print.time(nestedfor.time)
rm(nestedfor.time)

train = new.train
rm(new.train)
```

## Basic Training and Testing Sample
Please load the file - newData.RData for the remainder of our code.

```{r}
# Basic Trainng and Testing Set Sample ------------------------------------------
train = kddcup.data.ten.percent
test = kddcup.testdata.ten.percent
```

## Logistic Regression
Logistic Regression is a method for fitting a regression curve. It fits b0 and b1, the regression coefficients. The curve for this will not be linear. When we first ran this in October, the average was roughly 83%, after running this with the commented out predictors we gained an accuracy of 98%. We felt like that is too good to be true, possibly due to over-fitting. We ran again with the following predictors: src_bytes and logged_in and gained an accuracy of 94.2%. We believe that these two predictors help determine the optimum  access type as both the src bytes of the TCP connection and whether the user is logged in, go hand-in-hand with determining whether the connection is good or bad (an attack).

```{r}
# Logistic Regression -----------------------------------------------------
glm.fit.time <- proc.time()

#glm.fit = glm(access_type~
#              +flag
#              +src_bytes
#              +logged_in
#              +num_root
#              +num_file_creations
#              +count
#              +srv_count
#              +serror_rate
#              +srv_serror_rate
#              +rerror_rate
#              +srv_rerror_rate
#              +same_srv_rate
#              +diff_srv_rate
#              +srv_diff_host_rate
#              +dst_host_count
#              +dst_host_srv_count
#              +dst_host_same_srv_rate
#              +dst_host_diff_srv_rate
#              +dst_host_same_src_port_rate
#              +dst_host_srv_diff_host_rate
#              +dst_host_serror_rate
#              +dst_host_srv_serror_rate
#              +dst_host_rerror_rate
#              +dst_host_srv_rerror_rate
#              -access_type
#              -connection_type, data=train, family=binomial)
glm.fit = glm(access_type~src_bytes+logged_in, data=train, family=binomial)
summary(glm.fit)
glm.probs = predict(glm.fit, newdata=kddcup.data.ten.percent, type = "response")
glm.pred = ifelse(glm.probs > 0.5, 1, 0)
glm.pred.accesses = kddcup.data.ten.percent$access_type

table(glm.pred, glm.pred.accesses)
mean(glm.pred == glm.pred.accesses)

print.time(glm.fit.time)
rm(glm.fit.time)
```

## Linear Discriminant Analysis
```{r}
# Linear Discriminant Analysis (LDA) --------------------------------------
lda.fit.time <- proc.time()

#lda.fit = lda(access_type~
#              +flag
#              +src_bytes
#              +logged_in
#              +num_root
#              +num_file_creations
#              +count
#              +srv_count
#              +serror_rate
#              +srv_serror_rate
#              +rerror_rate
#              +srv_rerror_rate
#              +same_srv_rate
#              +diff_srv_rate
#              +srv_diff_host_rate
#              +dst_host_count
#              +dst_host_srv_count
#              +dst_host_same_srv_rate
#              +dst_host_diff_srv_rate
#              +dst_host_same_src_port_rate
#              +dst_host_srv_diff_host_rate
#              +dst_host_serror_rate
#              +dst_host_srv_serror_rate
#              +dst_host_rerror_rate
#              +dst_host_srv_rerror_rate
#              -access_type
#              -connection_type, data=train, family=binomial)

lda.fit = lda(access_type~src_bytes+logged_in, data=train, family=binomial)
summary(lda.fit)
lda.pred = predict(lda.fit, newdata=kddcup.data.ten.percent, type = "response")
table(lda.pred$class, kddcup.data.ten.percent$access_type)
mean(lda.pred$class == kddcup.data.ten.percent$access_type)

print.time(lda.fit.time)
rm(lda.fit.time)
```

## Quadratic Discriminant Analysis
```{r}
# Quadratic Discriminant Analysis (QDA) -----------------------------------

qda.fit.time <- proc.time()

#qda.fit = qda(access_type~
#              +flag
#              +src_bytes
#              +logged_in
#              +num_root
#              +num_file_creations
#              +count
#              +srv_count
#              +serror_rate
#              +srv_serror_rate
#              +rerror_rate
#              +srv_rerror_rate
#              +same_srv_rate
#              +diff_srv_rate
#              +srv_diff_host_rate
#              +dst_host_count
#              +dst_host_srv_count
#              +dst_host_same_srv_rate
#              +dst_host_diff_srv_rate
#              +dst_host_same_src_port_rate
#              +dst_host_srv_diff_host_rate
#              +dst_host_serror_rate
#              +dst_host_srv_serror_rate
#              +dst_host_rerror_rate
#              +dst_host_srv_rerror_rate
#              -access_type
#              -connection_type, data=train, family=binomial)
qda.fit = qda(access_type~src_bytes+logged_in, data=train, family=binomial)
summary(qda.fit)
qda.pred = predict(qda.fit, newdata=kddcup.data.ten.percent)
table(qda.pred$class, kddcup.data.ten.percent$access_type)
mean(qda.pred$class == kddcup.data.ten.percent$access_type)

print.time(qda.fit.time)
rm(qda.fit.time)
```

## Including Plots

Due to our large dataset, plotting would crash our systems

You can also embed plots, for example:

```{r pressure, echo=FALSE}
plot(glm.pred)

pairs(train)
```

Note that the `echo = FALSE` parameter was added to the code chunk to prevent printing of the R code that generated the plot.