knitr::opts_chunk$set(cache=TRUE, warning=FALSE, message=FALSE, fig.width=12)
library(lattice)
library(knitr)
library(caret)
## Loading required package: ggplot2
library(e1071)
library(randomForest)
## randomForest 4.6-12
## Type rfNews() to see new features/changes/bug fixes.

1. Introduction

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 is to predict the manner in which they did the exercise. This is the “classe” variable in the training set. You may use any of the other variables to predict with. You should create a report describing how you built your model, how you used cross validation, what you think the expected out of sample error is, and why you made the choices you did. You will also use your prediction model to predict 20 different test cases.

2. Data exploration and cleaning

2.1 Data load

In first place I load the train dataset and the submit_test will be uploaded for evaluation.

#download files
train_url <- "http://d396qusza40orc.cloudfront.net/predmachlearn/pml-training.csv"
download.file(url=train_url, destfile="training.csv")

test_url <- "http://d396qusza40orc.cloudfront.net/predmachlearn/pml-testing.csv"
download.file(url=test_url, destfile="testing.csv")

# Read in training and testing data

train <- read.csv("training.csv", na.strings=c("NA","#DIV/0!",""))
submit_test <- read.csv("testing.csv", na.strings=c("NA","#DIV/0!",""))

2.2 Data splitting and training set exploration

Now lets create a training and a testing data set, by splitting the train dataset in training and testing.

set.seed(1234)
in_train <- createDataPartition(y = train$classe,
                                p = 0.75,
                                list = FALSE)
training <- train[in_train,]
testing <- train[-in_train,]
dim(training); dim(testing); dim(submit_test)
## [1] 14718   160
## [1] 4904  160
## [1]  20 160

A rapid panoramic view of the classe variable that we will have to predict in submit_test. Its distribution suggests that a standardization should not be necessary.

histogram(~classe,  data = training,
          main = "Class of excercise")

summary(training$classe)
##    A    B    C    D    E 
## 4185 2848 2567 2412 2706
paste(round(prop.table(table(training$classe))*100, 2), "%", sep = "")
## [1] "28.43%" "19.35%" "17.44%" "16.39%" "18.39%"

3. Selecting the features

The train dataframe has many variables with a relevant amount of missing data. I will drop them, considering a really high cutpoint: 90%. This values depends by the particular distribution of missing data across the train dataframe.

na_count <-sapply(training, function(y) sum(length(which(is.na(y)))))
na_count <- data.frame(na_count)

The resulting dataset is thinner and more easily manageable. The submit_test dataset is reduced to the same dimensions, via column names comparison. An empty column classe is added.

summary(na_count)
##     na_count    
##  Min.   :    0  
##  1st Qu.:    0  
##  Median :14409  
##  Mean   : 9022  
##  3rd Qu.:14409  
##  Max.   :14718
# Eliminate variables with many NAs in dataset "train". Cutpoint 90%.

training_cut <- training[, ((na_count)/nrow(training)) < 0.9]

dim(training_cut)
## [1] 14718    60

Now, I am going to remove columns that are not predictors and with near zero variance:

training_cut2 <- training_cut[,8:length(training_cut)]

zero_var <- nearZeroVar(training_cut2, saveMetrics = TRUE)

4. Random Forest Model

Random forest model is appropriate for a classification problem, according to the lectures. Hence, I fit my model on training data first, and on the testing data used for cross validation.

set.seed(1234)
model_fit <- randomForest(classe~., data = training_cut2)
print(model_fit)
## 
## Call:
##  randomForest(formula = classe ~ ., data = training_cut2) 
##                Type of random forest: classification
##                      Number of trees: 500
## No. of variables tried at each split: 7
## 
##         OOB estimate of  error rate: 0.43%
## Confusion matrix:
##      A    B    C    D    E  class.error
## A 4182    2    0    0    1 0.0007168459
## B   13 2831    4    0    0 0.0059691011
## C    0   12 2554    1    0 0.0050642774
## D    0    0   20 2390    2 0.0091210614
## E    0    0    3    5 2698 0.0029563932

Cross validation testing for out-of-sample error.

prediction_1 <- predict(model_fit, testing, type = "class")
confusionMatrix(testing$classe, prediction_1)
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction    A    B    C    D    E
##          A 1395    0    0    0    0
##          B    2  945    2    0    0
##          C    0   10  844    1    0
##          D    0    0    8  796    0
##          E    0    0    0    0  901
## 
## Overall Statistics
##                                         
##                Accuracy : 0.9953        
##                  95% CI : (0.993, 0.997)
##     No Information Rate : 0.2849        
##     P-Value [Acc > NIR] : < 2.2e-16     
##                                         
##                   Kappa : 0.9941        
##  Mcnemar's Test P-Value : NA            
## 
## Statistics by Class:
## 
##                      Class: A Class: B Class: C Class: D Class: E
## Sensitivity            0.9986   0.9895   0.9883   0.9987   1.0000
## Specificity            1.0000   0.9990   0.9973   0.9981   1.0000
## Pos Pred Value         1.0000   0.9958   0.9871   0.9900   1.0000
## Neg Pred Value         0.9994   0.9975   0.9975   0.9998   1.0000
## Prevalence             0.2849   0.1947   0.1741   0.1625   0.1837
## Detection Rate         0.2845   0.1927   0.1721   0.1623   0.1837
## Detection Prevalence   0.2845   0.1935   0.1743   0.1639   0.1837
## Balanced Accuracy      0.9993   0.9943   0.9928   0.9984   1.0000

Cross validation for in-sample error

prediction_2 <- predict(model_fit, training_cut2, type = "class")
confusionMatrix(training_cut2$classe, prediction_2)
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction    A    B    C    D    E
##          A 4185    0    0    0    0
##          B    0 2848    0    0    0
##          C    0    0 2567    0    0
##          D    0    0    0 2412    0
##          E    0    0    0    0 2706
## 
## Overall Statistics
##                                      
##                Accuracy : 1          
##                  95% CI : (0.9997, 1)
##     No Information Rate : 0.2843     
##     P-Value [Acc > NIR] : < 2.2e-16  
##                                      
##                   Kappa : 1          
##  Mcnemar's Test P-Value : NA         
## 
## Statistics by Class:
## 
##                      Class: A Class: B Class: C Class: D Class: E
## Sensitivity            1.0000   1.0000   1.0000   1.0000   1.0000
## Specificity            1.0000   1.0000   1.0000   1.0000   1.0000
## Pos Pred Value         1.0000   1.0000   1.0000   1.0000   1.0000
## Neg Pred Value         1.0000   1.0000   1.0000   1.0000   1.0000
## Prevalence             0.2843   0.1935   0.1744   0.1639   0.1839
## Detection Rate         0.2843   0.1935   0.1744   0.1639   0.1839
## Detection Prevalence   0.2843   0.1935   0.1744   0.1639   0.1839
## Balanced Accuracy      1.0000   1.0000   1.0000   1.0000   1.0000

5. Prediction

Model applied to submit_test Upon submission all predictions were correct!

prediction_submit <- predict(model_fit, submit_test, type = "class")
print(prediction_submit)
##  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 
##  B  A  B  A  A  E  D  B  A  A  B  C  B  A  E  E  A  B  B  B 
## Levels: A B C D E
table(prediction_submit)
## prediction_submit
## A B C D E 
## 7 8 1 1 3