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("cobs-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('cobs') > > assign(".oldSearch", search(), env = .CheckExEnv) > assign(".oldNS", loadedNamespaces(), env = .CheckExEnv) > cleanEx(); ..nameEx <- "DublinWind" > > ### * DublinWind > > flush(stderr()); flush(stdout()) > > ### Name: DublinWind > ### Title: Daily Wind Speeds in Dublin > ### Aliases: DublinWind > ### Keywords: datasets > > ### ** Examples > > data(DublinWind) > str(DublinWind) `data.frame': 6574 obs. of 2 variables: $ speed: num 13.67 11.50 11.25 8.63 11.92 ... $ day : int 1 2 3 4 5 6 7 8 9 10 ... > plot(speed ~ day, data = DublinWind, type = "l") > attach(DublinWind) > plot(wSpeed <- ts(speed)) > > > > > cleanEx(); ..nameEx <- "USArmyRoofs" > > ### * USArmyRoofs > > flush(stderr()); flush(stdout()) > > ### Name: USArmyRoofs > ### Title: Roof Quality in US Army Bases > ### Aliases: USArmyRoofs > ### Keywords: datasets > > ### ** Examples > > data(USArmyRoofs) > plot(fci ~ age, data = USArmyRoofs, main = "US Army Roofs data") > > > > > cleanEx(); ..nameEx <- "cobs-methods" > > ### * cobs-methods > > flush(stderr()); flush(stdout()) > > ### Name: cobs-methods > ### Title: Methods for COBS Objects > ### Aliases: print.cobs summary.cobs fitted.cobs residuals.cobs > ### Keywords: print > > ### ** Examples > > example(cobs) cobs> x <- seq(-1, 1, , 50) cobs> y <- (f.true <- pnorm(2 * x)) + rnorm(50)/10 cobs> con <- rbind(c(1, min(x), 0), c(-1, max(x), 1), c(0, 0, 0.5)) cobs> Rbs <- cobs(x, y, constraint = "increase", pointwise = con) Performing general knot selection ... N[rq,L1]= (50, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,6, vars=5) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,7, vars=6) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,8, vars=7) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): Computing the final fit ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): cobs> Rbs COBS regression spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 3]: -1.0000030, 0.1836735, 1.0000030 cobs> plot(x, y) cobs> lines(predict(Rbs), col = 2, lwd = 1.5) cobs> lines(spline(x, f.true), col = "gray40") cobs> Rbsub <- cobs(x, y, constraint = "increase", pointwise = con, n.sub = 45) Performing general knot selection ... N[rq,L1]= (45, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,6, vars=5) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,7, vars=6) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,8, vars=7) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... N[rq,L1]= (45, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): Computing the final fit ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): cobs> summary(Rbsub) COBS regression spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", n.sub = 45, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 3]: -1.0000030, 0.1836735, 1.0000030 with 3 pointwise constraints coef : [1] 0.03313316 0.18255331 0.94257663 1.00000042 R^2 = 94.37% ; empirical tau (over all): 25/50 = 0.5 (target tau : 0.5) cobs> lines(predict(Rbsub), col = 4, lwd = 1) cobs> Sbs <- cobs(x, y, constraint = "increase", pointwise = con, lambda = -1) Searching for optimal lambda. This may take a while. While you are waiting, here is something you can consider to speed up the process: (a) Use a smaller number of knots; (b) Increase `factor' to 2 or 3; (c) Set lambda==0 to exclude the penalty term. N[rq,L1]= (50, 1); ([ei]qc= 1,60, vars=22) LAMBDA will be updated; (tau,lam)= (0.5, -1): cobs> Sbs COBS smoothing spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", lambda = -1, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 20]: -1.0000010, -0.9183673, -0.7959184, ... , 1.0000010 lambda = 0.4508836, selected via SIC, out of 115 ones. cobs> plot(Sbs$pp.lambda[-1], Sbs$sic[-1], log = "x", main = "SIC ~ lambda", xlab = expression(lambda), ylab = "SIC") cobs> axis(1, at = Sbs$lambda, label = expression(hat(lambda)), col.axis = 2, mgp = c(3, 0.5, 0)) cobs> Sb1 <- cobs(x, y, constraint = "increase", pointwise = con, lambda = -1, degree = 1) Searching for optimal lambda. This may take a while. While you are waiting, here is something you can consider to speed up the process: (a) Use a smaller number of knots; (b) Increase `factor' to 2 or 3; (c) Set lambda==0 to exclude the penalty term. N[rq,L1]= (50, 18); ([ei]qc= 1,21, vars=20) LAMBDA will be updated; (tau,lam)= (0.5, -1): cobs> summary(Sb1) COBS smoothing spline (degree = 1) from call: cobs(x = x, y = y, constraint = "increase", degree = 1, lambda = -1, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 5 (5) knots[1 .. 20]: -1.0000010, -0.9183673, -0.7959184, ... , 1.0000010 lambda = 0.2580284, selected via SIC, out of 44 ones. with 3 pointwise constraints coef : [1] 0.01355796 0.04681331 0.09669571 0.12995064 0.17983305 0.21308798 [7] 0.28194687 0.35080575 0.39671167 0.46557056 0.51147648 0.59368266 [13] 0.64848678 0.73069296 0.81289914 0.86770326 0.90378419 0.92783814 [19] 0.96391907 1.00000029 R^2 = 94.34% ; empirical tau (over all): 25/50 = 0.5 (target tau : 0.5) cobs> pxx <- predict(Sb1, xx <- seq(-1.2, 1.2, len = 201), interval = "both") cobs> plot(x, y, main = deparse(Sb1$call), xlim = range(xx), ylim = range(y, pxx[, "fit"])) cobs> lines(pxx, col = 2) cobs> rug(Sb1$knots, col = 4, lwd = 1.6) cobs> matlines(pxx[, 1], pxx[, -(1:2)], col = rep(c("green", "blue"), c(2, 2)), lty = 2) > summary(Sbs) COBS smoothing spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", lambda = -1, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 20]: -1.0000010, -0.9183673, -0.7959184, ... , 1.0000010 lambda = 0.4508836, selected via SIC, out of 115 ones. with 3 pointwise constraints coef : [1] 0.001090078 0.022607597 0.071968862 0.114681790 0.152962656 0.197891616 [7] 0.247252639 0.314463579 0.372686089 0.426476536 0.486915078 0.551785682 [13] 0.623304380 0.690391016 0.762917265 0.816707713 0.866066097 0.908776388 [19] 0.947054616 0.985010931 1.000000245 0.532069132 R^2 = 94.16% ; empirical tau (over all): 24/50 = 0.48 (target tau : 0.5) > > plot(fitted(Sbs), resid(Sbs), main = deparse(Sbs$call)) > > > > cleanEx(); ..nameEx <- "cobs" > > ### * cobs > > flush(stderr()); flush(stdout()) > > ### Name: cobs > ### Title: COnstrained B-Splines Nonparametric Regression Quantiles > ### Aliases: cobs n1000cut > ### Keywords: smooth regression > > ### ** Examples > > x <- seq(-1,1,,50) > y <- (f.true <- pnorm(2*x)) + rnorm(50)/10 > ## specify pointwise constraints (boundary conditions) > con <- rbind(c( 1,min(x),0), # f(min(x)) >= 0 + c(-1,max(x),1), # f(max(x)) <= 1 + c(0, 0, 0.5))# f(0) = 0.5 > > ## obtain the median regression B-spline using automatically selected knots > Rbs <- cobs(x,y,constraint="increase",pointwise=con) Performing general knot selection ... N[rq,L1]= (50, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,6, vars=5) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,7, vars=6) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,8, vars=7) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): Computing the final fit ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): > Rbs COBS regression spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 3]: -1.0000030, 0.1836735, 1.0000030 > > plot(x,y) > lines(predict(Rbs), col = 2, lwd = 1.5) > lines(spline(x,f.true), col = "gray40") > > Rbsub <- cobs(x,y,constraint="increase",pointwise=con, n.sub = 45) Performing general knot selection ... N[rq,L1]= (45, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,6, vars=5) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,7, vars=6) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,8, vars=7) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... N[rq,L1]= (45, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): Computing the final fit ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): > summary(Rbsub) COBS regression spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", n.sub = 45, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 3]: -1.0000030, 0.1836735, 1.0000030 with 3 pointwise constraints coef : [1] 0.03313316 0.18255331 0.94257663 1.00000042 R^2 = 94.37% ; empirical tau (over all): 25/50 = 0.5 (target tau : 0.5) > lines(predict(Rbsub), col = 4, lwd = 1) > > ## compute the median smoothing B-spline using automatically chosen lambda > Sbs <- cobs(x,y,constraint="increase",pointwise=con,lambda=-1) Searching for optimal lambda. This may take a while. While you are waiting, here is something you can consider to speed up the process: (a) Use a smaller number of knots; (b) Increase `factor' to 2 or 3; (c) Set lambda==0 to exclude the penalty term. N[rq,L1]= (50, 1); ([ei]qc= 1,60, vars=22) LAMBDA will be updated; (tau,lam)= (0.5, -1): > Sbs COBS smoothing spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", lambda = -1, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 20]: -1.0000010, -0.9183673, -0.7959184, ... , 1.0000010 lambda = 0.4508836, selected via SIC, out of 116 ones. > plot(Sbs$pp.lambda[-1], Sbs$sic[-1], log = "x", + main = "SIC ~ lambda", xlab = expression(lambda), ylab = "SIC") > axis(1, at = Sbs$lambda, label = expression(hat(lambda)), + col.axis = 2, mgp = c(3, 0.5, 0)) > > Sb1 <- cobs(x,y,constraint="increase",pointwise=con,lambda=-1, degree=1) Searching for optimal lambda. This may take a while. While you are waiting, here is something you can consider to speed up the process: (a) Use a smaller number of knots; (b) Increase `factor' to 2 or 3; (c) Set lambda==0 to exclude the penalty term. N[rq,L1]= (50, 18); ([ei]qc= 1,21, vars=20) LAMBDA will be updated; (tau,lam)= (0.5, -1): > summary(Sb1) COBS smoothing spline (degree = 1) from call: cobs(x = x, y = y, constraint = "increase", degree = 1, lambda = -1, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 5 (5) knots[1 .. 20]: -1.0000010, -0.9183673, -0.7959184, ... , 1.0000010 lambda = 0.2580284, selected via SIC, out of 44 ones. with 3 pointwise constraints coef : [1] 0.01355796 0.04681331 0.09669571 0.12995064 0.17983305 0.21308798 [7] 0.28194687 0.35080575 0.39671167 0.46557056 0.51147648 0.59368266 [13] 0.64848678 0.73069296 0.81289914 0.86770326 0.90378419 0.92783814 [19] 0.96391907 1.00000029 R^2 = 94.34% ; empirical tau (over all): 25/50 = 0.5 (target tau : 0.5) > pxx <- predict(Sb1, xx <- seq(-1.2, 1.2, len = 201), interval = "both") > plot(x,y, main = deparse(Sb1$call), + xlim = range(xx), ylim = range(y, pxx[,"fit"])) > lines(pxx, col = 2) > rug(Sb1$knots, col = 4, lwd = 1.6)# (too many knots) > matlines(pxx[,1], pxx[,-(1:2)], col= rep(c("green","blue"),c(2,2)), lty=2) > > > > cleanEx(); ..nameEx <- "cobsOld" > > ### * cobsOld > > flush(stderr()); flush(stdout()) > > ### Name: cobsOld > ### Title: COnstrained B-Splines Nonparametric Regression Quantiles > ### Aliases: cobsOld > ### Keywords: smooth regression > > ### ** Examples > > x <- seq(-1,1,,50) > y <- (f.true <- pnorm(2*x)) + rnorm(50)/10 > ## specify pointwise constraints (boundary conditions) > con <- rbind(c( 1,min(x),0), # f(min(x)) >= 0 + c(-1,max(x),1), # f(max(x)) <= 1 + c(0, 0, 0.5))# f(0) = 0.5 > > ## obtain the median regression B-spline using automatically selected knots > cobsOld(x,y,constraint="increase",pointwise=con)->cobs.o Performing general knot selection ... N[rq,L1]= (50, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,6, vars=5) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,7, vars=6) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,8, vars=7) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): Computing the final fit ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): > > plot(x,y) > lines(cobs.o$z, cobs.o$fit, col = 2, lwd = 1.5) > lines(spline(x,f.true), col = "gray40") > lines(cobs.o$z, cobs.o$cb.lo,lty=2, col = 3) > lines(cobs.o$z, cobs.o$cb.up,lty=2, col = 3) > > ## compute the median smoothing B-spline using automatically chosen lambda > cobsOld(x,y,constraint="increase",pointwise=con,lambda=-1)->cobs.oo Searching for optimal lambda. This may take a while. While you are waiting, here is something you can consider to speed up the process: (a) Use a smaller number of knots; (b) Increase `factor' to 2 or 3; (c) Set lambda==0 to exclude the penalty term. N[rq,L1]= (50, 1); ([ei]qc= 1,60, vars=22) LAMBDA will be updated; (tau,lam)= (0.5, -1): > plot(x,y, main = "COBS Median smoothing spline, automatical lambda") > lines(spline(x,f.true), col = "gray40") > lines(cobs.oo$z, cobs.oo$fit) > lines(cobs.oo$z, cobs.oo$cb.lo,lty=2) > lines(cobs.oo$z, cobs.oo$cb.up,lty=2) > > plot(cobs.oo$pp.lambda[-1], cobs.oo$sic[-1], log = "x", + main = "SIC ~ lambda", xlab = expression(lambda), ylab = "SIC") > axis(1, at = cobs.oo$lambda, label = expression(hat(lambda)), + col.axis = 2, mgp = c(3, 0.5, 0)) > > > > cleanEx(); ..nameEx <- "drqssbc" > > ### * drqssbc > > flush(stderr()); flush(stdout()) > > ### Name: drqssbc > ### Title: Regression Quantile Smoothing Spline with Constraints > ### Aliases: drqssbc > ### Keywords: smooth utilities > > ### ** Examples > > set.seed(1243) > x <- 1:32 > fx <- (x-5)*(x-15)^2*(x-21) > y <- fx + round(rnorm(x,s = 0.25),2) > ## FAILS drqssbc(x,y,nrq=32,lam=1,degree=1,knots=c(1,5,15,32)) > > > > cleanEx(); ..nameEx <- "exHe" > > ### * exHe > > flush(stderr()); flush(stdout()) > > ### Name: exHe > ### Title: Small Dataset Example of He > ### Aliases: exHe > ### Keywords: datasets > > ### ** Examples > > data(exHe) > plot(exHe, main = "He's 10 point example and cobs() fits") > tm <- tapply(exHe$y, exHe$x, mean) > lines(unique(exHe$x), tm, lty = 2) > ## 4 warnings : > cH <- cobsOld(exHe$x, exHe$y, constraint = "increase") WARNING! It looks like you are using 3 knots instead of the default number of 6 knots. Performing general knot selection ... N[rq,L1]= (10, 0); ([ei]qc= 0,2, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (10, 0); ([ei]qc= 0,3, vars=4) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... Computing the final fit ... N[rq,L1]= (10, 0); ([ei]qc= 0,2, vars=3) _fixed_ (tau,lam)= (0.5, 1): > lines(cH$z, cH$fit, col = 3) > cHn <- cobsOld(exHe$x, exHe$y, constraint = "none") WARNING! It looks like you are using 3 knots instead of the default number of 6 knots. Performing general knot selection ... N[rq,L1]= (10, 0); ([ei]qc= 0,0, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (10, 0); ([ei]qc= 0,0, vars=4) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... Computing the final fit ... N[rq,L1]= (10, 0); ([ei]qc= 0,0, vars=3) _fixed_ (tau,lam)= (0.5, 1): > lines(cHn$z, cHn$fit, col = 4) > cHd <- cobsOld(exHe$x, exHe$y, constraint = "decrease") WARNING! It looks like you are using 3 knots instead of the default number of 6 knots. Performing general knot selection ... N[rq,L1]= (10, 0); ([ei]qc= 0,2, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (10, 0); ([ei]qc= 0,3, vars=4) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... Computing the final fit ... N[rq,L1]= (10, 0); ([ei]qc= 0,2, vars=3) _fixed_ (tau,lam)= (0.5, 1): > lines(cHd$z, cHd$fit, col = 5) > > > > cleanEx(); ..nameEx <- "globtemp" > > ### * globtemp > > flush(stderr()); flush(stdout()) > > ### Name: globtemp > ### Title: Annual Average Global Surface Temperature > ### Aliases: globtemp > ### Keywords: datasets > > ### ** Examples > > data(globtemp) > plot(globtemp, main = "Annual Global Temperature Deviations") > > > > > cleanEx(); ..nameEx <- "predict.cobs" > > ### * predict.cobs > > flush(stderr()); flush(stdout()) > > ### Name: predict.cobs > ### Title: Predict method for COBS Fits > ### Aliases: predict.cobs > ### Keywords: regression > > ### ** Examples > > example(cobs) # continuing : cobs> x <- seq(-1, 1, , 50) cobs> y <- (f.true <- pnorm(2 * x)) + rnorm(50)/10 cobs> con <- rbind(c(1, min(x), 0), c(-1, max(x), 1), c(0, 0, 0.5)) cobs> Rbs <- cobs(x, y, constraint = "increase", pointwise = con) Performing general knot selection ... N[rq,L1]= (50, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,6, vars=5) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,7, vars=6) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,8, vars=7) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (50, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): Computing the final fit ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): cobs> Rbs COBS regression spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 3]: -1.0000030, 0.1836735, 1.0000030 cobs> plot(x, y) cobs> lines(predict(Rbs), col = 2, lwd = 1.5) cobs> lines(spline(x, f.true), col = "gray40") cobs> Rbsub <- cobs(x, y, constraint = "increase", pointwise = con, n.sub = 45) Performing general knot selection ... N[rq,L1]= (45, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,6, vars=5) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,7, vars=6) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,8, vars=7) _fixed_ (tau,lam)= (0.5, 1): Deleting unnecessary knots ... N[rq,L1]= (45, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): N[rq,L1]= (45, 0); ([ei]qc= 1,4, vars=3) _fixed_ (tau,lam)= (0.5, 1): Computing the final fit ... N[rq,L1]= (50, 0); ([ei]qc= 1,5, vars=4) _fixed_ (tau,lam)= (0.5, 1): cobs> summary(Rbsub) COBS regression spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", n.sub = 45, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 3]: -1.0000030, 0.1836735, 1.0000030 with 3 pointwise constraints coef : [1] 0.03313316 0.18255331 0.94257663 1.00000042 R^2 = 94.37% ; empirical tau (over all): 25/50 = 0.5 (target tau : 0.5) cobs> lines(predict(Rbsub), col = 4, lwd = 1) cobs> Sbs <- cobs(x, y, constraint = "increase", pointwise = con, lambda = -1) Searching for optimal lambda. This may take a while. While you are waiting, here is something you can consider to speed up the process: (a) Use a smaller number of knots; (b) Increase `factor' to 2 or 3; (c) Set lambda==0 to exclude the penalty term. N[rq,L1]= (50, 1); ([ei]qc= 1,60, vars=22) LAMBDA will be updated; (tau,lam)= (0.5, -1): cobs> Sbs COBS smoothing spline (degree = 2) from call: cobs(x = x, y = y, constraint = "increase", lambda = -1, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 4 (4) knots[1 .. 20]: -1.0000010, -0.9183673, -0.7959184, ... , 1.0000010 lambda = 0.4508836, selected via SIC, out of 113 ones. cobs> plot(Sbs$pp.lambda[-1], Sbs$sic[-1], log = "x", main = "SIC ~ lambda", xlab = expression(lambda), ylab = "SIC") cobs> axis(1, at = Sbs$lambda, label = expression(hat(lambda)), col.axis = 2, mgp = c(3, 0.5, 0)) cobs> Sb1 <- cobs(x, y, constraint = "increase", pointwise = con, lambda = -1, degree = 1) Searching for optimal lambda. This may take a while. While you are waiting, here is something you can consider to speed up the process: (a) Use a smaller number of knots; (b) Increase `factor' to 2 or 3; (c) Set lambda==0 to exclude the penalty term. N[rq,L1]= (50, 18); ([ei]qc= 1,21, vars=20) LAMBDA will be updated; (tau,lam)= (0.5, -1): cobs> summary(Sb1) COBS smoothing spline (degree = 1) from call: cobs(x = x, y = y, constraint = "increase", degree = 1, lambda = -1, pointwise = con) {tau=0.5}-quantile; dimensionality of fit: 5 (5) knots[1 .. 20]: -1.0000010, -0.9183673, -0.7959184, ... , 1.0000010 lambda = 0.2580284, selected via SIC, out of 44 ones. with 3 pointwise constraints coef : [1] 0.01355796 0.04681331 0.09669571 0.12995064 0.17983305 0.21308798 [7] 0.28194687 0.35080575 0.39671167 0.46557056 0.51147648 0.59368266 [13] 0.64848678 0.73069296 0.81289914 0.86770326 0.90378419 0.92783814 [19] 0.96391907 1.00000029 R^2 = 94.34% ; empirical tau (over all): 26/50 = 0.52 (target tau : 0.5) cobs> pxx <- predict(Sb1, xx <- seq(-1.2, 1.2, len = 201), interval = "both") cobs> plot(x, y, main = deparse(Sb1$call), xlim = range(xx), ylim = range(y, pxx[, "fit"])) cobs> lines(pxx, col = 2) cobs> rug(Sb1$knots, col = 4, lwd = 1.6) cobs> matlines(pxx[, 1], pxx[, -(1:2)], col = rep(c("green", "blue"), c(2, 2)), lty = 2) > (pRbs <- predict(Rbs)) z fit [1,] -1.00000288 0.03313319 [2,] -0.97980080 0.03832106 [3,] -0.95959872 0.04368394 [4,] -0.93939665 0.04922182 [5,] -0.91919457 0.05493470 [6,] -0.89899249 0.06082259 [7,] -0.87879041 0.06688547 [8,] -0.85858833 0.07312335 [9,] -0.83838625 0.07953624 [10,] -0.81818418 0.08612412 [11,] -0.79798210 0.09288701 [12,] -0.77778002 0.09982490 [13,] -0.75757794 0.10693779 [14,] -0.73737586 0.11422568 [15,] -0.71717378 0.12168857 [16,] -0.69697170 0.12932646 [17,] -0.67676963 0.13713935 [18,] -0.65656755 0.14512725 [19,] -0.63636547 0.15329014 [20,] -0.61616339 0.16162804 [21,] -0.59596131 0.17014094 [22,] -0.57575923 0.17882884 [23,] -0.55555716 0.18769173 [24,] -0.53535508 0.19672963 [25,] -0.51515300 0.20594254 [26,] -0.49495092 0.21533044 [27,] -0.47474884 0.22489334 [28,] -0.45454676 0.23463125 [29,] -0.43434469 0.24454415 [30,] -0.41414261 0.25463206 [31,] -0.39394053 0.26489497 [32,] -0.37373845 0.27533288 [33,] -0.35353637 0.28594579 [34,] -0.33333429 0.29673370 [35,] -0.31313222 0.30769661 [36,] -0.29293014 0.31883452 [37,] -0.27272806 0.33014744 [38,] -0.25252598 0.34163535 [39,] -0.23232390 0.35329827 [40,] -0.21212182 0.36513618 [41,] -0.19191974 0.37714910 [42,] -0.17171767 0.38933702 [43,] -0.15151559 0.40169994 [44,] -0.13131351 0.41423786 [45,] -0.11111143 0.42695078 [46,] -0.09090935 0.43983871 [47,] -0.07070727 0.45290163 [48,] -0.05050520 0.46613956 [49,] -0.03030312 0.47955248 [50,] -0.01010104 0.49314041 [51,] 0.01010104 0.50690334 [52,] 0.03030312 0.52084127 [53,] 0.05050520 0.53495420 [54,] 0.07070727 0.54924213 [55,] 0.09090935 0.56370506 [56,] 0.11111143 0.57834300 [57,] 0.13131351 0.59315593 [58,] 0.15151559 0.60814387 [59,] 0.17171767 0.62330680 [60,] 0.19191974 0.63860437 [61,] 0.21212182 0.65367716 [62,] 0.23232390 0.66844033 [63,] 0.25252598 0.68289385 [64,] 0.27272806 0.69703774 [65,] 0.29293014 0.71087200 [66,] 0.31313222 0.72439661 [67,] 0.33333429 0.73761160 [68,] 0.35353637 0.75051694 [69,] 0.37373845 0.76311265 [70,] 0.39394053 0.77539873 [71,] 0.41414261 0.78737517 [72,] 0.43434469 0.79904197 [73,] 0.45454676 0.81039914 [74,] 0.47474884 0.82144668 [75,] 0.49495092 0.83218457 [76,] 0.51515300 0.84261284 [77,] 0.53535508 0.85273146 [78,] 0.55555716 0.86254045 [79,] 0.57575923 0.87203981 [80,] 0.59596131 0.88122952 [81,] 0.61616339 0.89010961 [82,] 0.63636547 0.89868005 [83,] 0.65656755 0.90694087 [84,] 0.67676963 0.91489204 [85,] 0.69697170 0.92253358 [86,] 0.71717378 0.92986549 [87,] 0.73737586 0.93688776 [88,] 0.75757794 0.94360039 [89,] 0.77778002 0.95000339 [90,] 0.79798210 0.95609675 [91,] 0.81818418 0.96188048 [92,] 0.83838625 0.96735457 [93,] 0.85858833 0.97251902 [94,] 0.87879041 0.97737384 [95,] 0.89899249 0.98191902 [96,] 0.91919457 0.98615457 [97,] 0.93939665 0.99008048 [98,] 0.95959872 0.99369676 [99,] 0.97980080 0.99700340 [100,] 1.00000288 1.00000041 > str(pSbs <- predict(Sbs, xx, interval = "both")) num [1:201, 1:6] -1.20 -1.19 -1.18 -1.16 -1.15 ... - attr(*, "dimnames")=List of 2 ..$ : NULL ..$ : chr [1:6] "z" "fit" "cb.lo" "cb.up" ... > > plot(x,y, xlim = range(xx), ylim = range(y, pSbs[,2], finite = TRUE), + main = "COBS Median smoothing spline, automatical lambda") > lines(pSbs, col = "red") > lines(spline(x,f.true), col = "gray40") > matlines(pSbs[,1], pSbs[,-(1:2)], + col= rep(c("green","blue"),c(2,2)), lty=2) > > > > cleanEx(); ..nameEx <- "qbsks" > > ### * qbsks > > flush(stderr()); flush(stdout()) > > ### Name: qbsks > ### Title: Quantile B-Spline with Fixed Knots > ### Aliases: qbsks > ### Keywords: smooth utilities > > ### ** Examples > > > > > ### *