Computing part of Chapter 6

Today, I continued to finish Computing part of Chapter 6.

Codes (ordianry linear regression):
corthresh=0.9
toohigh=findCorrelation(cor(solTrainXtrans), corthresh)
corrpred=names(solTrainXtrans)[toohigh]
trainxfiltered=solTrainXtrans[, -toohigh]
testxfiltered=solTestXtrans[, -toohigh]
set.seed(100)
lmfiltered=train(x=trainxfiltered, y=solTrainY, method = "lm", trControl = ctr1)
lmfiltered
#solTrainY: 951*1; solTrainXtrans: 951*228; solTestXtrans: 316*228;
#trainxfiltered: 951*190; testxfiltered: 316*190
#the codes in P132 of Applied Predictive Modeling may be wrong
set.seed(100)
rlmpca=train(solTrainXtrans, solTrainY, method = "rlm", preProcess = "pca", trControl = ctr1)
#wait a minute
rlmpca

Codes (partial least squares):
library(pls)
plsfit=plsr(solubility~., data=trainingdata)
predict(plsfit, solTestXtrans[1:5,], ncomp = 1:2)
set.seed(100)
plstune=train(solTrainXtrans, solTrainY, method = "pls", tuneLength = 20, trControl = ctr1, preProcess = c("center", "scale"))
plstune

Codes (penalized regression models):
library(elasticnet)
#create ridge-regression models
#as lambda decreases to 0, we get the least squares solutions. lambda argument specifies the ridge-regression penalty
ridgemodel=enet(x=as.matrix(solTrainXtrans), y=solTrainY, lambda = 0.001)
ridgemodel
ridgepred=predict(ridgemodel, newx = as.matrix(solTestXtrans), s=1, model="fraction", type = "fit")
ridgepred
head(ridgepred$fit)
#define the candidate set of values
ridgegrid=data.frame(.lambda=seq(0, 0.1, length=15))
set.seed(100)
ridgeregfit=train(solTrainXtrans, solTrainY, method="ridge", tuneGrid=ridgegrid, trControl=ctr1, preProcess=c("center", "scale"))
ridgeregfit
enetmodel=enet(x=as.matrix(solTrainXtrans), y=solTrainY, lambda = 0.01, normalize = TRUE)
#normalize is for centering and scaling predictors prior to modeling
#lambda controls the ridge-regression penalty and, setting this value to 0, fits the lasso model
enetpred=predict(enetmodel, newx = as.matrix(solTestXtrans), s=0.1, mode="fraction", type="fit")
names(enetpred)
head(enetpred$fit)
enetcoef=predict(enetmodel, newx = as.matrix(solTestXtrans), s=0.1, mode="fraction", type="coefficients")
tail(enetcoef$coefficients)
#other packages: biglars (for large data sets), FLLat (for the fused lasso), grplasso (the group lasso), penalized, relaxo (the relaxed lasso)...
enetgrid=expand.grid(.lambda=c(0, 0.01, 0.1), .fraction=seq(0.05, 1, length=20))
set.seed(100)
enetTune=train(solTrainXtrans, solTrainY, method = "enet", tuneGrid = enetgrid, trControl = ctr1, preProcess = c("center", "scale"))
enetTune
plot(enetTune)

Tomorrow, I will do exercises of Chapter 6.

Chapter 6 Computing

Today, I did computing of Chapter 6.

library(AppliedPredictiveModeling)
data(solubility)
ls(pattern = "^solT")
set.seed(2)
sample(names(solTrainX), 8)
# 1. ordianry linear regression
trainingdata=solTrainXtrans
#add the solubility outcome
trainingdata$solubility=solTrainY
lmfitallpredictors=lm(solubility~., data=trainingdata)
summary(lmfitallpredictors)
lmpred1=predict(lmfitallpredictors, solTestXtrans)
head(lmpred1)
nrow(solTestXtrans)
lmvalues1=data.frame(obs=solTestY, pred=lmpred1)
library(caret)
defaultSummary(lmvalues1)
library(MASS)
rlmfitallpredictors=rlm(solubility~., data=trainingdata)
summary(rlmfitallpredictors)
rlmpred1=predict(rlmfitallpredictors, solTestXtrans)
head(rlmpred1)
rlmvalues1=data.frame(obs=solTestY, pred=rlmpred1)
defaultSummary(rlmvalues1)
#generate a resampling estimate of performance
ctr1=trainControl(method = "cv", number = 10)
set.seed(100)
lmfit1=train(x=solTrainXtrans, y=solTrainY, method = "lm", trControl = ctr1)
xyplot(solTrainY~predict(lmfit1),
       type=c("p","g"), xlab = "Predicted", ylab = "Observed")
xyplot(resid(lmfit1)~predict(lmfit1),
       type=c("p","g"), xlab = "Predicted", ylab = "Residuals")
#resid: generate the model residuals; predict: return the predicted values

I will try to finish computing at the weekend.

Chapter 6 (6.3 6.4)

Today, I finished reading 6.3 and 6.4 of Chapter 6.

Summary

Both PCR and PLS have similar predictive ability, but PLS does so with far fewer components.

The NIPALS algorithm works fairly efficiently for data sets of small-to-moderate size (< 2500 samples and < 30 predictors). When the number of samples and predictors climbs, the algorithm becomes inefficient.

Kernel approach: improve the speed of the algorithm

SIMPLS: deflate the covariance matrix between the predictors and the response

Covariance: .

GIFI approach: split each predictor into two or more bins for those predictors that are thought to have a nonlinear relationship with the response. Cut points for the bins are selected by the user and are based on either prior knowledge or characteristics of the data.

Penalized models:

A generalization of the lasso model is the elastic net:

This model will more effectively deal with groups of high correlated predictors.

Tomorrow, I will begin computing of Chapter 6.

Chapter 6 (6.1 6.2 6.3)

Today, I read some parts of Chapter 6.

Summary

QSAR: quantitative structure-activity relationship modeling

NIPALS: nonlinear iterative partial least squares algorithm

The objective of ordinary least squares linear regression is to find the plane that minimizes the sum-of-squared errors (SSE) between the observed and predicted response:

On many occasions, relationships among predictors can be complex and involve many predictors. In these cases, manual removal of specific predictors may not be possible and models that can tolerate collinearity may be more useful.

Box-Cox transformation can be applied to the continuous predictors in order to remove skewness.

If the correlation among predictors is high, then the ordinary least squares solution for multiple linear regression will have high variability and will become unstable.

If the number of predictors are greater than the number of observations, ordinary least squares in its usual form will be unable to find a unique set of regression coefficients that minimize the SSE.

Pre-processing predictors via PCA (dimension reduction) prior to performing regression is known as principal component regression (PCR). It has been widely applied in the context of problems with inherently highly correlated predictors or problems with more predictors than observations.

The author recommends using PLS when there are correlated predictors and a linear regression-type solution is desired.

PCA

PLS

While PCA linear combinations are chosen to maximally summarize predictor space variability, the PLS linear combinations of predictors are chosen to maximally summarize covariance with the response.

Prior to performing PLS, the predictors should be centered and scaled, especially if the predictors are on scales of differing magnitude.

PLS has one tuning parameter: the number of components to retain. Resampling techniques can be used to determine the optimal number of components.

Tomorrow, I will continue to read Chapter 6.

Chapter 5 Computing

Today, I finished Computing on Chapter 5.

observed=c(0.22, 0.83, -0.12, 0.89, -0.23, -1.30, -0.15, -1.4, 0.62, 0.99, -0.18, 0.32, 0.34,
           -0.3, 0.04, -0.87, 0.55, -1.3, -1.15, 0.2)
predicted=c(0.24, 0.78, -0.66, 0.53, 0.7, -0.75, -0.41, -0.43, 0.49, 0.79, -1.19, 0.06, 0.75,
           -0.07, 0.43, -0.42, -0.25, -0.64, -1.26, -0.07)
#in practice, the vector of predictions would be produced by the model function
residualvalues=observed-predicted
residualvalues
summary(residualvalues)
axisrange=extendrange(c(observed, predicted))
plot(observed, predicted, ylim=axisrange, xlim=axisrange)
abline(0, 1, col="darkgrey", lty=2)
plot(predicted, residualvalues, ylab = "residual")
abline(h=0, col="darkgrey", lty=2)
library(caret)
R2(predicted, observed)
RMSE(predicted, observed)
cor(predicted, observed)
cor(predicted, observed, method = "spearman")

Next week, I will read on Chapter 6.

Chapter 5 Reading

Today, I read Chapter 5 (5.1 5.2) of the book.

Summary

Quantitative Measures of Performance
RMSE: root mean squared error (take the square root of MSE so as to have the same units as the original data)
Coefficient of determination (R^2): the proportion of the information in the data that is explained by the model. It is a measure of correlation, not accuracy. It is dependent on the variation in the outcome. Same RMSE, larger variance, larger R^2.
Spearman’s rank correlation: the rank of the observed and predicted outcomes are obtained and the correlation coefficient between these ranks is calculated.

The Variance-Bias Trade-off
E[MSE]=σ^2+〖(Model Bias)〗^2+Model Variance
σ^2: residual variance, also called ‘irreducible noise’. ‘Model Bias’ reflects how close the functional form of the model can get to the true relationship between the predictors and the outcome.
It is generally true that more complex models can have very high variance, which leads to over-fitting.
The Variance-Bias Trade-off: increase the bias in the model to greatly reduce the model variance as a way to mitigate the problem of collinearity.

Tomorrow, I will do computing of Chapter 5.

Exercise 4.3 4.4

Today, I finished exercises in Chapter 4.

4.3
library(AppliedPredictiveModeling)
data(ChemicalManufacturingProcess)
#objective: find the number of PLS components that yields the optimal R2 value
#7 provides the most paprsimonious model using "one-standard error" method
exercise4.3
#2 is the best if tolerance should be less than 10%
colnames(exercise4.3)=c("Components", "Mean", "Tolerance, %", "Std.Error")
exercise4.3
names(exercise4.3)[names(exercise4.3)=="Tolerance, %"]="Tolerance"
exercise4.3
names(exercise4.3)[names(exercise4.3)=="Tolerance"]="Tolerance, %"
exercise4.3
#select Random Forests model that optimizes R2
#select SVM model with combined consideration of R2, prediction time and model complexity

4.4
library(caret)
data(oil)
oilType
table(oilType)
sam1=sample(oilType, 60, replace = F)
sam1
table(sam1)
set.seed(1)
sam2=createDataPartition(oilType, p=0.6, list=FALSE)
sam2
sam3=oilType[sam2]
sam3
table(sam3)
set.seed(1028)
sam2=createDataPartition(oilType, p=0.6, list=FALSE)
sam3=oilType[sam2]
table(sam3)
#obtain a confidence interval for the overall accuracy
binom.test(16, 20)
binom.test(15, 20)
#try different sample sizes and accuracy rates to understand the trade-off
#between the uncertainty in the results, the model performance and the test set size
#p-value: how extreme the observation is
#confidence interval: the probability of success in the interval

Tomorrow, I will start reading Chapter 5.