R : Copyright 2005, The R Foundation for Statistical Computing Version 2.1.1 (2005-06-20), ISBN 3-900051-07-0 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for a HTML browser interface to help. Type 'q()' to quit R. > ### *
> ### > attach(NULL, name = "CheckExEnv") > assign(".CheckExEnv", as.environment(2), pos = length(search())) # base > ## add some hooks to label plot pages for base and grid graphics > setHook("plot.new", ".newplot.hook") > setHook("persp", ".newplot.hook") > setHook("grid.newpage", ".gridplot.hook") > > assign("cleanEx", + function(env = .GlobalEnv) { + rm(list = ls(envir = env, all.names = TRUE), envir = env) + RNGkind("default", "default") + set.seed(1) + options(warn = 1) + delayedAssign("T", stop("T used instead of TRUE"), + assign.env = .CheckExEnv) + delayedAssign("F", stop("F used instead of FALSE"), + assign.env = .CheckExEnv) + sch <- search() + newitems <- sch[! sch %in% .oldSearch] + for(item in rev(newitems)) + eval(substitute(detach(item), list(item=item))) + missitems <- .oldSearch[! .oldSearch %in% sch] + if(length(missitems)) + warning("items ", paste(missitems, collapse=", "), + " have been removed from the search path") + }, + env = .CheckExEnv) > assign("..nameEx", "__{must remake R-ex/*.R}__", env = .CheckExEnv) # for now > assign("ptime", proc.time(), env = .CheckExEnv) > grDevices::postscript("ipred-Examples.ps") > assign("par.postscript", graphics::par(no.readonly = TRUE), env = .CheckExEnv) > options(contrasts = c(unordered = "contr.treatment", ordered = "contr.poly")) > options(warn = 1) > library('ipred') Loading required package: rpart Loading required package: MASS Loading required package: mlbench Loading required package: survival Loading required package: splines Loading required package: nnet Loading required package: class > > assign(".oldSearch", search(), env = .CheckExEnv) > assign(".oldNS", loadedNamespaces(), env = .CheckExEnv) > cleanEx(); ..nameEx <- "DLBCL" > > ### * DLBCL > > flush(stderr()); flush(stdout()) > > ### Name: DLBCL > ### Title: Diffuse Large B-Cell Lymphoma > ### Aliases: DLBCL > ### Keywords: datasets > > ### ** Examples > > data(DLBCL) > survfit(Surv(time, cens), data=DLBCL) Call: survfit(formula = Surv(time, cens), data = DLBCL) n events median 0.95LCL 0.95UCL 40.0 22.0 32.5 15.5 Inf > > > > > cleanEx(); ..nameEx <- "GBSG2" > > ### * GBSG2 > > flush(stderr()); flush(stdout()) > > ### Name: GBSG2 > ### Title: German Breast Cancer Study Group 2 > ### Aliases: GBSG2 > ### Keywords: datasets > > ### ** Examples > > data(GBSG2) > > thsum <- function(x) { + ret <- c(median(x), quantile(x, 0.25), quantile(x,0.75)) + names(ret)[1] <- "Median" + ret + } > > t(apply(GBSG2[,c("age", "tsize", "pnodes", + "progrec", "estrec")], 2, thsum)) Median 25% 75% age 53.0 46 61.00 tsize 25.0 20 35.00 pnodes 3.0 1 7.00 progrec 32.5 7 131.75 estrec 36.0 8 114.00 > > table(GBSG2$menostat) Pre Post 290 396 > table(GBSG2$tgrade) I II III 81 444 161 > table(GBSG2$horTh) no yes 440 246 > > # pooled Kaplan-Meier > > mod <- survfit(Surv(time, cens), data=GBSG2) > # integrated Brier score > sbrier(Surv(GBSG2$time, GBSG2$cens), mod) integrated Brier score 0.1919286 attr(,"time") [1] 8 2659 > # Brier score at 5 years > sbrier(Surv(GBSG2$time, GBSG2$cens), mod, btime=1825) Brier score 0.2499263 attr(,"time") [1] 1825 > > # Nottingham prognostic index > > GBSG2 <- GBSG2[order(GBSG2$time),] > > NPI <- 0.2*GBSG2$tsize/10 + 1 + as.integer(GBSG2$tgrade) > NPI[NPI < 3.4] <- 1 > NPI[NPI >= 3.4 & NPI <=5.4] <- 2 > NPI[NPI > 5.4] <- 3 > > mod <- survfit(Surv(time, cens) ~ NPI, data=GBSG2) > plot(mod) > > pred <- c() > survs <- c() > for (i in sort(unique(NPI))) + survs <- c(survs, getsurv(mod[i], 1825)) > > for (i in 1:nrow(GBSG2)) + pred <- c(pred, survs[NPI[i]]) > > # Brier score of NPI at t=5 years > sbrier(Surv(GBSG2$time, GBSG2$cens), pred, btime=1825) Brier score 0.2337553 attr(,"time") [1] 1825 > > > > > cleanEx(); ..nameEx <- "GlaucomaM" > > ### * GlaucomaM > > flush(stderr()); flush(stdout()) > > ### Name: GlaucomaM > ### Title: Glaucoma Database > ### Aliases: GlaucomaM > ### Keywords: datasets > > ### ** Examples > > data(GlaucomaM) > errorest(Class ~ ., data=GlaucomaM, model=rpart, + predict=function(obj, newdata) + predict(obj, newdata, type="class"), + control=rpart.control(xval=0)) Call: errorest.data.frame(formula = Class ~ ., data = GlaucomaM, model = rpart, predict = function(obj, newdata) predict(obj, newdata, type = "class"), control = rpart.control(xval = 0)) 10-fold cross-validation estimator of misclassification error Misclassification error: 0.2092 > glbagg <- bagging(Class ~ ., data=GlaucomaM, coob=TRUE) > glbagg Bagging classification trees with 25 bootstrap replications Call: bagging.data.frame(formula = Class ~ ., data = GlaucomaM, coob = TRUE) Out-of-bag estimate of misclassification error: 0.1888 > > > > > cleanEx(); ..nameEx <- "GlaucomaMVF" > > ### * GlaucomaMVF > > flush(stderr()); flush(stdout()) > > ### Name: GlaucomaMVF > ### Title: Glaucoma Database > ### Aliases: GlaucomaMVF > ### Keywords: datasets > > ### ** Examples > > ## Not run: > ##D data(GlaucomaMVF) > ##D > ##D response <- function (data) { > ##D attach(data) > ##D res <- ifelse((!is.na(clv) & !is.na(lora) & clv >= 5.1 & lora >= > ##D 49.23372) | (!is.na(clv) & !is.na(lora) & !is.na(cs) & > ##D clv < 5.1 & lora >= 58.55409 & cs < 1.405) | (is.na(clv) & > ##D !is.na(lora) & !is.na(cs) & lora >= 58.55409 & cs < 1.405) | > ##D (!is.na(clv) & is.na(lora) & cs < 1.405), 0, 1) > ##D detach(data) > ##D factor (res, labels = c("glaucoma", "normal")) > ##D } > ##D > ##D errorest(Class~clv+lora+cs~., data = GlaucomaMVF, model=inclass, > ##D estimator="cv", pFUN = list(list(model = rpart)), cFUN = response) > ## End(Not run) > > > > cleanEx(); ..nameEx <- "bagging" > > ### * bagging > > flush(stderr()); flush(stdout()) > > ### Name: bagging > ### Title: Bagging Classification, Regression and Survival Trees > ### Aliases: bagging ipredbagg ipredbagg.factor ipredbagg.integer > ### ipredbagg.numeric ipredbagg.Surv ipredbagg.default bagging.data.frame > ### bagging.default > ### Keywords: tree > > ### ** Examples > > > # Classification: Breast Cancer data > > data(BreastCancer) > > # Test set error bagging (nbagg = 50): 3.7% (Breiman, 1998, Table 5) > > mod <- bagging(Class ~ Cl.thickness + Cell.size + + Cell.shape + Marg.adhesion + + Epith.c.size + Bare.nuclei + + Bl.cromatin + Normal.nucleoli + + Mitoses, data=BreastCancer, coob=TRUE) > print(mod) Bagging classification trees with 25 bootstrap replications Call: bagging.data.frame(formula = Class ~ Cl.thickness + Cell.size + Cell.shape + Marg.adhesion + Epith.c.size + Bare.nuclei + Bl.cromatin + Normal.nucleoli + Mitoses, data = BreastCancer, coob = TRUE) Out-of-bag estimate of misclassification error: 0.0395 > > # Test set error bagging (nbagg=50): 7.9% (Breiman, 1996a, Table 2) > data(Ionosphere) > Ionosphere$V2 <- NULL # constant within groups > > bagging(Class ~ ., data=Ionosphere, coob=TRUE) Bagging classification trees with 25 bootstrap replications Call: bagging.data.frame(formula = Class ~ ., data = Ionosphere, coob = TRUE) Out-of-bag estimate of misclassification error: 0.0912 > > # Double-Bagging: combine LDA and classification trees > > # predict returns the linear discriminant values, i.e. linear combinations > # of the original predictors > > comb.lda <- list(list(model=lda, predict=function(obj, newdata) + predict(obj, newdata)$x)) > > # Note: out-of-bag estimator is not available in this situation, use > # errorest > > mod <- bagging(Class ~ ., data=Ionosphere, comb=comb.lda) > > predict(mod, Ionosphere[1:10,]) [1] good bad good bad good bad good bad good bad Levels: bad good > > # Regression: > > data(BostonHousing) > > # Test set error (nbagg=25, trees pruned): 3.41 (Breiman, 1996a, Table 8) > > mod <- bagging(medv ~ ., data=BostonHousing, coob=TRUE) > print(mod) Bagging regression trees with 25 bootstrap replications Call: bagging.data.frame(formula = medv ~ ., data = BostonHousing, coob = TRUE) Out-of-bag estimate of root mean squared error: 4.1106 > > learn <- as.data.frame(mlbench.friedman1(200)) > > # Test set error (nbagg=25, trees pruned): 2.47 (Breiman, 1996a, Table 8) > > mod <- bagging(y ~ ., data=learn, coob=TRUE) > print(mod) Bagging regression trees with 25 bootstrap replications Call: bagging.data.frame(formula = y ~ ., data = learn, coob = TRUE) Out-of-bag estimate of root mean squared error: 2.8989 > > # Survival data > > # Brier score for censored data estimated by > # 10 times 10-fold cross-validation: 0.2 (Hothorn et al, > # 2002) > > data(DLBCL) > mod <- bagging(Surv(time,cens) ~ MGEc.1 + MGEc.2 + MGEc.3 + MGEc.4 + MGEc.5 + + MGEc.6 + MGEc.7 + MGEc.8 + MGEc.9 + + MGEc.10 + IPI, data=DLBCL, coob=TRUE) > > print(mod) Bagging survival trees with 25 bootstrap replications Call: bagging.data.frame(formula = Surv(time, cens) ~ MGEc.1 + MGEc.2 + MGEc.3 + MGEc.4 + MGEc.5 + MGEc.6 + MGEc.7 + MGEc.8 + MGEc.9 + MGEc.10 + IPI, data = DLBCL, coob = TRUE) Out-of-bag estimate of Brier's score: 0.2119 > > > > > cleanEx(); ..nameEx <- "dystrophy" > > ### * dystrophy > > flush(stderr()); flush(stdout()) > > ### Name: dystrophy > ### Title: Detection of muscular dystrophy carriers. > ### Aliases: dystrophy > ### Keywords: datasets > > ### ** Examples > > ## Not run: > ##D data(dystrophy) > ##D errorest(Class~CK+H~AGE+PK+LD, data = dystrophy, model = inbagg, > ##D pFUN = list(list(model = lm, predict = mypredict.lm), list(model = rpart)), > ##D ns = 0.75, estimator = "cv") > ## End(Not run) > > > > cleanEx(); ..nameEx <- "errorest" > > ### * errorest > > flush(stderr()); flush(stdout()) > > ### Name: errorest > ### Title: Estimators of Prediction Error > ### Aliases: errorest errorest.data.frame errorest.default > ### Keywords: misc > > ### ** Examples > > > # Classification > > data(iris) > > # force predict to return class labels only > mypredict.lda <- function(object, newdata) + predict(object, newdata = newdata)$class > > # 10-fold cv of LDA for Iris data > errorest(Species ~ ., data=iris, model=lda, + estimator = "cv", predict= mypredict.lda) Call: errorest.data.frame(formula = Species ~ ., data = iris, model = lda, predict = mypredict.lda, estimator = "cv") 10-fold cross-validation estimator of misclassification error Misclassification error: 0.02 > > data(PimaIndiansDiabetes) > ## Not run: > ##D # 632+ bootstrap of LDA for Diabetes data > ##D errorest(diabetes ~ ., data=PimaIndiansDiabetes, model=lda, > ##D estimator = "632plus", predict= mypredict.lda) > ## End(Not run) > > #cv of a fixed partition of the data > list.tindx <- list(1:100, 101:200, 201:300, 301:400, 401:500, + 501:600, 601:700, 701:768) > > errorest(diabetes ~ ., data=PimaIndiansDiabetes, model=lda, + estimator = "cv", predict = mypredict.lda, + est.para = control.errorest(list.tindx = list.tindx)) Call: errorest.data.frame(formula = diabetes ~ ., data = PimaIndiansDiabetes, model = lda, predict = mypredict.lda, estimator = "cv", est.para = control.errorest(list.tindx = list.tindx)) 8-fold cross-validation estimator of misclassification error Misclassification error: 0.2227 > > ## Not run: > ##D #both bootstrap estimations based on fixed partitions > ##D > ##D list.tindx <- vector(mode = "list", length = 25) > ##D for(i in 1:25) { > ##D list.tindx[[i]] <- sample(1:768, 768, TRUE) > ##D } > ##D > ##D errorest(diabetes ~ ., data=PimaIndiansDiabetes, model=lda, > ##D estimator = c("boot", "632plus"), predict= mypredict.lda, > ##D est.para = control.errorest(list.tindx = list.tindx)) > ##D > ## End(Not run) > data(Glass) > > # LDA has cross-validated misclassification error of > # 38% (Ripley, 1996, page 98) > > # Pruned trees about 32% (Ripley, 1996, page 230) > > # use stratified sampling here, i.e. preserve the class proportions > errorest(Type ~ ., data=Glass, model=lda, + predict=mypredict.lda, est.para=control.errorest(strat=TRUE)) Call: errorest.data.frame(formula = Type ~ ., data = Glass, model = lda, predict = mypredict.lda, est.para = control.errorest(strat = TRUE)) 10-fold cross-validation estimator of misclassification error Misclassification error: 0.3785 > > # force predict to return class labels > mypredict.rpart <- function(object, newdata) + predict(object, newdata = newdata,type="class") > > pruneit <- function(formula, ...) + prune(rpart(formula, ...), cp =0.01) > > errorest(Type ~ ., data=Glass, model=pruneit, + predict=mypredict.rpart, est.para=control.errorest(strat=TRUE)) Call: errorest.data.frame(formula = Type ~ ., data = Glass, model = pruneit, predict = mypredict.rpart, est.para = control.errorest(strat = TRUE)) 10-fold cross-validation estimator of misclassification error Misclassification error: 0.3178 > > # compute sensitivity and specifity for stabilised LDA > > data(GlaucomaM) > > error <- errorest(Class ~ ., data=GlaucomaM, model=slda, + predict=mypredict.lda, est.para=control.errorest(predictions=TRUE)) > > # sensitivity > > mean(error$predictions[GlaucomaM$Class == "glaucoma"] == "glaucoma") [1] 0.8163265 > > # specifity > > mean(error$predictions[GlaucomaM$Class == "normal"] == "normal") [1] 0.8367347 > > # Indirect Classification: Smoking data > > data(Smoking) > # Set three groups of variables: > # 1) explanatory variables are: TarY, NicY, COY, Sex, Age > # 2) intermediate variables are: TVPS, BPNL, COHB > # 3) response (resp) is defined by: > > resp <- function(data){ + data <- data[, c("TVPS", "BPNL", "COHB")] + res <- t(t(data) > c(4438, 232.5, 58)) + res <- as.factor(ifelse(apply(res, 1, sum) > 2, 1, 0)) + res + } > > response <- resp(Smoking[ ,c("TVPS", "BPNL", "COHB")]) > smoking <- cbind(Smoking, response) > > formula <- response~TVPS+BPNL+COHB~TarY+NicY+COY+Sex+Age > > # Estimation per leave-one-out estimate for the misclassification is > # 36.36% (Hand et al., 2001), using indirect classification with > # linear models > ## Not run: > ##D errorest(formula, data = smoking, model = inclass,estimator = "cv", > ##D pFUN = list(list(model=lm, predict = mypredict.lm)), cFUN = resp, > ##D est.para=control.errorest(k=nrow(smoking))) > ## End(Not run) > > # Regression > > data(BostonHousing) > > # 10-fold cv of lm for Boston Housing data > errorest(medv ~ ., data=BostonHousing, model=lm, + est.para=control.errorest(random=FALSE)) Call: errorest.data.frame(formula = medv ~ ., data = BostonHousing, model = lm, est.para = control.errorest(random = FALSE)) 10-fold cross-validation estimator of root mean squared error Root mean squared error: 5.877 > > # the same, with "model" returning a function for prediction > # instead of an object of class "lm" > > mylm <- function(formula, data) { + mod <- lm(formula, data) + function(newdata) predict(mod, newdata) + } > > errorest(medv ~ ., data=BostonHousing, model=mylm, + est.para=control.errorest(random=FALSE)) Call: errorest.data.frame(formula = medv ~ ., data = BostonHousing, model = mylm, est.para = control.errorest(random = FALSE)) 10-fold cross-validation estimator of root mean squared error Root mean squared error: 5.877 > > # Survival data > > data(GBSG2) > > # prediction is fitted Kaplan-Meier > predict.survfit <- function(object, newdata) object > > # 5-fold cv of Kaplan-Meier for GBSG2 study > errorest(Surv(time, cens) ~ 1, data=GBSG2, model=survfit, + predict=predict.survfit, est.para=control.errorest(k=5)) Call: errorest.data.frame(formula = Surv(time, cens) ~ 1, data = GBSG2, model = survfit, predict = predict.survfit, est.para = control.errorest(k = 5)) 5-fold cross-validation estimator of Brier's score Brier's score: 0.1882 > > > > > cleanEx(); ..nameEx <- "inbagg" > > ### * inbagg > > flush(stderr()); flush(stdout()) > > ### Name: inbagg > ### Title: Indirect Bagging > ### Aliases: inbagg inbagg.default inbagg.data.frame > ### Keywords: misc > > ### ** Examples > > y <- as.factor(sample(1:2, 100, replace = TRUE)) > W <- mvrnorm(n = 200, mu = rep(0, 3), Sigma = diag(3)) > X <- mvrnorm(n = 200, mu = rep(2, 3), Sigma = diag(3)) > colnames(W) <- c("w1", "w2", "w3") > colnames(X) <- c("x1", "x2", "x3") > DATA <- data.frame(y, W, X) > > pFUN <- list(list(formula = w1~x1+x2, model = lm, predict = mypredict.lm), + list(model = rpart)) > > inbagg(y~w1+w2+w3~x1+x2+x3, data = DATA, pFUN = pFUN) Indirect bagging, with 25 bootstrap samples and intermediate variables: w1 w2 w3 > > > > cleanEx(); ..nameEx <- "inclass" > > ### * inclass > > flush(stderr()); flush(stdout()) > > ### Name: inclass > ### Title: Indirect Classification > ### Aliases: inclass inclass.default inclass.data.frame > ### Keywords: misc > > ### ** Examples > > data(Smoking) > # Set three groups of variables: > # 1) explanatory variables are: TarY, NicY, COY, Sex, Age > # 2) intermediate variables are: TVPS, BPNL, COHB > # 3) response (resp) is defined by: > > classify <- function(data){ + data <- data[,c("TVPS", "BPNL", "COHB")] + res <- t(t(data) > c(4438, 232.5, 58)) + res <- as.factor(ifelse(apply(res, 1, sum) > 2, 1, 0)) + res + } > > response <- classify(Smoking[ ,c("TVPS", "BPNL", "COHB")]) > smoking <- data.frame(Smoking, response) > > formula <- response~TVPS+BPNL+COHB~TarY+NicY+COY+Sex+Age > > inclass(formula, data = smoking, pFUN = list(list(model = lm, predict = + mypredict.lm)), cFUN = classify) Indirect classification, with 3 intermediate variables: TVPS BPNL COHB Predictive model per intermediate is lm > > > > > cleanEx(); ..nameEx <- "ipredknn" > > ### * ipredknn > > flush(stderr()); flush(stdout()) > > ### Name: ipredknn > ### Title: k-Nearest Neighbour Classification > ### Aliases: ipredknn > ### Keywords: multivariate > > ### ** Examples > > > learn <- as.data.frame(mlbench.twonorm(300)) > > mypredict.knn <- function(object, newdata) + predict.ipredknn(object, newdata, type="class") > > errorest(classes ~., data=learn, model=ipredknn, + predict=mypredict.knn) Call: errorest.data.frame(formula = classes ~ ., data = learn, model = ipredknn, predict = mypredict.knn) 10-fold cross-validation estimator of misclassification error Misclassification error: 0.0533 > > > > > cleanEx(); ..nameEx <- "kfoldcv" > > ### * kfoldcv > > flush(stderr()); flush(stdout()) > > ### Name: kfoldcv > ### Title: Subsamples for k-fold Cross-Validation > ### Aliases: kfoldcv > ### Keywords: misc > > ### ** Examples > > > # 10-fold CV with N = 91 > > kfoldcv(10, 91) [1] 10 9 9 9 9 9 9 9 9 9 > > ## Don't show: > k <- sample(5:15, 1) > k [1] 7 > N <- sample(50:150, 1) > N [1] 87 > stopifnot(sum(kfoldcv(k, N)) == N) > ## End Don't show > > > > > cleanEx(); ..nameEx <- "predict.bagging" > > ### * predict.bagging > > flush(stderr()); flush(stdout()) > > ### Name: predict.classbagg > ### Title: Predictions from Bagging Trees > ### Aliases: predict.classbagg predict.regbagg predict.survbagg > ### Keywords: tree > > ### ** Examples > > > data(Ionosphere) > Ionosphere$V2 <- NULL # constant within groups > > # nbagg = 10 for performance reasons here > mod <- bagging(Class ~ ., data=Ionosphere) > > # out-of-bag estimate > > mean(predict(mod) != Ionosphere$Class) [1] 0.08831909 > > # predictions for the first 10 observations > > predict(mod, newdata=Ionosphere[1:10,]) [1] good bad good bad good bad good bad good bad Levels: bad good > > predict(mod, newdata=Ionosphere[1:10,], type="prob") bad good [1,] 0.00 1.00 [2,] 1.00 0.00 [3,] 0.00 1.00 [4,] 0.64 0.36 [5,] 0.04 0.96 [6,] 1.00 0.00 [7,] 0.00 1.00 [8,] 1.00 0.00 [9,] 0.00 1.00 [10,] 1.00 0.00 > > > > > cleanEx(); ..nameEx <- "predict.inbagg" > > ### * predict.inbagg > > flush(stderr()); flush(stdout()) > > ### Name: predict.inbagg > ### Title: Predictions from an Inbagg Object > ### Aliases: predict.inbagg > ### Keywords: misc > > ### ** Examples > > library(mvtnorm) > y <- as.factor(sample(1:2, 100, replace = TRUE)) > W <- mvrnorm(n = 200, mu = rep(0, 3), Sigma = diag(3)) > X <- mvrnorm(n = 200, mu = rep(2, 3), Sigma = diag(3)) > colnames(W) <- c("w1", "w2", "w3") > colnames(X) <- c("x1", "x2", "x3") > DATA <- data.frame(y, W, X) > > pFUN <- list(list(formula = w1~x1+x2, model = lm), + list(model = rpart)) > > RES <- inbagg(y~w1+w2+w3~x1+x2+x3, data = DATA, pFUN = pFUN) > predict(RES, newdata = X) [1] 1 1 2 2 1 2 2 1 2 1 1 1 2 1 2 1 2 2 1 2 2 1 2 1 1 1 1 1 2 1 1 2 1 1 2 2 2 [38] 1 2 1 2 2 2 2 2 2 1 1 2 2 1 2 1 1 1 1 1 2 2 1 2 1 1 1 2 1 1 2 1 2 1 2 1 1 [75] 1 2 2 1 2 2 1 2 1 1 2 1 2 1 1 1 1 1 2 2 2 2 1 1 2 2 1 1 2 2 1 2 2 2 2 1 1 [112] 1 2 1 2 1 2 2 1 2 2 1 2 1 1 1 1 1 2 2 1 2 1 1 2 2 2 1 2 1 2 2 2 2 2 2 1 1 [149] 2 2 1 2 1 1 1 1 1 2 2 1 2 1 1 1 2 1 1 2 1 2 1 2 1 1 1 2 2 1 2 2 1 2 1 1 2 [186] 1 2 1 1 1 1 1 2 2 2 2 1 1 2 2 Levels: 1 2 > > > > cleanEx(); ..nameEx <- "predict.inclass" > > ### * predict.inclass > > flush(stderr()); flush(stdout()) > > ### Name: predict.inclass > ### Title: Predictions from an Inclass Object > ### Aliases: predict.inclass > ### Keywords: misc > > ### ** Examples > > ## Not run: > ##D # Simulation model, classification rule following Hand et al. (2001) > ##D > ##D theta90 <- varset(N = 1000, sigma = 0.1, theta = 90, threshold = 0) > ##D > ##D dataset <- as.data.frame(cbind(theta90$explanatory, theta90$intermediate)) > ##D names(dataset) <- c(colnames(theta90$explanatory), > ##D colnames(theta90$intermediate)) > ##D > ##D classify <- function(Y, threshold = 0) { > ##D Y <- Y[,c("y1", "y2")] > ##D z <- (Y > threshold) > ##D resp <- as.factor(ifelse((z[,1] + z[,2]) > 1, 1, 0)) > ##D return(resp) > ##D } > ##D > ##D formula <- response~y1+y2~x1+x2 > ##D > ##D fit <- inclass(formula, data = dataset, pFUN = list(list(model = lm)), > ##D cFUN = classify) > ##D > ##D predict(object = fit, newdata = dataset) > ##D > ##D data(Smoking) > ##D > ##D # explanatory variables are: TarY, NicY, COY, Sex, Age > ##D # intermediate variables are: TVPS, BPNL, COHB > ##D # reponse is defined by: > ##D > ##D classify <- function(data){ > ##D data <- data[,c("TVPS", "BPNL", "COHB")] > ##D res <- t(t(data) > c(4438, 232.5, 58)) > ##D res <- as.factor(ifelse(apply(res, 1, sum) > 2, 1, 0)) > ##D res > ##D } > ##D > ##D response <- classify(Smoking[ ,c("TVPS", "BPNL", "COHB")]) > ##D smoking <- cbind(Smoking, response) > ##D > ##D formula <- response~TVPS+BPNL+COHB~TarY+NicY+COY+Sex+Age > ##D > ##D fit <- inclass(formula, data = smoking, > ##D pFUN = list(list(model = lm)), cFUN = classify) > ##D > ##D predict(object = fit, newdata = smoking) > ## End(Not run) > > data(GlaucomaMVF) > glaucoma <- GlaucomaMVF[,(names(GlaucomaMVF) != "tension")] > # explanatory variables are derived by laser scanning image and intra occular pressure > # intermediate variables are: clv, cs, lora > # response is defined by > > classify <- function (data) { + attach(data) + res <- ifelse((!is.na(clv) & !is.na(lora) & clv >= 5.1 & lora >= + 49.23372) | (!is.na(clv) & !is.na(lora) & !is.na(cs) & + clv < 5.1 & lora >= 58.55409 & cs < 1.405) | (is.na(clv) & + !is.na(lora) & !is.na(cs) & lora >= 58.55409 & cs < 1.405) | + (!is.na(clv) & is.na(lora) & cs < 1.405), 0, 1) + detach(data) + factor (res, labels = c("glaucoma", "normal")) + } > > fit <- inclass(Class~clv+lora+cs~., data = glaucoma, + pFUN = list(list(model = rpart)), cFUN = classify) > > data(GlaucomaM) > predict(object = fit, newdata = GlaucomaM) [1] normal normal normal normal glaucoma glaucoma normal normal [9] normal normal normal normal glaucoma normal normal normal [17] normal normal normal glaucoma normal normal glaucoma normal [25] normal normal glaucoma normal glaucoma normal normal normal [33] normal normal normal normal glaucoma normal normal normal [41] normal normal glaucoma normal normal glaucoma normal normal [49] normal normal normal glaucoma glaucoma glaucoma glaucoma normal [57] glaucoma glaucoma normal normal glaucoma normal glaucoma normal [65] glaucoma normal normal normal normal normal glaucoma glaucoma [73] glaucoma normal normal normal glaucoma normal normal normal [81] glaucoma normal normal normal normal glaucoma glaucoma glaucoma [89] glaucoma normal normal normal glaucoma normal normal normal [97] normal normal glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma [105] normal glaucoma normal glaucoma glaucoma glaucoma glaucoma glaucoma [113] glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma [121] glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma [129] glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma [137] glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma [145] glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma [153] glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma [161] glaucoma glaucoma glaucoma normal glaucoma glaucoma glaucoma glaucoma [169] glaucoma glaucoma glaucoma glaucoma normal glaucoma glaucoma glaucoma [177] glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma glaucoma [185] glaucoma glaucoma normal glaucoma glaucoma glaucoma glaucoma glaucoma [193] glaucoma glaucoma glaucoma glaucoma Levels: glaucoma normal > > > > > cleanEx(); ..nameEx <- "prune.bagging" > > ### * prune.bagging > > flush(stderr()); flush(stdout()) > > ### Name: prune.classbagg > ### Title: Pruning for Bagging > ### Aliases: prune.classbagg prune.regbagg prune.survbagg > ### Keywords: tree > > ### ** Examples > > > data(Glass) > > mod <- bagging(Type ~ ., data=Glass, nbagg=10, coob=TRUE) > pmod <- prune(mod) Error in prune.rpart(tree$mtrees[[i]]$btree, cp = cp, ...) : subscript out of bounds Execution halted