library(dslabs)
library(tidyverse)
library(caret)
library(gridExtra)
The iris data is a widely used dataset in data analysis examples. It includes four botanical measurements related to three flower species: setosa, versicolor and virginica.
data(iris)
names(iris)
levels(iris$Species)
The following graphics show the three species and their distribution of the sepal and petal lengths and widths.
options(repr.plot.width=14, repr.plot.height=5)
grid.arrange(ggplot(data = iris, aes(Sepal.Length, Sepal.Width, col=Species)) + geom_point(),
ggplot(data = iris, aes(Petal.Length, Petal.Width, col=Species)) + geom_point(),
ncol = 2)
We will remove the setosa species and study the classification between versicolor and virginica species for a first approach. Later we will work on the entire set.
iris <- iris[-which(iris$Species=='setosa'),] %>%
mutate(Species = droplevels(Species))
y <- iris$Species
We create the train and test partitions.
test_index <- createDataPartition(y, times=1, p=0.4, list=FALSE)
test <- iris[test_index,]
train <- iris[-test_index,]
We start this study by classifying the species based on only one feature.
accuracies <- apply(train[,-5], 2, function(feature) {
cutoff <- seq(range(feature)[1], range(feature)[2], 0.1)
sapply(cutoff, function(x) {
y_hat <- ifelse(feature > x, "virginica", "versicolor")
mean(y_hat == train$Species)
})
})
sapply(accuracies, max)
Among the four features, the one that showed the best cutoff, with the highest accuracy, was the Petal Width.
So we will use the Petal Width to classify and calculate the accuracy over the test set.
feature <- which.max(sapply(accuracies, max))
cutoffs <- seq(range(train[,feature])[1], range(train[,feature])[2], 0.1)
best_cutoff <- cutoffs[which.max(accuracies[[feature]])]
best_cutoff
y_hat = factor(ifelse(test[,feature] > best_cutoff, "virginica", "versicolor"))
confusionMatrix(data = y_hat, reference = test$Species)
Confusion Matrix and Statistics
Reference
Prediction versicolor virginica
versicolor 18 2
virginica 2 18
Accuracy : 0.9
95% CI : (0.7634, 0.9721)
No Information Rate : 0.5
P-Value [Acc > NIR] : 9.285e-08
Kappa : 0.8
Mcnemar's Test P-Value : 1
Sensitivity : 0.90
Specificity : 0.90
Pos Pred Value : 0.90
Neg Pred Value : 0.90
Prevalence : 0.50
Detection Rate : 0.45
Detection Prevalence : 0.50
Balanced Accuracy : 0.90
'Positive' Class : versicolor
confusionMatrix(data = y_hat, reference = test$Species)$overall["Accuracy"]
cat("F1 Score:", F_meas(data = y_hat, reference = test$Species))
F1 Score: 0.9
The classification model we trained is in fact Linear Discriminant Analysis (LDA) applied over one predictor. By applying LDA, it has the same confusion matrix.
train_lda <- train(Species ~ Petal.Width, method = "lda", data = train)
train_lda$finalModel
Call:
lda(x, grouping = y)
Prior probabilities of groups:
versicolor virginica
0.5 0.5
Group means:
Petal.Width
versicolor 1.316667
virginica 2.030000
Coefficients of linear discriminants:
LD1
Petal.Width 4.613785
y_hat_lda <- predict(train_lda, test)
confusionMatrix(data = y_hat_lda, reference = test$Species)
Confusion Matrix and Statistics
Reference
Prediction versicolor virginica
versicolor 18 2
virginica 2 18
Accuracy : 0.9
95% CI : (0.7634, 0.9721)
No Information Rate : 0.5
P-Value [Acc > NIR] : 9.285e-08
Kappa : 0.8
Mcnemar's Test P-Value : 1
Sensitivity : 0.90
Specificity : 0.90
Pos Pred Value : 0.90
Neg Pred Value : 0.90
Prevalence : 0.50
Detection Rate : 0.45
Detection Prevalence : 0.50
Balanced Accuracy : 0.90
'Positive' Class : versicolor
lda_acc <- confusionMatrix(data = y_hat_lda, reference = test$Species)$overall["Accuracy"]
acc_results <- tibble(method = "2 classes: LDA: 1 feature", acc = lda_acc)
acc_results %>% knitr::kable()
|method | acc| |:-------------------------|---:| |2 classes: LDA: 1 feature | 0.9|
Using a single feature and optimizing the cutoff as we did on our training data can lead to overfitting.
We can explore the data, analysing each feature against the others.
options(repr.plot.width=10, repr.plot.height=10)
plot(iris,pch=21,bg=iris$Species)
By analyzing previous figure, we see that Petal.Length and Petal.Width in combination could potentially be more informative than either feature alone.
We will optimize the cutoffs separately for Petal.Length and Petal.Width in the train set and predict with Petal.Length and Petal.Width cutoffs together.
petal_length_cutoffs <- seq(range(train[,3])[1], range(train[,3])[2], 0.1)
petal_width_cutoffs <- seq(range(train[,4])[1], range(train[,4])[2], 0.1)
length_accuracies <- sapply(petal_length_cutoffs, function(x) {
y_hat <- ifelse(train$Petal.Length > x, "virginica", "versicolor")
mean(y_hat == train$Species)
})
cutoff_length <- petal_length_cutoffs[which.max(length_accuracies)]
width_accuracies <- sapply(petal_width_cutoffs, function(x) {
y_hat <- ifelse(train$Petal.Width > x, "virginica", "versicolor")
mean(y_hat == train$Species)
})
cutoff_width <- petal_width_cutoffs[which.max(width_accuracies)]
y_hat <- factor(ifelse(test$Petal.Length > cutoff_length |
test$Petal.Width > cutoff_width, "virginica", "versicolor"))
confusionMatrix(data = y_hat, reference = test$Species)
Confusion Matrix and Statistics
Reference
Prediction versicolor virginica
versicolor 17 0
virginica 3 20
Accuracy : 0.925
95% CI : (0.7961, 0.9843)
No Information Rate : 0.5
P-Value [Acc > NIR] : 9.733e-09
Kappa : 0.85
Mcnemar's Test P-Value : 0.2482
Sensitivity : 0.8500
Specificity : 1.0000
Pos Pred Value : 1.0000
Neg Pred Value : 0.8696
Prevalence : 0.5000
Detection Rate : 0.4250
Detection Prevalence : 0.4250
Balanced Accuracy : 0.9250
'Positive' Class : versicolor
confusionMatrix(data = y_hat, reference = test$Species)$overall["Accuracy"]
cat("F1 Score:", F_meas(data = y_hat, reference = test$Species))
F1 Score: 0.9189189
We can see that we achieved a better accuracy over the model. If we apply LDA over 2 predictors, the accuracy achieved is quite similar, but with a better distribution of sensitivity and specificity.
train_lda <- train(Species ~ Petal.Length + Petal.Width, method = "lda", data = train)
train_lda$finalModel
Call:
lda(x, grouping = y)
Prior probabilities of groups:
versicolor virginica
0.5 0.5
Group means:
Petal.Length Petal.Width
versicolor 4.27 1.316667
virginica 5.57 2.030000
Coefficients of linear discriminants:
LD1
Petal.Length 0.950446
Petal.Width 3.321780
y_hat_lda <- predict(train_lda, test)
confusionMatrix(data = y_hat_lda, reference = test$Species)
Confusion Matrix and Statistics
Reference
Prediction versicolor virginica
versicolor 18 2
virginica 2 18
Accuracy : 0.9
95% CI : (0.7634, 0.9721)
No Information Rate : 0.5
P-Value [Acc > NIR] : 9.285e-08
Kappa : 0.8
Mcnemar's Test P-Value : 1
Sensitivity : 0.90
Specificity : 0.90
Pos Pred Value : 0.90
Neg Pred Value : 0.90
Prevalence : 0.50
Detection Rate : 0.45
Detection Prevalence : 0.50
Balanced Accuracy : 0.90
'Positive' Class : versicolor
lda_acc <- confusionMatrix(data = y_hat_lda, reference = test$Species)$overall["Accuracy"]
acc_results <- bind_rows(acc_results,
tibble(method="2 classes: LDA: 2 features", acc = lda_acc))
acc_results %>% knitr::kable()
|method | acc| |:--------------------------|---:| |2 classes: LDA: 1 feature | 0.9| |2 classes: LDA: 2 features | 0.9|
We now apply Quadratic Discrimant Analysis, using all predictors, in search for a better accuracy.
train_qda <- train(Species ~ ., method = "qda", data = train)
train_qda$finalModel
Call:
qda(x, grouping = y)
Prior probabilities of groups:
versicolor virginica
0.5 0.5
Group means:
Sepal.Length Sepal.Width Petal.Length Petal.Width
versicolor 5.956667 2.763333 4.27 1.316667
virginica 6.616667 3.056667 5.57 2.030000
y_hat_qda <- predict(train_qda, test)
confusionMatrix(data = y_hat_qda, reference = test$Species)
Confusion Matrix and Statistics
Reference
Prediction versicolor virginica
versicolor 19 0
virginica 1 20
Accuracy : 0.975
95% CI : (0.8684, 0.9994)
No Information Rate : 0.5
P-Value [Acc > NIR] : 3.729e-11
Kappa : 0.95
Mcnemar's Test P-Value : 1
Sensitivity : 0.9500
Specificity : 1.0000
Pos Pred Value : 1.0000
Neg Pred Value : 0.9524
Prevalence : 0.5000
Detection Rate : 0.4750
Detection Prevalence : 0.4750
Balanced Accuracy : 0.9750
'Positive' Class : versicolor
qda_acc <- confusionMatrix(data = y_hat_qda, reference = test$Species)$overall["Accuracy"]
acc_results <- bind_rows(acc_results,
tibble(method="2 classes: QDA: all features", acc = qda_acc))
acc_results %>% knitr::kable()
|method | acc| |:----------------------------|-----:| |2 classes: LDA: 1 feature | 0.900| |2 classes: LDA: 2 features | 0.900| |2 classes: QDA: all features | 0.975|
We now proceed to apply the LDA and QDA models with the entire dataset and outcomes.
data(iris)
test_index <- createDataPartition(iris$Species, times=1, p=0.4, list=FALSE)
test <- iris[test_index,]
train <- iris[-test_index,]
train_lda <- train(Species ~ ., method = "lda", data = train)
train_lda$finalModel
Call:
lda(x, grouping = y)
Prior probabilities of groups:
setosa versicolor virginica
0.3333333 0.3333333 0.3333333
Group means:
Sepal.Length Sepal.Width Petal.Length Petal.Width
setosa 4.973333 3.453333 1.440000 0.2433333
versicolor 5.933333 2.743333 4.263333 1.3066667
virginica 6.466667 2.966667 5.426667 2.0100000
Coefficients of linear discriminants:
LD1 LD2
Sepal.Length 1.018446 -0.6747386
Sepal.Width 1.349276 2.3200825
Petal.Length -2.395358 -1.2192081
Petal.Width -2.452692 4.2048620
Proportion of trace:
LD1 LD2
0.983 0.017
y_hat_lda <- predict(train_lda, test)
confusionMatrix(data = y_hat_lda, reference = test$Species)
Confusion Matrix and Statistics
Reference
Prediction setosa versicolor virginica
setosa 20 0 0
versicolor 0 18 0
virginica 0 2 20
Overall Statistics
Accuracy : 0.9667
95% CI : (0.8847, 0.9959)
No Information Rate : 0.3333
P-Value [Acc > NIR] : < 2.2e-16
Kappa : 0.95
Mcnemar's Test P-Value : NA
Statistics by Class:
Class: setosa Class: versicolor Class: virginica
Sensitivity 1.0000 0.9000 1.0000
Specificity 1.0000 1.0000 0.9500
Pos Pred Value 1.0000 1.0000 0.9091
Neg Pred Value 1.0000 0.9524 1.0000
Prevalence 0.3333 0.3333 0.3333
Detection Rate 0.3333 0.3000 0.3333
Detection Prevalence 0.3333 0.3000 0.3667
Balanced Accuracy 1.0000 0.9500 0.9750
lda_acc <- confusionMatrix(data = y_hat_lda, reference = test$Species)$overall["Accuracy"]
acc_results <- bind_rows(acc_results,
tibble(method="3 classes: LDA", acc = lda_acc))
acc_results %>% knitr::kable()
|method | acc| |:----------------------------|---------:| |2 classes: LDA: 1 feature | 0.9000000| |2 classes: LDA: 2 features | 0.9000000| |2 classes: QDA: all features | 0.9750000| |3 classes: LDA | 0.9666667|
train_qda <- train(Species ~ ., method = "qda", data = train)
train_qda$finalModel
Call:
qda(x, grouping = y)
Prior probabilities of groups:
setosa versicolor virginica
0.3333333 0.3333333 0.3333333
Group means:
Sepal.Length Sepal.Width Petal.Length Petal.Width
setosa 4.973333 3.453333 1.440000 0.2433333
versicolor 5.933333 2.743333 4.263333 1.3066667
virginica 6.466667 2.966667 5.426667 2.0100000
y_hat_qda <- predict(train_qda, test)
confusionMatrix(data = y_hat_qda, reference = test$Species)
Confusion Matrix and Statistics
Reference
Prediction setosa versicolor virginica
setosa 20 0 0
versicolor 0 18 0
virginica 0 2 20
Overall Statistics
Accuracy : 0.9667
95% CI : (0.8847, 0.9959)
No Information Rate : 0.3333
P-Value [Acc > NIR] : < 2.2e-16
Kappa : 0.95
Mcnemar's Test P-Value : NA
Statistics by Class:
Class: setosa Class: versicolor Class: virginica
Sensitivity 1.0000 0.9000 1.0000
Specificity 1.0000 1.0000 0.9500
Pos Pred Value 1.0000 1.0000 0.9091
Neg Pred Value 1.0000 0.9524 1.0000
Prevalence 0.3333 0.3333 0.3333
Detection Rate 0.3333 0.3000 0.3333
Detection Prevalence 0.3333 0.3000 0.3667
Balanced Accuracy 1.0000 0.9500 0.9750
qda_acc <- confusionMatrix(data = y_hat_qda, reference = test$Species)$overall["Accuracy"]
acc_results <- bind_rows(acc_results,
tibble(method="3 classes: QDA", acc = qda_acc))
acc_results %>% knitr::kable()
|method | acc| |:----------------------------|---------:| |2 classes: LDA: 1 feature | 0.9000000| |2 classes: LDA: 2 features | 0.9000000| |2 classes: QDA: all features | 0.9750000| |3 classes: LDA | 0.9666667| |3 classes: QDA | 0.9666667|
Since this dataset is small and has few predictors, LDA and QDA was fast and showed high accuracy. The predictors also looked to follow multivariate normal distribution, which is required for these models.