APM Chapter 7 HW: Nonlinear Regression

Author

Alex Ptacek

library(tidyverse)
library(mlbench)
library(caret)
library(janitor)
library(earth)
library(kernlab)
library(VIM)

set.seed(624)

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

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:

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)
## or other methods.

## 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:

knnModel <- train(x = trainingData$x,
y = trainingData$y,
method = "knn",
preProc = c("center", "scale"),
tuneLength = 10)

knnModel

knnPred <- predict(knnModel, newdata = testData$x)
## The function 'postResample' can be used to get the test set
## perforamnce values
postResample(pred = knnPred, obs = testData$y)

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

Answer: The MARS model has the best performance because it has the lowest RMSE and highest R-squared. The tuned MARS model has incorporates predictors X1-X5 and excludes the five non-informative variables, and it explains 93.8% of the test set variance. Overall, it is a great model for this data.

mars_grid <- expand.grid(.degree = 1:2, .nprune = 2:10)
ctrl <- trainControl(method = "cv")

mars_tuned <- train(trainingData$x, trainingData$y,
                    method = "earth",
                    tuneGrid = mars_grid,
                    trControl = ctrl)

mars_pred <- predict(mars_tuned, newdata = testData$x)

postResample(pred = mars_pred, obs = testData$y)

     RMSE  Rsquared       MAE 
1.3738355 0.9277961 1.0903503

summary(mars_tuned)

Call: earth(x=data.frame[200,10], y=c(10.41,16.21,1...), keepxy=TRUE, degree=2,
            nprune=10)

                                coefficients
(Intercept)                         2.536313
h(X1-0.274555)                     11.969090
h(X2-0.334452)                     27.428282
h(X3-0.295314)                     11.984873
h(0.668056-X3)                     13.365415
h(X4-0.183828)                     10.828397
h(0.784898-X5)                     -4.669924
h(X1-0.139921) * h(0.334452-X2)   -31.677802
h(X1-0.274555) * h(X2-0.4577)     -67.791291
h(0.607225-X1) * h(X2-0.334452)   -50.249746

Selected 10 of 20 terms, and 5 of 10 predictors (nprune=10)
Termination condition: Reached nk 21
Importance: X4, X1, X2, X5, X3, X6-unused, X7-unused, X8-unused, X9-unused, ...
Number of terms at each degree of interaction: 1 6 3
GCV 1.782412    RSS 277.6306    GRSq 0.9235662    RSq 0.9398731

svm_tuned <- train(trainingData$x, trainingData$y,
                   method = "svmRadial",
                   preProcess = c("center", "scale"),
                   tuneLength = 14,
                   trControl = ctrl)

svm_pred <- predict(svm_tuned, newdata = testData$x)

postResample(pred = svm_pred, obs = testData$y)

     RMSE  Rsquared       MAE 
2.1112981 0.8302345 1.6246801

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.

**ALL DATA PROCESSING STEPS FROM EXERCISE 6.3**

library(AppliedPredictiveModeling)
data(ChemicalManufacturingProcess)

chem_manu_process <- ChemicalManufacturingProcess |> 
  as_tibble() |> 
  clean_names()

na_counts <- sapply(chem_manu_process, \(x) sum(is.na(x)))

na_counts <- tibble(vars = names(na_counts),
       na_counts = as.integer(na_counts)) |> 
  filter(na_counts > 0)

missing_chem <- chem_manu_process |> 
  select(all_of(na_counts$vars))

imp_values <- kNN(missing_chem, imp_var = FALSE) |> as_tibble()

chem_data_wip <- chem_manu_process |> 
  select(-all_of(na_counts$vars)) |> 
  bind_cols(imp_values)

chem_data_wip <- chem_data_wip |> 
  mutate(index = row_number())

chem_train <- slice_sample(chem_data_wip, prop = 0.8)

chem_test <- chem_data_wip |> 
  filter(!index %in% chem_train$index) |> 
  select(-index)

chem_train <- chem_train |> 
  select(-index)

chem_train_x <- chem_train |> 
  select(-yield)

chem_train_y <- chem_train |> 
  pull(yield) |> 
  as.vector()

chem_test_x <- chem_test |> 
  select(-yield)

chem_test_y <- chem_test |> 
  select(yield)

Which nonlinear regression model gives the optimal resampling and test set performance?

Answer: The MARS model has the best resampling and test set performance.

mars_grid <- expand.grid(.degree = 1:2, .nprune = 2:30)

mars_tuned <- train(chem_train_x, chem_train_y,
                    method = "earth",
                    tuneGrid = mars_grid,
                    trControl = ctrl)

mars_pred <- predict(mars_tuned, newdata = chem_test_x)

postResample(pred = mars_pred, obs = chem_test_y)

     RMSE  Rsquared       MAE 
0.7381151 0.7514460 0.6406255

svm_tuned <- train(chem_train_x, chem_train_y,
                   method = "svmRadial",
                   preProcess = c("center", "scale"),
                   tuneLength = 14,
                   trControl = ctrl)

svm_pred <- predict(svm_tuned, newdata = chem_test_x)

postResample(pred = svm_pred, obs = as_vector(chem_test_y))

     RMSE  Rsquared       MAE 
0.8834718 0.6409164 0.7154885

knn_tuned <- train(chem_train_x, chem_train_y,
                  method = "knn",
                  preProc = c("center", "scale"),
                  tuneLength = 20,
                  trControl = ctrl)

knn_pred <- predict(knn_tuned, newdata = chem_test_x)

postResample(pred = knn_pred, obs = as_vector(chem_test_y))

     RMSE  Rsquared       MAE 
1.0730033 0.4712348 0.9021605

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?

Answer: The MARS model only has 4 predictors, which are all hinge functions of process variables. Of these, `manufacturing_process32` and `manufacturing_process09` are the most important. This result is similar to my linear model, because the linear model was also dominated by process variables. Furthermore, the most important variables in the MARS model were also in the top ten in my linear model.

varImp(mars_tuned)

earth variable importance

                        Overall
manufacturing_process32  100.00
manufacturing_process09   58.21
manufacturing_process13   16.23
manufacturing_process39    0.00

summary(mars_tuned)

Call: earth(x=tbl_df[140,57], y=c(40.59,39.84,4...), keepxy=TRUE, degree=1,
            nprune=6)

                                 coefficients
(Intercept)                         38.255344
h(43.92-manufacturing_process09)    -0.716177
h(manufacturing_process09-43.92)     0.406495
h(33.1-manufacturing_process13)      2.794455
h(manufacturing_process32-154)       0.239567
h(6.9-manufacturing_process39)      -0.157776

Selected 6 of 21 terms, and 4 of 57 predictors (nprune=6)
Termination condition: RSq changed by less than 0.001 at 21 terms
Importance: manufacturing_process32, manufacturing_process09, ...
Number of terms at each degree of interaction: 1 5 (additive model)
GCV 1.483282    RSS 176.3092    GRSq 0.6038459    RSq 0.6587961

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?

Answer: These plots provide further evidence of a strong relationship between these predictors and `yield`. Both predictors appear to have a positive correlation with yield. The plot of `manufacturing_process32` also appears to have a sinusoidal pattern.

ggplot(chem_manu_process, aes(x = manufacturing_process32, y = yield)) +
  geom_point()

ggplot(chem_manu_process, aes(x = manufacturing_process09, y = yield)) +
  geom_point()

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

Answer: The MARS model has the best performance because it has the lowest RMSE and highest R-squared. The tuned MARS model has incorporates predictors X1-X5 and excludes the five non-informative variables, and it explains 93.8% of the test set variance. Overall, it is a great model for this data.

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.

Answer: The MARS model has the best resampling and test set performance.

Answer: These plots provide further evidence of a strong relationship between these predictors and yield. Both predictors appear to have a positive correlation with yield. The plot of manufacturing_process32 also appears to have a sinusoidal pattern.

Answer: These plots provide further evidence of a strong relationship between these predictors and `yield`. Both predictors appear to have a positive correlation with yield. The plot of `manufacturing_process32` also appears to have a sinusoidal pattern.