Applied_Predictive_Model_11
来源：互联网 发布：spark时间序列算法 编辑：程序博客网 时间：2024/06/11 05:24
一直在用Python/R/Shell做后台开发，但缺少交流。第一次尝试把自己学习中东西贴到博客中，希望得到高手们的指点让自己也能得到一些阳光普照。
library(AppliedPredictiveModeling)data(concrete)library(caret)library(plyr)head(concrete,2)featurePlot(concrete[, -9], concrete$CompressiveStrength, between = list(x = 1, y = 1), type = c("g", "p", "smooth"))
averaged <- ddply(mixtures, .(Cement, BlastFurnaceSlag, FlyAsh, Water, Superplasticizer, CoarseAggregate, FineAggregate, Age), function(x) c(CompressiveStrength = mean(x$CompressiveStrength)))head(mixtures,2)set.seed(975)inTrain <- createDataPartition(averaged$CompressiveStrength, p = 3/4)[[1]]training <- averaged[ inTrain,]testing <- averaged[-inTrain,]### Repeated Training/Test Splitsctrl <- trainControl(method = "repeatedcv", repeats = 5, number = 10)### Create a model formula that can be used repeatedlymodForm <- paste("CompressiveStrength ~ (.)^2 + I(Cement^2) + I(BlastFurnaceSlag^2) +", "I(FlyAsh^2) + I(Water^2) + I(Superplasticizer^2) +", "I(CoarseAggregate^2) + I(FineAggregate^2) + I(Age^2)")modForm <- as.formula(modForm)### lm & pls: Partial Least Squares set.seed(669)lmFit <- train(modForm, data = training, method = "lm", trControl = ctrl)set.seed(669)plsFit <- train(modForm, data = training, method = "pls", preProc = c("center", "scale"), tuneLength = 15, trControl = ctrl)Attaching package: 'pls'The following object is masked from 'package:caret': R2The following object is masked from 'package:stats': loadings### ridge regressionlassoGrid <- expand.grid(lambda = c(0, .001, .01, .1), fraction = seq(0.05, 1, length = 20))set.seed(669)lassoFit <- train(modForm, data = training, method = "enet", preProc = c("center", "scale"), tuneGrid = lassoGrid, trControl = ctrl)Loading required package: elasticnet### Multivariate Adaptive Regression Splinesset.seed(669)earthFit <- train(CompressiveStrength ~ ., data = training, method = "earth", tuneGrid = expand.grid(degree = 1, nprune = 2:25), trControl = ctrl)Loading required package: earthLoading required package: plotmoLoading required package: plotrixLoading required package: TeachingDemosWarning message:"package 'TeachingDemos' was built under R version 3.2.5"### SVMset.seed(669)svmRFit <- train(CompressiveStrength ~ ., data = training, method = "svmRadial", tuneLength = 15, preProc = c("center", "scale"), trControl = ctrl)Loading required package: kernlabAttaching package: 'kernlab'The following object is masked from 'package:ggplot2': alpha### neural network models### don't run it due that it'll take more times!nnetGrid <- expand.grid(decay = c(0.001, .01, .1), size = seq(1, 27, by = 2), bag = FALSE)set.seed(669)nnetFit <- train(CompressiveStrength ~ ., data = training, method = "avNNet", tuneGrid = nnetGrid, preProc = c("center", "scale"), linout = TRUE, trace = FALSE, maxit = 50, allowParallel = FALSE, trControl = ctrl)Loading required package: nnetset.seed(669)rpartFit <- train(CompressiveStrength ~ ., data = training, method = "rpart", tuneLength = 30, trControl = ctrl)Loading required package: rpartset.seed(669)treebagFit <- train(CompressiveStrength ~ ., data = training, method = "treebag", trControl = ctrl)Loading required package: ipredWarning message:"package 'ipred' was built under R version 3.2.5"Loading required package: e1071set.seed(669)ctreeFit <- train(CompressiveStrength ~ ., data = training, method = "ctree", tuneLength = 10, trControl = ctrl)Loading required package: partyWarning message:"package 'party' was built under R version 3.2.5"Loading required package: gridLoading required package: mvtnormWarning message:"package 'mvtnorm' was built under R version 3.3.0"Loading required package: modeltoolsWarning message:"package 'modeltools' was built under R version 3.2.5"Loading required package: stats4Attaching package: 'modeltools'The following object is masked from 'package:kernlab': priorThe following object is masked from 'package:plyr': emptyLoading required package: strucchangeWarning message:"package 'strucchange' was built under R version 3.2.5"Loading required package: zooAttaching package: 'zoo'The following objects are masked from 'package:base': as.Date, as.Date.numericLoading required package: sandwichWarning message:"package 'sandwich' was built under R version 3.2.5"set.seed(669)rfFit <- train(CompressiveStrength ~ ., data = training, method = "rf", tuneLength = 10, ntrees = 1000, importance = TRUE, trControl = ctrl)Loading required package: randomForestrandomForest 4.6-12Type rfNews() to see new features/changes/bug fixes.Attaching package: 'randomForest'The following object is masked from 'package:ggplot2': marginnote: only 7 unique complexity parameters in default grid. Truncating the grid to 7 .cbGrid <- expand.grid(committees = c(1, 5, 10, 50, 75, 100), neighbors = c(0, 1, 3, 5, 7, 9))set.seed(669)cbFit <- train(CompressiveStrength ~ ., data = training, method = "cubist", tuneGrid = cbGrid, trControl = ctrl)Loading required package: CubistWarning message:"package 'Cubist' was built under R version 3.2.5"rs <- resamples(list("Linear Reg" = lmFit, " PLS" = plsFit, "Elastic Net" = lassoFit, MARS = earthFit, SVM = svmRFit, CART = rpartFit, "Cond Inf Tree" = ctreeFit, "Bagged Tree" = treebagFit, "Random Forest" = rfFit, Cubist = cbFit))parallelplot(rs)
parallelplot(rs, metric = "Rsquared")
library(proxy)Attaching package: 'proxy'The following objects are masked from 'package:stats': as.dist, distThe following object is masked from 'package:base': as.matrix### Create a function to maximize compressive strength* while keeping### the predictor values as mixtures. Water (in x[7]) is used as the ### 'slack variable'. ### * We are actually minimizing the negative compressive strengthmodelPrediction <- function(x, mod, limit = 2500){ if(x[1] < 0 | x[1] > 1) return(10^38) if(x[2] < 0 | x[2] > 1) return(10^38) if(x[3] < 0 | x[3] > 1) return(10^38) if(x[4] < 0 | x[4] > 1) return(10^38) if(x[5] < 0 | x[5] > 1) return(10^38) if(x[6] < 0 | x[6] > 1) return(10^38) x <- c(x, 1 - sum(x)) if(x[7] < 0.05) return(10^38) tmp <- as.data.frame(t(x)) names(tmp) <- c('Cement','BlastFurnaceSlag','FlyAsh', 'Superplasticizer','CoarseAggregate', 'FineAggregate', 'Water') tmp$Age <- 28 -predict(mod, tmp)}### Get mixtures at 28 days subTrain <- subset(training, Age == 28)### Center and scale the data to use dissimilarity samplingpp1 <- preProcess(subTrain[, -(8:9)], c("center", "scale"))scaledTrain <- predict(pp1, subTrain[, 1:7])scaledTrain### Randomly select a few mixtures as a starting poolset.seed(91)startMixture <- sample(1:nrow(subTrain), 1)starters <- scaledTrain[startMixture, 1:7]pool <- scaledTrainindex <- maxDissim(starters, pool, 14)startPoints <- c(startMixture, index)starters <- subTrain[startPoints,1:7]startingValues <- starters[, -4]### For each starting mixture, optimize the Cubist model using### a simplex search routinecbResults <- startingValuescbResults$Water <- NAcbResults$Prediction <- NAfor(i in 1:nrow(cbResults)){ results <- optim(unlist(cbResults[i,1:6]), modelPrediction, method = "Nelder-Mead", control=list(maxit=5000), mod = cbFit) cbResults$Prediction[i] <- -results$value cbResults[i,1:6] <- results$par}cbResults$Water <- 1 - apply(cbResults[,1:6], 1, sum)cbResults <- subset(cbResults, Prediction > 0 & Water > .02)cbResults <- cbResults[order(-cbResults$Prediction),][1:3,]cbResults$Model <- "Cubist"### Do the same for the neural network modelnnetResults <- startingValuesnnetResults$Water <- NAnnetResults$Prediction <- NAfor(i in 1:nrow(nnetResults)){ results <- optim(unlist(nnetResults[i, 1:6,]), modelPrediction, method = "Nelder-Mead", control=list(maxit=5000), mod = nnetFit) nnetResults$Prediction[i] <- -results$value nnetResults[i,1:6] <- results$par}nnetResults$Water <- 1 - apply(nnetResults[,1:6], 1, sum)nnetResults <- subset(nnetResults, Prediction > 0 & Water > .02)nnetResults <- nnetResults[order(-nnetResults$Prediction),][1:3,]nnetResults$Model <- "NNet"nnetResults### Convert the predicted mixtures to PCA space and plotpp2 <- preProcess(subTrain[, 1:7], "pca")pca1 <- predict(pp2, subTrain[, 1:7])pca1$Data <- "Training Set"pca1$Data[startPoints] <- "Starting Values"pca3 <- predict(pp2, cbResults[, names(subTrain[, 1:7])])pca3$Data <- "Cubist"pca4 <- predict(pp2, nnetResults[, names(subTrain[, 1:7])])pca4$Data <- "Neural Network"pcaData <- rbind(pca1, pca3, pca4)pcaData$Data <- factor(pcaData$Data, levels = c("Training Set","Starting Values", "Cubist","Neural Network"))lim <- extendrange(pcaData[, 1:2])xyplot(PC2 ~ PC1, data = pcaData, groups = Data, auto.key = list(columns = 2), xlim = lim, ylim = lim, type = c("g", "p"))
xyplot(PC2 ~ PC1, data = pcaData, groups = Data, auto.key = list(columns = 2), par.settings = simpleTheme(pch = 20:25), xlim = lim, ylim = lim, type = c("g", "p"))
Section 11.1 Class Predictions
Basic: Simulate some two class data with two predictors
library(AppliedPredictiveModeling)### Simulate some two class data with two predictorsset.seed(975)training <- quadBoundaryFunc(500)testing <- quadBoundaryFunc(1000)testing$class2 <- ifelse(testing$class == "Class1", 1, 0)testing$ID <- 1:nrow(testing)### Fit modelslibrary(MASS)qdaFit <- qda(class ~ X1 + X2, data = training)library(randomForest)rfFit <- randomForest(class ~ X1 + X2, data = training, ntree = 2000)### Predict the test settesting$qda <- predict(qdaFit, testing)$posterior[,1]testing$rf <- predict(rfFit, testing, type = "prob")[,1]One way to assess the quality of the class probabilities is using a calibration plot.
For a given set of data, this plot shows some measure of the observed
probability of an event versus the predicted class probability. One approach
for creating this visualization is to score a collection of samples with known
outcomes (preferably a test set) using a classification model. The next step
is to bin the data into groups based on their class probabilities. For example,
a set of bins might be [0, 10%], (10%, 20%], …, (90%, 100%]. For each
bin, determine the observed event rate. Suppose that 50 samples fell into
the bin for class probabilities less than 10% and there was a single event.
The midpoint of the bin is 5% and the observed event rate would be 2 %.
The calibration plot would display the midpoint of the bin on the x-axis and
the observed event rate on the y-axis. If the points fall along a 45◦ line, the
model has produced well-calibrated probabilities.
### Generate the calibration analysislibrary(caret)calData1 <- calibration(class ~ qda + rf, data = testing, cuts = 10)### Plot the curvexyplot(calData1, auto.key = list(columns = 2))
### To calibrate the data, treat the probabilities as inputs into the### modeltrainProbs <- trainingtrainProbs$qda <- predict(qdaFit)$posterior[,1]### These models take the probabilities as inputs and, based on the### true class, re-calibrate them.library(klaR)nbCal <- NaiveBayes(class ~ qda, data = trainProbs, usekernel = TRUE)### We use relevel() here because glm() models the probability of the### second factor level.lrCal <- glm(relevel(class, "Class2") ~ qda, data = trainProbs, family = binomial)Warning message:"package 'klaR' was built under R version 3.2.5"Attaching package: 'klaR'The following object is masked from 'package:TeachingDemos': triplot### Now re-predict the test set using the modified class probability### estimatestesting$qda2 <- predict(nbCal, testing[, "qda", drop = FALSE])$posterior[,1]testing$qda3 <- predict(lrCal, testing[, "qda", drop = FALSE], type = "response")### Manipulate the data a bit for pretty plottingsimulatedProbs <- testing[, c("class", "rf", "qda3")]names(simulatedProbs) <- c("TrueClass", "RandomForestProb", "QDACalibrated")simulatedProbs$RandomForestClass <- predict(rfFit, testing)Manipulate the data a bit for pretty plotting
calData2 <- calibration(class ~ qda + qda2 + qda3, data = testing)calData2$data$calibModelVar <- as.character(calData2$data$calibModelVar)calData2$data$calibModelVar <- ifelse(calData2$data$calibModelVar == "qda", "QDA", calData2$data$calibModelVar)calData2$data$calibModelVar <- ifelse(calData2$data$calibModelVar == "qda2", "Bayesian Calibration", calData2$data$calibModelVar)calData2$data$calibModelVar <- ifelse(calData2$data$calibModelVar == "qda3", "Sigmoidal Calibration", calData2$data$calibModelVar)calData2$data$calibModelVar <- factor(calData2$data$calibModelVar, levels = c("QDA", "Bayesian Calibration", "Sigmoidal Calibration"))QDA is re-calibrated use Bayesian & Sigmoidal
xyplot(calData2, auto.key = list(columns = 1))
Presenting Class Probabilities
### Recreate the model used in the over-fitting chapterlibrary(caret)data(GermanCredit)head(GermanCredit,2)## First, remove near-zero variance predictors then get rid of a few predictors ## that duplicate values. For example, there are two possible values for the ## housing variable: "Rent", "Own" and "ForFree". So that we don't have linear## dependencies, we get rid of one of the levels (e.g. "ForFree")GermanCredit <- GermanCredit[, -nearZeroVar(GermanCredit)]GermanCredit$CheckingAccountStatus.lt.0 <- NULLGermanCredit$SavingsAccountBonds.lt.100 <- NULLGermanCredit$EmploymentDuration.lt.1 <- NULLGermanCredit$EmploymentDuration.Unemployed <- NULLGermanCredit$Personal.Male.Married.Widowed <- NULLGermanCredit$Property.Unknown <- NULLGermanCredit$Housing.ForFree <- NULLnames(GermanCredit)[ncol(GermanCredit) ]<- 'Class'## Split the data into training (80%) and test sets (20%)set.seed(100)inTrain <- createDataPartition(GermanCredit$Class, p = .8)[[1]]GermanCreditTrain <- GermanCredit[ inTrain, ]GermanCreditTest <- GermanCredit[-inTrain, ]set.seed(1056)logisticReg <- train(Class ~ ., data = GermanCreditTrain, method = "glm", trControl = trainControl(method = "repeatedcv", repeats = 5))Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"logisticRegGeneralized Linear Model 800 samples 19 predictor 2 classes: 'good', 'bad' No pre-processingResampling: Cross-Validated (10 fold, repeated 5 times) Summary of sample sizes: 720, 720, 720, 720, 720, 720, ... Resampling results: Accuracy Kappa 0.75425 0.3797266### Predict the test setcreditResults <- data.frame(obs = GermanCreditTest$Class)creditResults$prob <- predict(logisticReg, GermanCreditTest, type = "prob")[, "bad"]creditResults$pred <- predict(logisticReg, GermanCreditTest)creditResults$Label <- ifelse(creditResults$obs == "bad", "True Outcome: Bad Credit", "True Outcome: Good Credit")Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"Warning message in predict.lm(object, newdata, se.fit, scale = 1, type = ifelse(type == :"prediction from a rank-deficient fit may be misleading"### Plot the probability of bad credithistogram(~prob|Label, data = creditResults, layout = c(2, 1), nint = 20, xlab = "Probability of Bad Credit", type = "count")### Calculate and plot the calibration curvecreditCalib <- calibration(obs ~ prob, data = creditResults)xyplot(creditCalib)
Top: Histograms for a set of probabilities associated with bad credit.
The two panels split the customers by their true class. Bottom: A calibration
plot for these probabilities
The top panel shows histograms of
the test set probabilities for the logistic regression model (the panels indicate
the true credit status). The probability of bad credit for the customers with
good credit shows a skewed distribution where most customers’ probabilities
are quite low. In contrast, the probabilities for the customers with bad
credit are flat (or uniformly distributed), reflecting the model’s inability to
distinguish bad credit cases.
### Create the confusion matrix from the test set.confusionMatrix(data = creditResults$pred, reference = creditResults$obs)Confusion Matrix and Statistics ReferencePrediction good bad good 123 31 bad 17 29 Accuracy : 0.76 95% CI : (0.6947, 0.8174) No Information Rate : 0.7 P-Value [Acc > NIR] : 0.03595 Kappa : 0.3878 Mcnemar's Test P-Value : 0.06060 Sensitivity : 0.8786 Specificity : 0.4833 Pos Pred Value : 0.7987 Neg Pred Value : 0.6304 Prevalence : 0.7000 Detection Rate : 0.6150 Detection Prevalence : 0.7700 Balanced Accuracy : 0.6810 'Positive' Class : good The Kappa statistic (also known as
Cohen’s Kappa) was originally designed to assess the agreement between two
raters (Cohen 1960). Kappa takes into account the accuracy that would be
generated simply by chance.
ROC: Receiver Operating Characteristic Curves
### ROC curves:### Like glm(), roc() treats the last level of the factor as the event### of interest so we use relevel() to change the observed class datalibrary(pROC)creditROC <- roc(relevel(creditResults$obs, "good"), creditResults$prob)coords(creditROC, "all")[,1:3]auc(creditROC)ci.auc(creditROC)0.78202380952381
- 0.71286003739189
- 0.78202380952381
- 0.851187581655729
### Note the x-axis is reversedplot(creditROC)
Lift Charts
Lift charts (Ling and Li 1998) are a visualization tool for assessing the ability
of a model to detect events in a data set with two classes.
Lift charts do just this: rank the samples by
their scores and determine the cumulative event rate as more samples are
evaluated.In the optimal case, the M highest-ranked samples would contain
all M events. When the model is non-informative, the highest-ranked X %
of the data would contain, on average, X events.
creditLift <- lift(obs ~ prob, data = creditResults)xyplot(creditLift)
