[R] Classification Problem : LDA(Linear Discriminant Analysis)

1. Discriminant Analysis

• Bayes Theorem : $Pr(Y=k|X=x)=\frac{{\pi }_k{f}_k(x)}{\sum _{l=1}^K{\pi }_l{f}_l(x)}$
• The density for $X$ in class $k$ : $f_k(x)=\frac{1}{\sqrt{2\pi }{\sigma }_k}{e}^{-\frac{(x-{\mu }_k{)}^2}{2{\sigma }_k^2}}$
• The prior probability for class $k$ : ${\pi }_k=Pr(Y=k)$
• If the distribution of the $X$ is approximately normal, LDA and QDA is more stable.

 

1.1 LDA(Linear Discriminant Analysis)

• Assume that $\sigma ={\sigma }_1={\sigma }_2=...{\sigma }_k$
• $f_k(x)$
  • p = 1 : $f_k(x)\to N({\mu }_k,\sigma )$
  • p > 1 : $f_k(x)\to N({\mu }_k,\sum )$
•  $p_k(x)=\frac{{\pi }_k\frac{1}{\sqrt{2\pi }\sigma }{e}^{-\frac{(x-{\mu }_k{)}^2}{2{\sigma }^2}}}{\sum _{l=1}^K{\pi }_l\frac{1}{\sqrt{2\pi }\sigma }{e}^{-\frac{(x-{\mu }_l{)}^2}{2{\sigma }^2}}}$
• To classify at the value $X=x$, see which of the $p_k(x)$ is largest.
• Because below term of $p_k(x)$ is same, we need to consider upper terms.
  • Discriminant score :
    • p = 1 : ${\sigma }_k(x)=x\frac{{\mu }_k}{{\sigma }^2}-\frac{{\mu }_k^2}{2{\sigma }^2}+log({\pi }_k)$
    • p > 1 : ${\sigma }_k(x)=x^T\sum ^{-1}{\mu }_k-\frac{1}{2}{\mu }_k^T\sum ^{-1}{\mu }_k+log({\pi }_k)$
  • ${\sigma }_k(x)$ is a linear function of $x$. 
•  Decision boundary at $x=\frac{\widehat{{\mu }_1}+\widehat{{\mu }_2}}{2}$
• Need to estimate : $({\mu }_1,...,{\mu }_k),({\pi }_1,...,{\pi }_k),\sum$

 

1.2 [Ex] Iris Data

Training model with LDA

 

# open the iris dataset 
data(iris)
str(iris)
summary(iris)
plot(iris[, -5], col=as.numeric(iris$Species) + 1)

# Apply LDA for iris data
library(MASS)
g <- lda(Species ~., data=iris)
plot(g)
plot(g, dimen=1)

 

 

Compute Missclassification Rate for training sets

 

# Compute misclassification error for training sets
pred <- predict(g)
table(pred$class, iris$Species)
mean(pred$class!=iris$Species)

 

 

We can get inferred probability 𝑝̂ 𝑘(𝑥) : predict(g, iris)$posterior[-tran, ].

 

Performance comparison between LDA vs Multinomial Regression

 

library(nnet)
set.seed(1234)
K <- 100
RES <- array(0, c(K, 2))
for (i in 1:K) {
    tran.num <- sample(nrow(iris), size=floor(nrow(iris)*2/3))
    tran <- as.logical(rep(0, nrow(iris)))
    tran[tran.num] <- TRUE
    g1 <- lda(Species ~., data=iris, subset=tran)
    g2 <- multinom(Species ~., data=iris, subset=tran, trace=FALSE)
    pred1 <- predict(g1, iris[!tran,])$class
    pred2 <- predict(g2, iris[!tran,])
    RES[i, 1] <- mean(pred1!=iris$Species[!tran])
    RES[i, 2] <- mean(pred2!=iris$Species[!tran])
}
apply(RES, 2, mean)

# 0.0254 0.0436

 

1.3 [Ex] Default Dataset

# Importing library 
library(ISLR)
data(Default)
attach(Default)
library(MASS)

# Training model 
g <- lda(default~., data=Default)
pred <- predict(g, default)
table(pred$class, default)
mean(pred$class!=default)

 

 

• False Positive rate : The fraction of negative examples that are classified as positive(0.22%)
• False Negative rate : The fraction of positive examples that are classified as negative(76.2%)
• If we classified the prior - always to class "No", then we would make 333/10000 errors, which is only 3.33%.