In Kuhn and Johnson do problems 7.2 and 7.5

Ex 7.2

7.2. Friedman (1991) introduced several benchmark data sets create by simulation. One of these simulations used the following nonlinear equation to create data:

y = 10 sin(πx1x2) + 20(x3 − 0.5)2 + 10x4 + 5x5 + N(0, σ2)

where the x values are random variables uniformly distributed between [0, 1] (there are also 5 other non-informative variables also created in the simulation). The package mlbench contains a function called mlbench.friedman1 that simulates these data:

library(AppliedPredictiveModeling)
library(mlbench)
library(caret)

set.seed(200)
trainingData <- mlbench.friedman1(200, sd = 1)

## We convert the 'x' data from a matrix to a data frame
## One reason is that this will give the columns names.
trainingData$x <- data.frame(trainingData$x)

## Look at the data using
featurePlot(trainingData$x, trainingData$y)

#This creates a list with a vector 'y' and a matrix of predictors 'x'. Also simulate a large test set to estimate the true error rate with good precision:

testData <- mlbench.friedman1(5000, sd = 1)
testData$x <- data.frame(testData$x)

#Tune several models on these data. For example:
library(caret)
knnModel <- train(x = trainingData$x, y = trainingData$y, method = "knn", preProc = c("center", "scale"), tuneLength = 10)

knnModel
## k-Nearest Neighbors 
## 
## 200 samples
##  10 predictor
## 
## Pre-processing: centered (10), scaled (10) 
## Resampling: Bootstrapped (25 reps) 
## Summary of sample sizes: 200, 200, 200, 200, 200, 200, ... 
## Resampling results across tuning parameters:
## 
##   k   RMSE      Rsquared   MAE     
##    5  3.466085  0.5121775  2.816838
##    7  3.349428  0.5452823  2.727410
##    9  3.264276  0.5785990  2.660026
##   11  3.214216  0.6024244  2.603767
##   13  3.196510  0.6176570  2.591935
##   15  3.184173  0.6305506  2.577482
##   17  3.183130  0.6425367  2.567787
##   19  3.198752  0.6483184  2.592683
##   21  3.188993  0.6611428  2.588787
##   23  3.200458  0.6638353  2.604529
## 
## RMSE was used to select the optimal model using the smallest value.
## The final value used for the model was k = 17.
knnPred <- predict(knnModel, newdata = testData$x)

## The function 'postResample' can be used to get the test set performance values
postResample(pred = knnPred, obs = testData$y)
##      RMSE  Rsquared       MAE 
## 3.2040595 0.6819919 2.5683461

Which models appear to give the best performance? Does MARS select the informative predictors (those named X1–X5)?

#Try MARS model
library(earth)
marsGrid <- expand.grid(.degree = 1:2, .nprune = 2:12)
# Fix the seed so that the results can be reproduced
set.seed(1111)
marsTuned <- train(x=trainingData$x, y=trainingData$y, method = "earth", tuneGrid = marsGrid, trControl = trainControl(method = "cv"))

marsPred<-predict(marsTuned, testData$x)

#Try Support Vector Machine model
library(kernlab)

svmRTuned <- train(x=trainingData$x, y=trainingData$y, method = "svmRadial", preProc =c("center", "scale"), tuneLength = 14, trControl = trainControl(method = "cv"))

svmRPred<-predict(svmRTuned, testData$x)

#Try Neural Network model
library(nnet)

nnetGrid <- expand.grid(.decay = c(0, 0.01, .1), .size = c(1:10), .bag = FALSE)

set.seed(999)
crtl<-trainControl(method = "cv")

nnetTuned <- train(x=trainingData$x, y=trainingData$y, method = "avNNet", tuneGrid = nnetGrid, trControl = crtl, preProc = c("center", "scale"), linout = TRUE, trace = FALSE, MaxNWts = 10 * (ncol(trainingData$x) + 1) + 10 + 1, maxit = 500)

nnetPred<-predict(nnetTuned, testData$x)

#Create a summary for all models:
rbind(
  "knn" = postResample(pred = knnPred, obs = testData$y),
  "mars" = postResample(pred = marsPred, obs = testData$y),
  "svm" = postResample(pred = svmRPred, obs = testData$y),
  "nnet" = postResample(pred = nnetPred, obs = testData$y)
)
##          RMSE  Rsquared      MAE
## knn  3.204059 0.6819919 2.568346
## mars 1.280306 0.9335241 1.016867
## svm  2.077675 0.8250305 1.578314
## nnet 2.290661 0.7943951 1.723649

From the value of RMSE and R2, the best model is MARS.

varImp(marsTuned)
## earth variable importance
## 
##    Overall
## X1  100.00
## X4   75.33
## X2   48.88
## X5   15.63
## X3    0.00

MARS selected the informative predictors named X1–X5 as shown above.

Ex 7.5

Exercise 6.3 describes data for a chemical manufacturing process. Use the same data imputation, data splitting, and pre-processing steps as before and train several nonlinear regression models.

#Prepare the dataframe
data(ChemicalManufacturingProcess)
ChemicalManufacturingProcess<-as.data.frame(ChemicalManufacturingProcess)

#The Impute package is not available, I can only use preProcess with KNN for imputation
preProcValues<-preProcess(ChemicalManufacturingProcess[,-1], method="knnImpute", k=5, knnSummary=mean)

impute_CMP<-predict(preProcValues, ChemicalManufacturingProcess[,-1], na.action=na.pass)

#The impute_CMP looks strange and I found out it's normalized, we need to de-normalize and get the original data back
procNames <- data.frame(col = names(preProcValues$mean), mean = preProcValues$mean, sd = preProcValues$std)
for(i in procNames$col){
 impute_CMP[i] <- impute_CMP[i]*preProcValues$std[i]+preProcValues$mean[i] 
}

set.seed(2744)
fold2 <- ChemicalManufacturingProcess$Yield|> createDataPartition(p = 0.9, list = FALSE, times = 1)

#Create training and testing set
bio_manu_train<-impute_CMP[fold2,]
bio_manu_test<-impute_CMP[-fold2,]
yield_train<-ChemicalManufacturingProcess[fold2,1]
yield_test<-ChemicalManufacturingProcess[-fold2,1]
  1. Which nonlinear regression model gives the optimal resampling and test set performance? Let try KNN, MARS SVM and Neural Network model
#Try KNN
knnModel2 <- train(x = bio_manu_train, y = yield_train, method = "knn", preProc = c("center", "scale"), tuneLength = 10)
knnPred2 <- predict(knnModel2, newdata = bio_manu_test)

#Try MARS model
marsGrid <- expand.grid(.degree = 1:2, .nprune = 2:12)
# Fix the seed so that the results can be reproduced
set.seed(1111)
marsTuned2 <- train(x=bio_manu_train, y=yield_train, method = "earth", tuneGrid = marsGrid, trControl = trainControl(method = "cv"))
marsPred2<-predict(marsTuned2, bio_manu_test)

#Try Support Vector Machine model
svmRTuned2 <- train(x=bio_manu_train, y=yield_train, method = "svmRadial", preProc =c("center", "scale"), tuneLength = 14, trControl = trainControl(method = "cv"))
svmRPred2<-predict(svmRTuned2, bio_manu_test)

#Try Neural Network model
nnetGrid <- expand.grid(.decay = c(0, 0.01, .1), .size = c(1:10), .bag = FALSE)
set.seed(999)
crtl<-trainControl(method = "cv")
nnetTuned2 <- train(x=bio_manu_train, y=yield_train, method = "avNNet", tuneGrid = nnetGrid, trControl = crtl, preProc = c("center", "scale"), linout = TRUE, trace = FALSE, MaxNWts = 10 * (ncol(trainingData$x) + 1) + 10 + 1, maxit = 500)
nnetPred2<-predict(nnetTuned2, bio_manu_test)

#Create a summary for all models:
rbind(
  "knn" = postResample(pred = knnPred2, obs = yield_test),
  "mars" = postResample(pred = marsPred2, obs = yield_test),
  "svm" = postResample(pred = svmRPred2, obs = yield_test),
  "nnet" = postResample(pred = nnetPred2, obs = yield_test)
)
##          RMSE  Rsquared       MAE
## knn  1.366693 0.4979899 1.0942361
## mars 1.303821 0.4802854 1.0066770
## svm  0.932722 0.7305057 0.7397594
## nnet 1.259355 0.5491511 0.9441607

According to the RMSE and R2 value, the SVM model is the best among these four modesl.

  1. Which predictors are most important in the optimal nonlinear regression model? Do either the biological or process variables dominate the list? How do the top ten important predictors compare to the top ten predictors from the optimal linear model?
SVMRImp<-varImp(svmRTuned2)

#Redo the previous linear model (PLS)
bio_manu_pls <- train(x = bio_manu_train, y = yield_train, method = "pls", tuneLength = 20, trControl = trainControl(method = "cv"), preProcess = c("center", "scale"))

PLSImp<-varImp(bio_manu_pls)

#Plot both Important data in parallel
p1<-plot(PLSImp, top=10, main='Linear Model: PLS')
p2<-plot(SVMRImp, top=10, main='Non-linear Model: SVM')
gridExtra::grid.arrange(p1, p2, ncol = 2)

Only the top 10 variables are shown above. There are both biological and process variables on the list. Comparing to the top ten predictors from the optimal linear model (PLS) which has a dominated list of manufacturing process with 2 biological material in top 10, this list from non-linear model makes more sense as it includes 4 biological material and 6 manufacturing process.

  1. Explore the relationships between the top predictors and the response for the predictors that are unique to the optimal nonlinear regression model. Do these plots reveal intuition about the biological or process predictors and their relationship with yield?
#Pick the top 3 predictors for the analysis
top3 <- ChemicalManufacturingProcess|>
  dplyr::select(Yield, ManufacturingProcess32, BiologicalMaterial06, ManufacturingProcess13)

ggplot(top3, aes(x=ManufacturingProcess32,y=Yield)) + geom_point()

ggplot(top3, aes(x=BiologicalMaterial06,y=Yield)) + geom_point()

ggplot(top3, aes(x=ManufacturingProcess13,y=Yield)) + geom_point()

cor_res_t3 <- cor(top3, use = "na.or.complete")
corrplot::corrplot(cor_res_t3,
  type = "lower",
  order = "original",
  tl.col = "Blue",
  tl.srt = 45,
  tl.cex = 0.5
  )

Both manufacturing process predictors are quite linear, but the biological material shows some non-linear relationship with yield. Moreover, there is correlation between the biological materials and manufacturing process.