In Kuhn and Johnson do problems 6.2 and 6.3. There are only two but they consist of many parts. Please submit a link to your Rpubs and submit the .rmd file as well.

Exercise 6.2

Developing a model to predict permeability (see Sect. 1.4) could save significant resources for a pharmaceutical company, while at the same time more rapidly identifying molecules that have a sufficient permeability to become a drug:

a.

library(AppliedPredictiveModeling)
library(tidyverse)
library(caret)
library(RANN)
data(permeability)

b.

The fingerprint predictors indicate the presence or absence of substructures of a molecule and are often sparse meaning that relatively few of the molecules contain each substructure. Filter out the predictors that have low frequencies using the nearZeroVar function from the caret package. How many predictors are left for modeling?

The nearZeroVar function filtered out 719 predictors, leaving 388 predictors,

nearZero <- nearZeroVar(fingerprints)
lowFreqPred <- fingerprints[, -nearZero]
dim(lowFreqPred)
## [1] 165 388

Set dataframe before split

fingerprints_df <- fingerprints %>% 
    as_tibble() %>% 
    mutate(Target = permeability)

c.

Split the data into a training and a test set, pre-process the data, and tune a PLS model. How many latent variables are optimal and what is the corresponding resampled estimate of R2?

Rsquared is highest where ncomp = 8 (0.5030588).

set.seed(02180)

train_partition <- createDataPartition(fingerprints_df$Target, times = 1, p = 0.8, list = F)
train_set <- fingerprints_df[train_partition, ]
test_set <- fingerprints_df[-train_partition, ]

## Set up Partial Least Squares Model

pls_mod <- train(Target ~ ., data = train_set, metric = "Rsquared", method = "pls", center = T, trControl = trainControl("cv", number = 10), tuneLength = 25)
pls_mod
## Partial Least Squares 
## 
##  133 samples
## 1107 predictors
## 
## No pre-processing
## Resampling: Cross-Validated (10 fold) 
## Summary of sample sizes: 118, 119, 119, 121, 120, 120, ... 
## Resampling results across tuning parameters:
## 
##   ncomp  RMSE      Rsquared   MAE      
##    1     13.48756  0.3117614  10.506544
##    2     11.80448  0.4500985   8.639619
##    3     11.64713  0.4885465   8.702974
##    4     11.84474  0.4753360   8.973997
##    5     12.23946  0.4588232   9.292188
##    6     11.91100  0.4771585   9.026839
##    7     11.74528  0.4849933   9.046796
##    8     11.64133  0.5030588   8.867180
##    9     11.83744  0.4934233   8.916343
##   10     11.90220  0.4891102   8.964149
##   11     12.06976  0.4796015   9.135313
##   12     12.16789  0.4824971   9.133353
##   13     12.43952  0.4579997   9.362038
##   14     12.61017  0.4575988   9.470925
##   15     12.90751  0.4399179   9.703233
##   16     13.19751  0.4334757   9.797438
##   17     13.30235  0.4334355   9.932949
##   18     13.49177  0.4234612  10.053385
##   19     13.52269  0.4222782  10.028409
##   20     13.51128  0.4239817   9.961870
##   21     13.60855  0.4228055  10.021615
##   22     13.71529  0.4234043  10.115173
##   23     13.82370  0.4198877  10.157318
##   24     13.92916  0.4146172  10.247913
##   25     14.04928  0.4126800  10.358906
## 
## Rsquared was used to select the optimal model using the largest value.
## The final value used for the model was ncomp = 8.
summary(pls_mod)
## Data:    X dimension: 133 1107 
##  Y dimension: 133 1
## Fit method: oscorespls
## Number of components considered: 8
## TRAINING: % variance explained
##           1 comps  2 comps  3 comps  4 comps  5 comps  6 comps  7 comps
## X           25.12    37.00    43.80    47.45    52.77    58.86    61.30
## .outcome    30.17    53.74    62.11    71.34    77.34    80.68    83.59
##           8 comps
## X           62.80
## .outcome    87.29

d.

Predict the response for the test set. What is the test set estimate of R2?

An Rsquared of 0.3182895 for the test dataset is comparatively worse than the training data set

# Set predictor test
pred_test <- predict(pls_mod, test_set)

pred_data <- data.frame(obs = test_set$Target, pred = pred_test)

# Set colnames for default summary
colnames(pred_data) <- c("obs", "pred")

defaultSummary(pred_data)
##       RMSE   Rsquared        MAE 
## 11.3829626  0.3182895  7.9950484

e.

Try building other models discussed in this chapter. Do any have better predictive performance?

I will use PCR to try and find a better fitting model The PCR model yielded a lower Rsquared score compared to my original model which was sitting around .50 R2

set.seed(02180)

PCR_mod <-  train(Target ~ ., data = train_set, method = "pcr",
                  center = T, trControl = trainControl("cv", number = 10),
                  tuneLength = 25)

pred_pcr <- predict(PCR_mod, test_set)

pred_pcr_data <- data.frame(obs = test_set$Target, pred = pred_pcr)

colnames(pred_pcr_data) <- c("obs","pred")

defaultSummary(pred_pcr_data)
##       RMSE   Rsquared        MAE 
## 10.7883822  0.3865399  6.4985218

f.

Would you recommend any of your models to replace the permeability laboratory experiment?

The PCR model in this instance performed much worse than my original PLS model, thus I would likely continue using the PLS model.

Exercise 6.3

a.

Start R and use these commands to load the data:

data("ChemicalManufacturingProcess")
dim(ChemicalManufacturingProcess)
## [1] 176  58

b.

A small percentage of cells in the predictor set contain missing values. Use an imputation function to fill in these missing values (e.g., see Sect. 3.8).

set.seed(02180)
missing_val <- preProcess(ChemicalManufacturingProcess, method = c('knnImpute'))
missing_df <- predict(missing_val, ChemicalManufacturingProcess)
chem_pre_process <- missing_df[, -nearZeroVar(missing_df)]

c.

Split the data into a training and a test set, pre-process the data, and tune a model of your choice from this chapter. What is the optimal value of the performance metric?

Ill start with PCA: Ncomp = 18 is where the Rsquared is maximized and the RMSE was minimized. With an R2 of 0.5501713 and an RMSE of 0.6662796

dp_chem <- createDataPartition(chem_pre_process$Yield, p = 0.8, list = F)
train_chem <- chem_pre_process[dp_chem,]
test_chem <- chem_pre_process[-dp_chem,]
pcr_chem <- train(Yield ~., data = train_chem, method = "pcr", center = T, trControl = trainControl("cv", number = 10), tuneLength = 25)
pcr_chem
## Principal Component Analysis 
## 
## 144 samples
##  56 predictor
## 
## No pre-processing
## Resampling: Cross-Validated (10 fold) 
## Summary of sample sizes: 131, 128, 129, 130, 128, 129, ... 
## Resampling results across tuning parameters:
## 
##   ncomp  RMSE       Rsquared   MAE      
##    1     0.8323741  0.3199209  0.6818622
##    2     0.8328189  0.3184820  0.6823891
##    3     0.8119696  0.3688657  0.6523541
##    4     0.8286502  0.3501283  0.6627367
##    5     0.8001197  0.3732916  0.6541014
##    6     0.7902794  0.3847506  0.6517525
##    7     0.7790191  0.4109751  0.6354220
##    8     0.7555472  0.4320210  0.6147314
##    9     0.7142023  0.5002865  0.5755154
##   10     0.6675721  0.5438919  0.5398642
##   11     0.6688872  0.5428337  0.5394437
##   12     0.6737817  0.5423294  0.5459377
##   13     0.6788423  0.5365625  0.5546770
##   14     0.6749455  0.5334829  0.5552679
##   15     0.6749398  0.5374176  0.5530005
##   16     0.6888797  0.5201496  0.5634165
##   17     0.6808903  0.5407401  0.5574088
##   18     0.6948702  0.5143322  0.5714891
##   19     0.6962055  0.5174285  0.5696421
##   20     0.7073315  0.5121826  0.5766562
##   21     0.7004638  0.5279897  0.5673212
##   22     0.6974623  0.5330592  0.5675422
##   23     0.7016080  0.5285056  0.5703760
##   24     0.7032819  0.5318976  0.5694746
##   25     0.6988398  0.5366181  0.5672114
## 
## RMSE was used to select the optimal model using the smallest value.
## The final value used for the model was ncomp = 10.

d. 

Predict the response for the test set.What is the value of the performance metric and how does this compare with the resampled performance metric on the training set?

The RMSE lowered to 0.6359710 and Rsquared jumped to 0.6128665

pred_chem <- predict(pcr_chem, test_chem)
chem_pcr_test <- data.frame(obs = test_chem$Yield, pred = pred_chem)
colnames(chem_pcr_test) <- c("obs", "pred")
defaultSummary(chem_pcr_test)
##      RMSE  Rsquared       MAE 
## 0.7183884 0.4983845 0.6123158

e.

Which predictors are most important in the model you have trained? Do either the biological or process predictors dominate the list?

I can use varImp from caret to check out which predictors are the most important in the model. ManufacturingProcess and Biological Material are a mix of the top 10 most important variable however BiologicalMaterial06 & ManufacturingProcess13 seem to be the highest. A corr plot

varImp(pcr_chem)
## loess r-squared variable importance
## 
##   only 20 most important variables shown (out of 56)
## 
##                        Overall
## ManufacturingProcess32  100.00
## BiologicalMaterial06     95.82
## ManufacturingProcess13   95.51
## ManufacturingProcess17   84.82
## ManufacturingProcess31   77.70
## BiologicalMaterial03     73.77
## BiologicalMaterial12     73.32
## ManufacturingProcess36   68.91
## BiologicalMaterial02     67.38
## ManufacturingProcess09   65.95
## BiologicalMaterial04     59.46
## ManufacturingProcess06   58.67
## ManufacturingProcess02   51.27
## ManufacturingProcess33   51.26
## BiologicalMaterial11     49.56
## ManufacturingProcess11   47.71
## ManufacturingProcess29   47.00
## BiologicalMaterial08     43.37
## BiologicalMaterial01     42.41
## BiologicalMaterial09     34.82

f.

Explore the relationships between each of the top predictors and the response. How could this information be helpful in improving yield in future runs of the manufacturing process?

This corr plot shows that the ManufacturingProcess predictors are each other, while the BiologicalMaterials are also related to each other. The Yield is highly related to a few ManufacturingProcess and BiologicalMaterial predictors (both negatively and positively) while not hardly related to others. In this case the manufacturer might want to maximize the amount of MP32, BM06, BM03 and BM12 used while trying to minimize the amount of MP13, MP17 and MP31 used (negative or no correlation to yield)

correlation <- cor(select(chem_pre_process, 'ManufacturingProcess32','BiologicalMaterial06','ManufacturingProcess13', 'ManufacturingProcess17', 'ManufacturingProcess31','BiologicalMaterial03','BiologicalMaterial12', 'Yield'))
corrplot::corrplot(correlation, method='square', type="upper")