Pactical Machine Learning Project
Background
Using devices such as Jawbone Up, Nike FuelBand, and Fitbit it is now possible to collect a large amount of data about personal activity relatively inexpensively. These type of devices are part of the quantified self movement â a group of enthusiasts who take measurements about themselves regularly to improve their health, to find patterns in their behavior, or because they are tech geeks. One thing that people regularly do is quantify how much of a particular activity they do, but they rarely quantify how well they do it. In this project, the goal will be to use data from accelerometers on the belt, forearm, arm, and dumbell of 6 participants. They were asked to perform barbell lifts correctly and incorrectly in 5 different ways. More information is available from the website here: http://groupware.les.inf.puc-rio.br/har (see the section on the Weight Lifting Exercise Dataset).
The goal of this project is to predict the manner in which they did the exercise by building a machine learning algorithm. This is the “classe” variable in the training set.
Data Preparation
The following libraries shall be used for this task:
library(caret)## Loading required package: lattice## Loading required package: ggplot2library(rattle)## Rattle: A free graphical interface for data mining with R.
## Version 4.1.0 Copyright (c) 2006-2015 Togaware Pty Ltd.
## Type 'rattle()' to shake, rattle, and roll your data.library(rpart)For reproducibility we shall set the seed for the experiment.
set.seed(9876)The training data for this project are available here:
pml.training <- "http://d396qusza40orc.cloudfront.net/predmachlearn/pml-training.csv"The test data are available here:
pml.testing <- "http://d396qusza40orc.cloudfront.net/predmachlearn/pml-testing.csv"Read the training and testing data and specify the NA strings
training <- read.csv(url(pml.training), na.strings=c("NA","#DIV/0!",""))
testing <- read.csv(url(pml.testing), na.strings=c("NA","#DIV/0!",""))Partition Data
Partition the training data into 60:40 ratios
inTrain <- createDataPartition(y=training$classe, p=0.6, list=FALSE)
theTraining <- training[inTrain, ]
dim(theTraining)## [1] 11776 160theTesting <- training[-inTrain, ]
dim(theTesting)## [1] 7846 160Clean Data
We need to reduce the number of features by removing variables with nearly zero variance, the first ID column which is not useful for prediction, as well as variables that are almost always NA.
Firstly, clean the partitioned training data.
## remove near zero variables
nzv <- nearZeroVar(theTraining)
theTraining <- theTraining[, -nzv]
## exclude the first ID column from the training
theTraining<- theTraining[c(-1)]
## remove variables containing mostly NAs (90%)
theNAs <- sapply(theTraining, function(x) mean(is.na(x))) > 0.90
theTraining <- theTraining[, theNAs==F]Secondly, clean the partitioned testing data
## remove near zero variables
nzv <- nearZeroVar(theTesting)
theTesting <- theTesting[, -nzv]
## exclude the ID column from the testing
theTesting<- theTesting[c(-1)]
## remove variables containing mostly NAs (90%)
theNAs <- sapply(theTesting, function(x) mean(is.na(x))) > 0.90
theTesting <- theTesting[, theNAs==F]Make sure the data are of the same types
## coerce data into same data type
for (i in 1:length(theTesting) ) {
for(j in 1:length(theTraining)) {
if( length( grep(names(theTraining[i]), names(theTesting)[j]) ) ==1) {
class(theTesting[j]) <- class(theTraining[i])
}
}
}Modeling and Prediction
We shall try 2 machine learning algorithms, namely the classification tree and random forests, and then choose the one that gives a better accuracy.
Model using classification tree and then predict:
## train data and create model using rpart, then print fancy plot
modelFit <- train(classe ~ ., method="rpart", data=theTraining)
print(modelFit$finalModel)## n= 11776
##
## node), split, n, loss, yval, (yprob)
## * denotes terminal node
##
## 1) root 11776 8428 A (0.28 0.19 0.17 0.16 0.18)
## 2) roll_belt< 129.5 10708 7393 A (0.31 0.21 0.19 0.18 0.11)
## 4) pitch_forearm< -33.95 981 8 A (0.99 0.0082 0 0 0) *
## 5) pitch_forearm>=-33.95 9727 7385 A (0.24 0.23 0.21 0.2 0.12)
## 10) cvtd_timestamp02/12/2011 13:33>=0.5 765 201 A (0.74 0.26 0 0 0) *
## 11) cvtd_timestamp02/12/2011 13:33< 0.5 8962 6892 B (0.2 0.23 0.23 0.22 0.13)
## 22) magnet_dumbbell_z< -34.5 2869 1750 A (0.39 0.31 0.08 0.19 0.038)
## 44) raw_timestamp_part_1< 1.322838e+09 696 31 A (0.96 0.04 0.0043 0 0) *
## 45) raw_timestamp_part_1>=1.322838e+09 2173 1321 B (0.21 0.39 0.1 0.24 0.05) *
## 23) magnet_dumbbell_z>=-34.5 6093 4269 C (0.11 0.2 0.3 0.23 0.17) *
## 3) roll_belt>=129.5 1068 33 E (0.031 0 0 0 0.97) *fancyRpartPlot(modelFit$finalModel)
## perform prediction and show confusion matrix
prediction1 <- predict(modelFit, theTesting, type = "raw")
confusionMatrix(prediction1, theTesting$classe)## Confusion Matrix and Statistics
##
## Reference
## Prediction A B C D E
## A 1422 152 1 0 0
## B 299 584 180 386 83
## C 470 782 1187 900 695
## D 0 0 0 0 0
## E 41 0 0 0 664
##
## Overall Statistics
##
## Accuracy : 0.4916
## 95% CI : (0.4805, 0.5027)
## No Information Rate : 0.2845
## P-Value [Acc > NIR] : < 2.2e-16
##
## Kappa : 0.3637
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: A Class: B Class: C Class: D Class: E
## Sensitivity 0.6371 0.38472 0.8677 0.0000 0.46047
## Specificity 0.9727 0.85019 0.5605 1.0000 0.99360
## Pos Pred Value 0.9029 0.38120 0.2942 NaN 0.94184
## Neg Pred Value 0.8708 0.85207 0.9525 0.8361 0.89105
## Prevalence 0.2845 0.19347 0.1744 0.1639 0.18379
## Detection Rate 0.1812 0.07443 0.1513 0.0000 0.08463
## Detection Prevalence 0.2007 0.19526 0.5141 0.0000 0.08985
## Balanced Accuracy 0.8049 0.61745 0.7141 0.5000 0.72703The accuracy is only 49% and the split is not balanced. Therefore we need to try another algorithm.
Model using random forests with 10 cross validation and then predict:
## 10 cross validation to train data and create model using random forest
fitControl <- trainControl(method="cv", number=10, verboseIter=F)
fit <- train(classe ~ ., method ="rf", trControl=fitControl, data=theTraining)## Loading required package: randomForest## randomForest 4.6-12## Type rfNews() to see new features/changes/bug fixes.##
## Attaching package: 'randomForest'## The following object is masked from 'package:ggplot2':
##
## marginprint(fit$finalModel)##
## Call:
## randomForest(x = x, y = y, mtry = param$mtry)
## Type of random forest: classification
## Number of trees: 500
## No. of variables tried at each split: 40
##
## OOB estimate of error rate: 0.08%
## Confusion matrix:
## A B C D E class.error
## A 3348 0 0 0 0 0.0000000000
## B 1 2277 1 0 0 0.0008775779
## C 0 1 2050 3 0 0.0019474197
## D 0 0 2 1926 2 0.0020725389
## E 0 0 0 0 2165 0.0000000000## perform prediction and show confusion matrix
prediction2 <- predict(fit, theTesting, type = "raw")
confusionMatrix(prediction2, theTesting$classe)## Confusion Matrix and Statistics
##
## Reference
## Prediction A B C D E
## A 2232 0 0 0 0
## B 0 1518 2 0 0
## C 0 0 1364 4 0
## D 0 0 2 1281 0
## E 0 0 0 1 1442
##
## Overall Statistics
##
## Accuracy : 0.9989
## 95% CI : (0.9978, 0.9995)
## No Information Rate : 0.2845
## P-Value [Acc > NIR] : < 2.2e-16
##
## Kappa : 0.9985
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: A Class: B Class: C Class: D Class: E
## Sensitivity 1.0000 1.0000 0.9971 0.9961 1.0000
## Specificity 1.0000 0.9997 0.9994 0.9997 0.9998
## Pos Pred Value 1.0000 0.9987 0.9971 0.9984 0.9993
## Neg Pred Value 1.0000 1.0000 0.9994 0.9992 1.0000
## Prevalence 0.2845 0.1935 0.1744 0.1639 0.1838
## Detection Rate 0.2845 0.1935 0.1738 0.1633 0.1838
## Detection Prevalence 0.2845 0.1937 0.1744 0.1635 0.1839
## Balanced Accuracy 1.0000 0.9998 0.9982 0.9979 0.9999The accuracy is 99.89%, thus the out-of-sample error is 0.08%. Since randomForest gives better accuracy, write the output of the 20 test cases to files.
Validation
thePrediction <- predict(fit, testing, type = "raw")
pml_write_files = function(x){
n = length(x)
for(i in 1:n){
filename = paste0("test_",i,".txt")
write.table(x[i],file=filename,quote=FALSE,row.names=FALSE,col.names=FALSE)
}
}
pml_write_files(thePrediction)