--- title: 'Classification: Naive Bayes and SVM' author: "Torben" date: "August 24, 2018" output: pdf_document: default html_document: default --- ```{r setup, include=FALSE} knitr::opts_chunk$set(echo = TRUE, message = FALSE) library(knitr) library(tidyverse) theme_set(theme_bw()) ``` # Classification Today we will discuss different types of classification methods. One is based on a probabilistic argument, the other on separable *hyperplanes* (which is just a fancy word for a plane that divides the feature space). In classification we are set in a supervised learning situation. That is, we have some training data for which we have some features $X$ and a response $Y$. In classification we always has that $Y$ is categorical - i.e. it has a certain number of levels. They may be ordered, but ordinal situations is out of our scope. Typically, we begin with the case of assuming that $Y$ is binary -- it one has two levels, e.g. $\{0, 1\}$ or $\{-1, 1\}$. Most classification methods have been derived from this simplest case, and is then extended to deal with more than two classes. ## Model based approach For the model based approaches, the typical procedure is to compute a *posterior probability* for a class given the features. That is, $P(Y = y\mid X = x)$, where $Y$ is the *stochastic variable* and $y$ is *a* state, e.g $y = 0$ or $y = 1$ in the binary case. Similarly, $X$ is the random quantity and for a specific observation $X$ has the value(s) $x$. A posterior probability is the probability of $Y$ *given* we have seen the feature information $X$. The *a priori probability*, or short *prior*, reflects the believe we have about $Y$ *before* we see the features. ### Bayes classifier A classifier that assigns the class to be the most probable is called a *Bayes classifier*: $$\hat{k} = \arg\max_k P(Y=k \mid X=x_0),$$ for some features $x_0$. ### Binary case First, we observe that $P(Y = 0 \mid X) + P(Y = 1 \mid X) = 1$, which implies $P(Y = 0 \mid X) = 1 - P(Y = 1 \mid X)$. Hence, we can focus just on $P(Y = 1 \mid X)$, say. $$P(Y = 1 \mid X) = \frac{P(Y = 1, X)}{P(X)} = \frac{P(Y = 1, X)}{P(X)}\frac{P(Y = 1)}{P(Y = 1)} = \frac{P(X\mid Y = 1)P(Y = 1)}{P(X)}$$ For the binary case, we have that a Bayes classifier will assign the class to be $1$ if $P(Y = 1\mid X)>P(Y=0\mid X)$, because we know they sum to 1 (there are only two outcomes). Furthermore, we can write $$\frac{P(Y = 1\mid X)}{P(Y=0\mid X)}>1$$ And since the above factorisation of $P(Y=1\mid X)$ also holds for $P(Y=0\mid X)$ we have $$\frac{P(Y = 1\mid X)}{P(Y=0\mid X)} = \frac{P(X\mid Y = 1)}{P(X\mid Y = 0)}\frac{P(Y = 1)}{P(Y = 0)}$$ Hence, we just need a model for $P(X\mid Y = y)$ and then assign a *prior* probability to $P(Y = 1)$. Often this is (maybe confusingly) denoted $\pi = P(Y = 1)$ with $P(Y=0) = 1-\pi$. Note, $\pi$ is in this case **not** `r pi`. ## Naive Bayes Modelling $P(X\mid Y = y)$ may not be easy as $X$ can high (extremely) high-dimensional. However, one model assumption that simplifies this dramatically is that of (conditional) independence: $$P(X\mid Y = y) = P(X_1\mid Y=y)P(X_2\mid Y=y)\cdots P(X_p\mid Y=y) = \prod_{i=1}^p P(X_i\mid Y=y)$$ This turns a complicated model for $P(X\mid Y=y)$ into a product of simple models -- one for each $X_i$. Including this in the expression above yields for $X = x_0$, $$\frac{P(Y = 1\mid X = x_0)}{P(Y=0\mid X = x_0)} = \frac{\pi}{1-\pi}\prod_{i=1}^p\frac{P(X_i = x_{0i}\mid Y = 1)}{P(X_i = x_{0i}\mid Y = 0)}$$ Hence, if the ratio above is greater than 1, we classify a new observation $x_0$ as $\hat{Y} = 1$. Otherwise, classify $\hat{Y} = 0$. The model has several advantages: * Simple to estimate parameters (no need for iterative procedures) * Insensitive to missing data (the term just disappears) * Works for $n\ll p$ ### `naiveBayes` in R The package `e1071` contains several methodologies, including `naiveBayes` and `svm` (which we will use here). The `naiveBayes` assumes that for numerical features, $x_i$ follows a normal distribution, with group specific mean and variance: $\mu_y$ and $\sigma_y^2$. For categorical cases we simply tabulate and estimate the associated probabilities from the counts. #### Example: `Titanic` A small example with the survival of Titanic passengers. ```{r} Titanic_tbl <- Titanic %>% as_tibble() %>% mutate(n = as.integer(n)) %>% mutate_if(is.character, factor) %>% mutate_at(c("Age", "Survived"), funs(fct_rev)) Titanic_tbl %>% kable() ``` Fitting a `naiveBayes` (we use the tabular form of the data here) ```{r} library(e1071) nb_titanic <- naiveBayes(Survived ~ ., data = Titanic) ``` The estimates of *prior* and $P(x_i\mid y)$ ```{r} nb_titanic$apriori nb_titanic$tables ``` Which variables are more informative to predict survival? ```{r} ## Class: p_class <- nb_titanic$tables$Class p_class["Yes",]/p_class["No",] # Sex p_sex <- nb_titanic$tables$Sex p_sex["Yes",]/p_sex["No",] # Age p_age <- nb_titanic$tables$Age p_age["Yes",]/p_age["No",] ``` #### Exercise: ```{r} iris %>% ggplot(aes(x = Petal.Length, y = Petal.Width, colour = Species)) + geom_point() ``` * Use `naiveBayes` to classify the `iris` flowers into their three classes. * What information does the `$tables` from the fit contain? * Try to identify the most informative variables for predicting the classes. * Are the assumptions about normality satisfied for each of $x_i\mid y$? * Can you visualise the decision boundaries in the `Petal`-plane (i.e `x = Petal.Length, y = Petal.Width`)? *Tip:* Create a grid (cf below), make the prediction for *each* point in the grid and plot this. ```{r} iris_petal_grid <- iris %>% expand( Petal.Length = seq(min(Petal.Length), max(Petal.Length), len = 100), Petal.Width = seq(min(Petal.Width), max(Petal.Width), len = 100) ) %>% mutate(Sepal.Length = NA, Sepal.Width = NA) ``` # Support Vector Machines: SVM See slides `day-5-SVM.pdf` for introduction to SVMs ```{r} iris_svm_linear <- svm(Species ~ ., data = iris, cost = 1, kerner = "linear") iris_svm_linear summary(iris_svm_linear) plot(iris_svm_linear, formula = Petal.Length ~ Petal.Width, data = iris, slice = list(Sepal.Width = 3, Sepal.Length = 4)) iris_svm_radial <- svm(Species ~ ., data = iris, cost = 1, kerner = "radial") plot(iris_svm_radial, formula = Petal.Length ~ Petal.Width, data = iris, slice = list(Sepal.Width = 3, Sepal.Length = 4)) iris_svm_poly_3 <- svm(Species ~ ., data = iris, cost = 100, kerner = "polynomial", degree = 5) plot(iris_svm_poly_3, formula = Petal.Length ~ Petal.Width, data = iris, slice = list(Sepal.Width = 3, Sepal.Length = 4)) ``` ## Tuning ```{r} ncol_data <- ncol(iris) iris_radial_tune <- tune.svm(Species ~ ., data = iris, kerner = "radial", cost = 10^(0:3), gamma = 1/(ncol_data*c(0.5,1,2))) plot(iris_radial_tune) summary(iris_radial_tune) ``` # Topics not covered * ROC curves - Used to decide on an optimal threshold value. That is, the threshold of 0.5 may not be optimal * $k$-Nearest Neighbouhrs - A 'simple' technique where a *test sample* is classified based by a majority vote among its $k$ closest data points in the *training data* - This is called on *online* or *lazy* learner as it does not fit a model to data, but uses all the training data for each new classification task. - Typically cross-validation is used to decide on $k$ * Imbalanced training data case - When samples of one type is much more dominant in the training data. One approach is to use weights for methods that allows/incorporates this. * Linear (LDA) and quadratic discriminant analysis (QDA) - The predecessors of SVM - Relies on an assumption of multivariate normality of the data (given the class) * Ordinal classification - Situations where the levels of $Y$ has an ordering to them, e.g. *low*, *mid* and *high* - See the `ordinal` package for regression methods to deal with this type of analysis. - See the `rpartOrdinal` or `rpartScore` for extentions to `rpart` for classification trees with ordinal responses. - See `glmnetcr` for a `glmnet` like approach to ordinal response prediction