diff --git a/rpanel/hist_button.R b/rpanel/hist_button.R
new file mode 100644
index 0000000000000000000000000000000000000000..d7529cc182c96d5dfef73e872e052d441c238793
--- /dev/null
+++ b/rpanel/hist_button.R
@@ -0,0 +1,20 @@
+## Botão de ação (rp.button)
+
+require(rpanel)
+
+x <- precip
+ht <- hist(x)
+
+hist.reactive <- function(input){
+ col <- sample(colors(), size=1)
+ plot(ht, main=NULL,
+ ylab="Frequência absoluta", xlab="Precipitação",
+ col=col, sub=col)
+ return(input)
+}
+
+panel <- rp.control(title="Histograma")
+rp.button(panel=panel,
+ title="Nova cor!",
+ action=hist.reactive)
+
diff --git a/rpanel/hist_button.gif b/rpanel/hist_button.gif
new file mode 100644
index 0000000000000000000000000000000000000000..e4036099cb3b36cf691ec983b24ba6e45cf5aa28
Binary files /dev/null and b/rpanel/hist_button.gif differ
diff --git a/rpanel/hist_checkbox.R b/rpanel/hist_checkbox.R
new file mode 100644
index 0000000000000000000000000000000000000000..8c7cd54214b99ad7af376505cf20a00feeb3984b
--- /dev/null
+++ b/rpanel/hist_checkbox.R
@@ -0,0 +1,30 @@
+## Caixa de seleção (rp.checkbox)
+
+require(rpanel)
+
+x <- precip
+ht <- hist(x)
+col <- rep("#3366CC", length(ht$counts))
+
+hist.reactive <- function(input){
+ if(input$modal){
+ col[which.max(ht$counts)] <- "#142952"
+ }
+ plot(ht, col=col, main=NULL,
+ ylab="Frequência absoluta",
+ xlab="Precipitação")
+ if(input$rg){
+ rug(x)
+ }
+ return(input)
+}
+
+panel <- rp.control(title="Histograma")
+rp.checkbox(panel=panel, variable=rg,
+ title="Marcar sobre eixo com os valores?",
+ initval=FALSE,
+ action=hist.reactive)
+rp.checkbox(panel=panel, variable=modal,
+ title="Destacal a classe modal?",
+ initval=FALSE,
+ action=hist.reactive)
diff --git a/rpanel/hist_checkbox.gif b/rpanel/hist_checkbox.gif
new file mode 100644
index 0000000000000000000000000000000000000000..ba98b98cd8781fe9d3480472e8c01d628052851f
Binary files /dev/null and b/rpanel/hist_checkbox.gif differ
diff --git a/rpanel/hist_checkboxgroup.R b/rpanel/hist_checkboxgroup.R
new file mode 100644
index 0000000000000000000000000000000000000000..f2c7776a4f63b4971a75dd851275f0c436fd81bd
--- /dev/null
+++ b/rpanel/hist_checkboxgroup.R
@@ -0,0 +1,28 @@
+## Caixa de seleção múltipla (rp.checkbox)
+
+require(rpanel)
+
+x <- precip
+ht <- hist(x)
+nc <- length(ht$counts)
+
+cols <- c(Vermelho="#F81D54", Amarelo="#FF9F1E",
+ Azul="#2791E1", Verde="#72F51D")
+cols2 <- c(cols, rev(cols))
+
+hist.reactive <- function(input){
+ seqcol <- colorRampPalette(cols2[input$colors])
+ plot(ht, col=seqcol(nc),
+ main=NULL,
+ ylab="Frequência absoluta",
+ xlab="Precipitação")
+ return(input)
+}
+
+panel <- rp.control(title="Histograma")
+rp.checkbox(panel=panel, variable=colors,
+ title="Escolha as cores para interpolar:",
+ labels=names(cols2),
+ initval=c(TRUE, is.na(cols2)[-1]),
+ action=hist.reactive)
+
diff --git a/rpanel/hist_checkboxgroup.gif b/rpanel/hist_checkboxgroup.gif
new file mode 100644
index 0000000000000000000000000000000000000000..b393512494530ede568b632368c8942d53d8a6af
Binary files /dev/null and b/rpanel/hist_checkboxgroup.gif differ
diff --git a/rpanel/hist_numeric.R b/rpanel/hist_numeric.R
new file mode 100644
index 0000000000000000000000000000000000000000..5a72707bddf0f41d3dc2b01cee95e7342bbaf921
--- /dev/null
+++ b/rpanel/hist_numeric.R
@@ -0,0 +1,35 @@
+## Entrada numérica (rp.numeric)
+
+require(rpanel)
+
+x <- precip
+ht <- hist(x)
+
+hist.reactive <- function(input){
+ m <- input$mar
+ par(mar=c(m, m, 1, 1))
+ plot(ht, col="#660066",
+ main=NULL, axes=FALSE, ann=FALSE,
+ xaxt="n", yaxt="n")
+ box(bty="L")
+ axis(side=1, cex.axis=input$cexaxis)
+ axis(side=2, cex.axis=input$cexaxis)
+ title(ylab="Frequência absoluta",
+ xlab="Precipitação",
+ line=input$line)
+ return(input)
+}
+
+panel <- rp.control(title="Histograma")
+rp.doublebutton(panel=panel, variable=mar,
+ title="Tamanho das margens:",
+ initval=5, range=c(3, 7), step=0.5,
+ action=hist.reactive)
+rp.doublebutton(panel=panel, variable=cexaxis,
+ title="Tamanho do texto dos eixos:",
+ initval=1, range=c(0.5, 2), step=0.1,
+ action=hist.reactive)
+rp.doublebutton(panel=panel, variable=line,
+ title="Distância dos rótulos dos eixos:",
+ initval=3, range=c(1, 4), step=0.1,
+ action=hist.reactive)
diff --git a/rpanel/hist_numeric.gif b/rpanel/hist_numeric.gif
new file mode 100644
index 0000000000000000000000000000000000000000..ea558f25350d526873642f958bbcbc56d8774c67
Binary files /dev/null and b/rpanel/hist_numeric.gif differ
diff --git a/rpanel/hist_radio.R b/rpanel/hist_radio.R
new file mode 100644
index 0000000000000000000000000000000000000000..ad01a71bb57ca269b4aedda3981320f3c13a6e34
--- /dev/null
+++ b/rpanel/hist_radio.R
@@ -0,0 +1,28 @@
+## Múltipla escolha (rp.radiogroup)
+
+require(rpanel)
+
+x <- precip
+ht <- hist(x)
+
+hist.reactive <- function(input){
+ plot(ht,
+ col=input$col,
+ main=NULL,
+ ylab="Frequência absoluta",
+ xlab="Precipitação")
+ return(input)
+}
+
+choices <- c(Turquesa="#00CC99",
+ Azul="#0066FF",
+ Rosa="#FF3399",
+ Laranja="#FF6600",
+ Roxo="#660066",
+ "Verde limão"="#99FF33")
+
+panel <- rp.control(title="Histograma")
+rp.radiogroup(panel=panel, variable=col,
+ title="Escolha a cor para as barras:",
+ vals=choices, labels=names(choices),
+ action=hist.reactive)
diff --git a/rpanel/hist_radio.gif b/rpanel/hist_radio.gif
new file mode 100644
index 0000000000000000000000000000000000000000..93467479a0d6817348d577f6d4dd0126b19f40f7
Binary files /dev/null and b/rpanel/hist_radio.gif differ
diff --git a/rpanel/hist_select.R b/rpanel/hist_select.R
new file mode 100644
index 0000000000000000000000000000000000000000..adb198662287ca01ca65389692deada2ac646cec
--- /dev/null
+++ b/rpanel/hist_select.R
@@ -0,0 +1,44 @@
+## Caixa de seleção (rp.listbox)
+
+require(rpanel)
+
+nclass <- c("Sturges", "Scott", "Freedman-Diaconis")
+obj <- c("precip","rivers","islands")
+
+hist.reactive <- function(input){
+ L <- switch(input$obj,
+ precip=list(x=precip,
+ xlab="Precipitação anual média (polegadas)"),
+ rivers=list(x=rivers,
+ xlab="Comprimento dos rios (milhas)"),
+ islands=list(x=islands,
+ xlab="Ãrea de ilhas (1000 milhas quadradas)"))
+ hist(L$x,
+ breaks=input$nclass,
+ col="#8F0047",
+ main=NULL,
+ ylab="Frequência absoluta",
+ xlab=L$xlab)
+ rug(L$x)
+ return(input)
+}
+
+panel <- rp.control(title="Histograma")
+rp.combo(panel=panel, variable=obj,
+ prompt="Escolha o conjunto de dados:",
+ vals=obj, initval=obj[1],
+ action=hist.reactive)
+rp.combo(panel=panel, variable=nclass,
+ prompt="Escolha a regra para número de classes:",
+ vals=nclass, initval=nclass[1],
+ action=hist.reactive)
+
+panel <- rp.control(title="Histograma")
+rp.listbox(panel=panel, variable=obj,
+ title="Escolha o conjunto de dados:",
+ vals=obj, initval=obj[1],
+ action=hist.reactive)
+rp.listbox(panel=panel, variable=nclass,
+ title="Escolha a regra para número de classes:",
+ vals=nclass, initval=nclass[1],
+ action=hist.reactive)
diff --git a/rpanel/hist_select.gif b/rpanel/hist_select.gif
new file mode 100644
index 0000000000000000000000000000000000000000..5a4f841af273a1612252d00df5431c129f8865f1
Binary files /dev/null and b/rpanel/hist_select.gif differ
diff --git a/rpanel/hist_select2.R b/rpanel/hist_select2.R
new file mode 100644
index 0000000000000000000000000000000000000000..8c7d63b02ac24b11aac6687db12698f185d37fee
--- /dev/null
+++ b/rpanel/hist_select2.R
@@ -0,0 +1,32 @@
+## Caixas de seleção (rp.listbox e rp.radiogroup)
+
+require(rpanel)
+
+fml <- names(X11Fonts())
+fnt <- c("plain"=1, "bold"=2, "italic"=3, "bold-italic"=4)
+
+x <- precip
+ht <- hist(x)
+
+hist.reactive <- function(input){
+ f <- as.integer(input$fnt)
+ plot(ht,
+ family=input$fml,
+ font=as.integer(input$fnt),
+ col="#FF9200",
+ main=NULL,
+ ylab="Frequência absoluta",
+ xlab="Precipitação")
+ return(input)
+}
+
+panel <- rp.control(title="Histograma")
+rp.listbox(panel=panel, variable=fml,
+ title="Escolha o tipo de fonte:",
+ vals=fml, initval=fml[1],
+ action=hist.reactive)
+rp.radiogroup(panel=panel, variable=fnt,
+ title="Escolha o estilo de fonte:",
+ vals=fnt, initval=fnt[1],
+ labels=names(fnt),
+ action=hist.reactive)
diff --git a/rpanel/hist_select2.gif b/rpanel/hist_select2.gif
new file mode 100644
index 0000000000000000000000000000000000000000..9cfce5b4a60b5ed9437a399a6f5d96e1e1a9a777
Binary files /dev/null and b/rpanel/hist_select2.gif differ
diff --git a/rpanel/hist_slider.R b/rpanel/hist_slider.R
new file mode 100644
index 0000000000000000000000000000000000000000..2cf6ccf10e5a5d12f8dcde5deaff092a77d4bb71
--- /dev/null
+++ b/rpanel/hist_slider.R
@@ -0,0 +1,25 @@
+## Deslizador (rp.slider)
+
+require(rpanel)
+
+x <- precip
+
+## Extremos com amplitude estendida em 5%.
+a <- extendrange(x, f=0.05)
+
+hist.reactive <- function(input){
+ bks <- seq(a[1], a[2], length.out=input$nclass+1)
+ hist(x,
+ breaks=bks,
+ main=NULL,
+ col="#008A8A",
+ ylab="Frequência absoluta",
+ xlab="Precipitação")
+ return(input)
+}
+
+panel <- rp.control(title="Histograma")
+rp.slider(panel=panel, variable=nclass,
+ title="Escolha o número de classes:",
+ from=1, to=30, resolution=1, initval=10,
+ action=hist.reactive)
diff --git a/rpanel/hist_slider.gif b/rpanel/hist_slider.gif
new file mode 100644
index 0000000000000000000000000000000000000000..97b5081dd47cdf8a059172120abd4c69ae65bc45
Binary files /dev/null and b/rpanel/hist_slider.gif differ
diff --git a/rpanel/hist_text.R b/rpanel/hist_text.R
new file mode 100644
index 0000000000000000000000000000000000000000..7e72561048e9aa7a2c2837655e9f78261b7060c1
--- /dev/null
+++ b/rpanel/hist_text.R
@@ -0,0 +1,26 @@
+## Entrada de texto (rp.textentry)
+
+require(rpanel)
+
+x <- precip
+ht <- hist(x)
+
+hist.reactive <- function(input){
+ plot(ht, col="#006666",
+ ylab="Frequência absoluta",
+ xlab="Precipitação",
+ main=input$main,
+ sub=input$sub)
+ return(input)
+}
+
+panel <- rp.control(title="Histograma")
+rp.textentry(panel=panel, variable=main,
+ labels="Texto para o tÃtulo:",
+ initval="",
+ action=hist.reactive)
+rp.textentry(panel=panel, variable=sub,
+ labels="Texto para o subtÃtulo:",
+ initval="",
+ action=hist.reactive)
+
diff --git a/rpanel/hist_text.gif b/rpanel/hist_text.gif
new file mode 100644
index 0000000000000000000000000000000000000000..7a0f85fff9fe9e99a9f90f2a335cddca87f2e95c
Binary files /dev/null and b/rpanel/hist_text.gif differ
diff --git a/rpanel/lm.R b/rpanel/lm.R
new file mode 100644
index 0000000000000000000000000000000000000000..60f3cf404209b806c1e153e44caeb18e948634db
--- /dev/null
+++ b/rpanel/lm.R
@@ -0,0 +1,99 @@
+require(rpanel)
+require(plyr)
+require(latticeExtra)
+## require(car)
+require(alr3)
+require(wzRfun)
+
+##----------------------------------------------------------------------
+
+url <- "http://westfall.ba.ttu.edu/isqs5349/Rdata/turtles.txt"
+tur <- read.table(url, header=TRUE, sep="\t")
+str(tur)
+
+xtabs(~gender, tur)
+
+tur$Gender <- factor(tur$gender)
+
+form <- c("null"=.~1,
+ "gend"=.~Gender,
+ "leng"=.~length,
+ "gend+leng"=.~Gender+length,
+ "gend*leng"=.~Gender*length)
+
+m1 <- lm(height~Gender*length, data=tur)
+
+p1 <- xyplot(height~length, groups=Gender, data=tur,
+ ylab="Altura", xlab="Comprimento")
+p1
+
+pred <- ddply(tur, .(Gender), summarise,
+ length=seq(min(tur$length), max(tur$length), l=20))
+
+draw.model <- function(panel){
+ m3 <- update(m1, as.formula(form[[panel$f]]))
+ pred <- cbind(pred, predict(m3, newdata=pred, interval="confidence"))
+ p2 <- xyplot(fit~length, groups=Gender, data=pred, type="l",
+ ly=pred$lwr, uy=pred$upr, cty="bands", alpha=0.5,
+ prepanel=prepanel.cbH, panel=panel.superpose,
+ panel.groups=panel.cbH)
+ print(p1+as.layer(p2, under=TRUE))
+ panel
+}
+
+panel <- rp.control()
+rp.listbox(panel, variable=f, vals=names(form), action=draw.model)
+
+##----------------------------------------------------------------------
+
+str(sleep1)
+## help(sleep1, h="html")
+
+## Danger Index.
+sleep1$D <- factor(sleep1$D)
+sleep1$lbw <- log(sleep1$BodyWt)
+
+xyplot(TS~lbw|D, data=sleep1)
+xyplot(TS~lbw|D, data=sleep1, type=c("p","r"))
+
+##-----------------------------------------------------------------------------
+## Ajuste do modelo.
+
+m1 <- lm(log(TS)~D*lbw, data=sleep1)
+
+## Diagnóstico.
+par(mfrow=c(2,2)); plot(m1); layout(1)
+
+##-----------------------------------------------------------------------------
+## Estimativas por grupo.
+
+## R² com SQT corrigida para média.
+summary(m1)
+
+form <- c("null"=.~1,
+ "D"=.~D,
+ "lbw"=.~lbw,
+ "D+lbw"=.~D+lbw,
+ "D*lbw"=.~D*lbw)
+
+pred <- ddply(sleep1, .(D), summarise,
+ lbw=seq(min(lbw), max(lbw), l=20)
+ ## lbw=extendseq(lbw, l=12)
+ )
+
+draw.model <- function(panel){
+ m3 <- update(m1, as.formula(form[[panel$f]]))
+ pred <- cbind(pred,
+ predict(m3, newdata=pred, interval="confidence"))
+ p1 <- xyplot(log(TS)~lbw|D, data=sleep1)
+ p2 <- xyplot(fit~lbw|D, data=pred, type="l",
+ ly=pred$lwr, uy=pred$upr, cty="bands", alpha=0.25,
+ prepanel=prepanel.cbH,
+ panel=panel.cbH)
+ print(p1+as.layer(p2, under=TRUE))
+ panel
+}
+
+panel <- rp.control()
+rp.listbox(panel, variable=f, vals=names(form), action=draw.model,
+ title="Select a model")
diff --git a/rpanel/rp.density.R b/rpanel/rp.density.R
new file mode 100644
index 0000000000000000000000000000000000000000..b664de8be0f1e8866700a7059977a677a0d18884
--- /dev/null
+++ b/rpanel/rp.density.R
@@ -0,0 +1,67 @@
+require(rpanel)
+
+rp.density <- function(x){
+ draw.density <- function(panel){
+ ## Plot density curve and add the rugs.
+ n <- length(na.omit(panel$x))
+ aux <- density(panel$x, width=panel$width, kernel=panel$kernel)
+ plot(aux, main=NA)
+ rug(panel$x)
+ ## The interval.
+ arrows(panel$xc-0.5*panel$width, 0,
+ panel$xc+0.5*panel$width, 0,
+ length=0.1, code=3, angle=90, col=2)
+ ## Line at center.
+ yc <- approx(aux$x, aux$y, xout=panel$xc)
+ arrows(panel$xc, 0, panel$xc, yc$y, length=0.1, col=2)
+ ## Density for a single point.
+ d <- density(panel$xc, width=panel$width,
+ kernel=panel$kernel, n=128)
+ lines(d$x, d$y/n, col=2)
+ switch(panel$points,
+ "points"={
+ i <- findInterval(panel$x, vec=d$x[c(1,128)])
+ ## i <- findInterval(panel$x, vec=range(d$x))
+ xin <- panel$x[i==1L]
+ yin <- approx(d$x, d$y/n, xout=xin)
+ points(yin$x, yin$y, pch=19, cex=0.5, col=2)
+ },
+ "segments"={
+ i <- findInterval(panel$x, vec=d$x[c(1,128)])
+ xin <- panel$x[i==1L]
+ yin <- approx(d$x, d$y/n, xout=xin)
+ segments(yin$x, 0, yin$x, yin$y, col=2)
+ },
+ "none"={
+ NULL
+ }
+ )
+ panel
+ }
+ ## Initialization.
+ init <- density(x)
+ init$er <- extendrange(x, f=0.1)
+ init$kernels <- eval(formals(density.default)$kernel)
+ init$w <- 4*init$bw ## For gaussian.
+ ## Panels.
+ panel <- rp.control(x=x)
+ rp.slider(panel, variable=width,
+ from=init$w*0.05, to=init$w*3, initval=init$w,
+ showvalue=TRUE, action=draw.density, title="Width")
+ rp.slider(panel, variable=xc,
+ from=init$er[1], to=init$er[2],
+ initval=mean(x, na.rm=TRUE),
+ showvalue=TRUE, action=draw.density, title="Center")
+ rp.radiogroup(panel, variable=kernel, vals=init$kernels,
+ action=draw.density, title="Kernel function")
+ rp.radiogroup(panel, variable=points,
+ vals=c("points","segments","none"),
+ action=draw.density, title="Density on points")
+ rp.do(panel, action=draw.density)
+}
+
+x <- rnorm(30)
+
+x11()
+rp.density(precip)
+
diff --git a/rpanel/rp.semiparreg.R b/rpanel/rp.semiparreg.R
new file mode 100644
index 0000000000000000000000000000000000000000..92e28cd63daaf206a8694e32567f63c5a1b16c34
--- /dev/null
+++ b/rpanel/rp.semiparreg.R
@@ -0,0 +1,455 @@
+##----------------------------------------------------------------------
+## Definições da sessão.
+
+require(rpanel)
+require(splines)
+require(mgcv)
+
+##----------------------------------------------------------------------
+
+sui <- read.table("http://www.leg.ufpr.br/~walmes/data/preabate.txt",
+ header=TRUE, dec=",")
+sui$trat <- factor(-1*sui$trat)
+levels(sui$trat) <- c("Sem aspersão","Com aspersão")
+str(sui)
+
+## xyplot(temp~hora, groups=trat, data=sui,
+## xlab="Minutos à partir das 07:00 da manhã",
+## ylab="Temperatura corporal (ºC)",
+## auto.key=TRUE)+
+## glayer(panel.smoother(..., span=0.4))
+
+## xyplot(temp~hora|trat, groups=asp, data=sui, type=c("p","smooth"))
+
+str(sui)
+ca <- subset(sui, trat=="Com aspersão")
+
+##----------------------------------------------------------------------
+## (Global) Polynomial regression.
+
+rp.poly <- function(y, x, er=0){
+ annotations <- function(m0){
+ mtext(side=3, adj=0, line=1.5,
+ text=sprintf("X rank: %i", m0$rank))
+ mtext(side=3, adj=0, line=0.5,
+ text=sprintf("Residual DF: %i", m0$df.residual))
+ press <- sum((residuals(m0)/(1-hatvalues(m0)))^2)
+ mtext(side=3, adj=1, line=0.5,
+ text=sprintf("PRESS: %0.2f", press))
+ sm0 <- summary(m0)
+ lastcoef <- sprintf("High term p-value: %0.5f",
+ sm0$coeff[length(coef(m0)), 4])
+ mtext(side=3, adj=1, line=2.5, text=lastcoef)
+ mtext(side=3, adj=1, line=1.5,
+ text=sprintf("R² (adj. R²): %0.2f (%0.2f)",
+ 100*sm0$r.squared, 100*sm0$adj.r.squared))
+ }
+ draw.poly <- function(panel){
+ with(panel,
+ {
+ m0 <- lm(y~poly(x, degree=degree))
+ xr <- extendrange(x, f=er)
+ yr <- extendrange(y, f=er)
+ xx <- seq(xr[1], xr[2], length.out=200)
+ cb <- predict(m0, newdata=data.frame(x=xx),
+ interval="confidence")
+ plot(y~x, xlim=xr, ylim=yr)
+ matlines(xx, cb, lty=c(1,2,2), col=1)
+ annotations(m0)
+ })
+ panel
+ }
+ maxd <- length(unique(x))-1
+ panel <- rp.control(x=x, y=y, er=er)
+ rp.doublebutton(panel, variable=degree,
+ action=draw.poly, showvalue=TRUE,
+ step=1, initval=1,
+ range=c(1, maxd),
+ title="Degree of polynomial")
+ rp.do(panel, action=draw.poly)
+}
+
+rp.poly(x=rock$area, y=rock$peri, er=0.3)
+rp.poly(x=cars$speed, y=cars$dist, er=0.3)
+rp.poly(x=faithful$eruptions, y=faithful$waiting, er=0.3)
+rp.poly(x=ca$hora, y=ca$temp, er=0.3)
+
+##-----------------------------------------------------------------------------
+## Local polynomial regression.
+
+rp.loess <- function(y, x, n=20, er=0, ...){
+ annotations <- function(m0, y){
+ mtext(side=3, adj=0, line=1.5,
+ text=sprintf("H trace: %0.2f", m0$trace.hat))
+ mtext(side=3, adj=0, line=0.5,
+ text=sprintf("Equivalent num. param.: %0.2f", m0$enp))
+ r2 <- 100*cor(fitted(m0), y)^2
+ mtext(side=3, adj=1, line=1.5,
+ text=sprintf("R²: %0.2f", r2))
+ }
+ draw.loess <- function(panel){
+ with(panel,
+ {
+ m0 <- loess(y~x, span=span, degree=degree,
+ family="gaussian")
+ erx <- extendrange(x, f=er)
+ xx <- seq(erx[1], erx[2], length.out=200)
+ yy <- predict(m0, newdata=data.frame(x=xx))
+ a <- abs(x-x0)
+ if(span < 1){
+ q <- as.integer(span*length(x))
+ d <- sort(a)[q]
+ } else {
+ q <- length(x)
+ d <- max(abs(x-x0))*sqrt(span)
+ }
+ w <- rep(0, length(x))
+ s <- a <= d
+ w[s] <- (1-(a[s]/d)^3)^3
+ i <- as.integer(s)
+ xl <- range(x)
+ xl[1] <- ifelse(x0 - d>xl[1], x0-d, xl[1])
+ xl[2] <- ifelse(x0 + d<xl[2], x0+d, xl[2])
+ f0 <- function(...){
+ y.pred <- sum(w*y)/sum(w)
+ segments(xl[1], y.pred, xl[2], y.pred, col=2,
+ ...)
+ }
+ f1 <- function(...){
+ m <- lm(y ~ poly(x, degree=1), weights=w)
+ y.pred <- predict(m, newdata=list(x=xl))
+ segments(xl[1], y.pred[1], xl[2], y.pred[2], col=2,
+ ...)
+ }
+ f2 <- function(...){
+ m <- lm(y~poly(x, degree=as.integer(degree)),
+ weights=w)
+ x.pred <- seq(xl[1], xl[2], length.out=n)
+ y.pred <- predict(m, newdata=list(x=x.pred))
+ lines(x.pred, y.pred, col=2, ...)
+ }
+ plot(x, y, pch=2*(!s)+1, cex=i*3*w+1,
+ xlim=extendrange(x, f=er),
+ ylim=extendrange(y, f=er), ...)
+ lines(yy~xx)
+ abline(v=c(x0, xl), lty=c(2,3,3))
+ annotations(m0, y)
+ mtext(side=3, adj=1, line=0.5,
+ text=sprintf("Number of obs. used/total: %i/%i",
+ sum(w>0), length(y)))
+ switch(findInterval(degree, c(-Inf, 1, 2, Inf)),
+ "0"=f0(), "1"=f1(), "2"=f2())
+ })
+ panel
+ }
+ xr <- extendrange(x, f=er)
+ panel <- rp.control(x=x, y=y, n=n, er=er)
+ rp.doublebutton(panel, variable=degree, action=draw.loess,
+ showvalue=TRUE, step=1, initval=0, range=c(0, 2),
+ title="Degree")
+ rp.slider(panel, variable=x0, action=draw.loess,
+ from=xr[1], to=xr[2], initval=median(range(x)),
+ showvalue=TRUE, title="x-value")
+ rp.slider(panel, variable=span, action=draw.loess,
+ from=0, to=1.5, initval=0.75,
+ showvalue=TRUE, title="span")
+ rp.do(panel, action=draw.loess)
+}
+
+rp.loess(x=cars$speed, y=cars$dist, n=25, er=0.2)
+rp.loess(x=rock$area, y=rock$peri, n=25)
+rp.loess(x=faithful$eruptions, y=faithful$waiting, er=0.3)
+rp.loess(x=ca$hora, y=ca$temp, er=0.3)
+
+##-----------------------------------------------------------------------------
+## Regression splines.
+
+rp.spline <- function(x, y, er=0, ...){
+ plotfit <- function(x, y, m0, er){
+ xr <- extendrange(x, f=er)
+ xx <- seq(xr[1], xr[2], length.out=200)
+ cb <- predict(m0, newdata=data.frame(x=xx),
+ interval="confidence")
+ matlines(xx, cb, lty=c(1,2,2), col=1)
+ }
+ annotations <- function(m0){
+ mtext(side=3, adj=0, line=1.5,
+ text=sprintf("X rank: %i", m0$rank))
+ mtext(side=3, adj=0, line=0.5,
+ text=sprintf("Residual DF: %i", m0$df.residual))
+ press <- sum((residuals(m0)/(1-hatvalues(m0)))^2)
+ mtext(side=3, adj=1, line=0.5,
+ text=sprintf("PRESS: %0.2f", press))
+ sm0 <- summary(m0)
+ mtext(side=3, adj=1, line=1.5,
+ text=sprintf("R² (adj. R²): %0.2f (%0.2f)",
+ 100*sm0$r.squared, 100*sm0$adj.r.squared))
+ }
+ draw.df.degree <- function(panel){
+ with(panel, {
+ xr <- extendrange(x, f=er)
+ yr <- extendrange(y, f=er)
+ plot(y~x, xlim=xr, ylim=yr, ...)
+ nk <- as.integer(df1)-as.integer(degree1)
+ nk <- max(c(0, nk))
+ if(quantile){
+ s <- seq(0, 1, l=2+nk)
+ q <- quantile(x, probs=s[-c(1,nk+2)])
+ abline(v=q, col="gray50", lty=2)
+ }
+ switch(base1,
+ "bs"={
+ m0 <- lm(y~bs(x, degree=as.integer(degree1),
+ df=as.integer(df1)))
+ plotfit(x, y, m0, er)
+ annotations(m0)
+ mtext(side=3, line=0.5,
+ text=sprintf("Number of internal nodes: %i", nk))
+ },
+ "ns"={
+ m0 <- lm(y~ns(x, df=as.integer(df1)))
+ plotfit(x, y, m0, er)
+ annotations(m0)
+ mtext(side=3, line=0.5,
+ text=sprintf("Number of internal nodes: %i", nk))
+ })
+ })
+ panel
+ }
+ draw.degree.knots <- function(panel){
+ with(panel, {
+ xr <- extendrange(x, f=er)
+ yr <- extendrange(y, f=er)
+ plot(y~x, xlim=xr, ylim=yr, ...)
+ abline(v=k, col=2, lty=2)
+ switch(base2,
+ "bs"={
+ m0 <- lm(y~bs(x, degree=as.integer(degree2),
+ knots=k))
+ plotfit(x, y, m0, er)
+ annotations(m0)
+ },
+ "ns"={
+ if(is.na(k)){
+ m0 <- lm(y~ns(x, df=as.integer(degree2)))
+ } else {
+ m0 <- lm(y~ns(x, knots=k))
+ }
+ plotfit(x, y, m0, er)
+ annotations(m0)
+ })
+ })
+ panel
+ }
+ reset.draw <- function(panel){
+ k <- locator(type="p", pch=19, col=2)$x
+ if(length(k)==0){
+ panel$k <- NA
+ } else {
+ panel$k <- k
+ abline(v=k, col=2, lty=2)
+ }
+ panel$degree2 <- 1
+ rp.do(panel, draw.degree.knots)
+ panel
+ }
+ panel <- rp.control(x=x, y=y, k=NA, ...)
+ rp.notebook(panel, tabnames=c("DF", "Knots"),
+ tabs=c("Control DF", "Choose knots"),
+ width=280, height=200)
+ rp.radiogroup(panel, variable=base1, vals=c("bs","ns"),
+ action=draw.df.degree, title="Base spline",
+ parentname="DF")
+ rp.checkbox(panel, variable=quantile, action=draw.df.degree,
+ title="Show quantiles of x", parentname="DF")
+ rp.doublebutton(panel, variable=df1, action=draw.df.degree,
+ showvalue=TRUE, step=1, initval=3, range=c(1, 10),
+ title="Degrees of freedom", parentname="DF")
+ rp.doublebutton(panel, variable=degree1, action=draw.df.degree,
+ showvalue=TRUE, step=1, initval=3, range=c(1, 10),
+ title="Polynomial degree", parentname="DF")
+ rp.radiogroup(panel, variable=base2, vals=c("bs","ns"),
+ action=draw.degree.knots, title="Base spline",
+ parentname="Knots")
+ rp.doublebutton(panel, variable=degree2, action=draw.degree.knots,
+ showvalue=TRUE, step=1, initval=3, range=c(1, 10),
+ title="Polynomial degree", parentname="Knots")
+ rp.button(panel, action=reset.draw, title="Click to select knots",
+ parentname="Knots")
+ rp.do(panel, action=draw.df.degree)
+}
+
+rp.spline(x=cars$speed, y=cars$dist, er=0.3,
+ xlab="Velocidade", ylab="Distância")
+rp.spline(x=rock$area, y=rock$peri, er=0.2)
+rp.spline(x=faithful$eruptions, y=faithful$waiting, er=0.3)
+rp.spline(x=ca$hora, y=ca$temp, er=0.3)
+
+##-------------------------------------------
+## Base do spline.
+
+rp.base <- function(x){
+ draw <- function(panel){
+ with(panel, {
+ switch(base,
+ "bs"={
+ X <- model.matrix(~bs(x, df=df, degree=degree))
+ matplot(x=x, y=X[,-1], type="l")
+ },
+ "ns"={
+ X <- model.matrix(~ns(x, df=df))
+ matplot(x=x, y=X[,-1], type="l")
+ })})
+ panel
+ }
+ panel <- rp.control(x=x)
+ rp.radiogroup(panel, variable=base, vals=c("bs","ns"),
+ action=draw, title="Base spline")
+ rp.doublebutton(panel, variable=df, action=draw,
+ showvalue=TRUE, step=1, initval=3, range=c(1, 10),
+ title="Degrees of freedom")
+ rp.doublebutton(panel, variable=degree, action=draw,
+ showvalue=TRUE, step=1, initval=3, range=c(1, 10),
+ title="Polynomial degree")
+}
+
+x <- seq(0,10,0.1)
+rp.base(x)
+
+## matplot(x=x, y=X[,-1], type="l")
+
+##-----------------------------------------------------------------------------
+## Smoothing splines.
+
+rp.smooth.spline <- function(x, y, er=0, ...){
+ plotfit <- function(x, y, m0, er, ...){
+ xr <- extendrange(x, f=er)
+ yr <- extendrange(y, f=er)
+ xx <- seq(xr[1], xr[2], length.out=200)
+ plot(y~x, xlim=xr, ylim=yr, ...)
+ cb <- predict(m0, x=xx)
+ lines(cb$x, cb$y, lty=1, col=1)
+ }
+ annotations <- function(m0, y){
+ mtext(side=3, adj=0, line=1.5,
+ text=sprintf("Equiv. DF: %0.1f", m0$df))
+ mtext(side=3, adj=0, line=0.5,
+ text=sprintf("Minimized crit.: %0.3f", m0$crit))
+ mtext(side=3, adj=1, line=0.5,
+ text=sprintf("Smoothing parameter: %0.2f", m0$spar))
+ r2 <- 100*cor(fitted(m0), y)^2
+ mtext(side=3, adj=1, line=1.5,
+ text=sprintf("R²: %0.2f", r2))
+ }
+ draw.sspline <- function(panel){
+ with(panel,
+ {
+ switch(control,
+ "df"={
+ m0 <- smooth.spline(x=x, y=y, df=df)
+ },
+ "spar"={
+ m0 <- smooth.spline(x=x, y=y, spar=spar)
+ })
+ plotfit(x, y, m0, er, ...)
+ annotations(m0, y)
+ })
+ panel
+ }
+ draw.default <- function(panel){
+ m0 <- smooth.spline(x=x, y=y)
+ plotfit(x, y, m0, er, ...)
+ annotations(m0, y)
+ mtext(side=3, line=-1.5, text="Default fit", col=4)
+ panel
+ }
+ panel <- rp.control(x=x, y=y, er=er)
+ rp.do(panel, action=draw.default)
+ rp.radiogroup(panel, variable=control, vals=c("spar","df"),
+ labels=c("Smoothing parameter",
+ "Equivalent degrees of freedom"),
+ title="Argument in control", action=draw.sspline)
+ rp.slider(panel, variable=spar, from=0, to=1, initval=0.5,
+ action=draw.sspline, title="Smoothing parameter",
+ resolution=0.05, showvalue=TRUE)
+ rp.doublebutton(panel, variable=df,
+ range=c(1.5, length(x)-1),
+ initval=1.5, step=0.5, showvalue=TRUE,
+ action=draw.sspline,
+ title="Equivalent DF")
+ rp.button(panel, action=draw.default, title="Default fit")
+}
+
+rp.smooth.spline(cars$speed, cars$dist, er=0.3,
+ xlab="Velocidade", ylab="Comprimento")
+
+rp.smooth.spline(x=ca$hora, y=ca$temp, er=0.3)
+
+##-----------------------------------------------------------------------------
+## GAM.
+
+rp.gam <- function(y, x, er=0, ...){
+ annotations <- function(m0){
+ sm0 <- summary(m0)
+ mtext(side=3, adj=0, line=1.5,
+ text=sprintf("Equivalent DF: %0.2f", sm0$edf))
+ mtext(side=3, adj=0, line=0.5,
+ text=sprintf("Residual DF: %0.1f", sm0$residual.df))
+ ## mtext(side=3, adj=1, line=0.5,
+ ## text=sprintf("AIC: %0.2f", m0$aic))
+ mtext(side=3, adj=1, line=0.5,
+ text=sprintf("Smoothing parameter: %0.2f",
+ m0$smooth[[1]]$sp))
+ mtext(side=3, adj=1, line=1.5,
+ text=sprintf("R²: %0.2f",
+ 100*sm0$r.sq))
+ }
+ plotfit <- function(x, y, m0, er, ...){
+ b0 <- coef(m0)[1]
+ xr <- extendrange(x, f=er)
+ yr <- extendrange(y, f=er)
+ plot(m0, shift=b0, xlim=xr, ylim=yr-b0, ...)
+ points(y~x)
+ }
+ draw.gam <- function(panel){
+ with(panel,
+ {
+ m0 <- gam(y~s(x, k=k, sp=sp))
+ plotfit(x, y, m0, er, ...)
+ annotations(m0)
+ })
+ panel
+ }
+ draw.default <- function(panel){
+ with(panel,
+ {
+ m0 <- gam(y~s(x))
+ plotfit(x, y, m0, er, ...)
+ annotations(m0)
+ mtext(side=3, line=-1.5, text="Default fit", col=4)
+ })
+ panel
+ }
+ panel <- rp.control(x=x, y=y, er=er, ...)
+ maxd <- length(unique(x))-1
+ rp.doublebutton(panel, variable=k,
+ action=draw.gam, showvalue=TRUE,
+ step=1, initval=-1,
+ range=c(-1, maxd),
+ title="Degree of polynomial")
+ rp.slider(panel, variable=sp,
+ from=0, to=1, initval=0.1, action=draw.gam,
+ title="Smoothing parameter", showvalue=TRUE)
+ rp.button(panel, action=draw.default, title="Default fit")
+ rp.do(panel, action=draw.default)
+}
+
+rp.gam(x=cars$speed, y=cars$dist, er=0.3,
+ xlab="Velocidade", ylab="Comprimento")
+
+rp.gam(x=rock$peri, y=rock$perm, er=0.3)
+
+rp.gam(x=faithful$eruptions, y=faithful$waiting, er=0.3)
+rp.gam(x=ca$hora, y=ca$temp, er=0.3)
+
+##----------------------------------------------------------------------