Difference between revisions of "Bsplines.R"

From Organic Design wiki
(library(splines))
(Seems to be a bug in variable naming from 'fit')
Line 1: Line 1:
 
# {{R}}
 
# {{R}}
 
library(splines)
 
library(splines)
 
+
 
# p38
 
# p38
 
pos <- function(x) ifelse(x>0, x, 0)
 
pos <- function(x) ifelse(x>0, x, 0)
Line 8: Line 8:
 
knot <- 4
 
knot <- 4
 
plot(x,y)
 
plot(x,y)
X <- cbind("int"=1,"x"=x,"x2"=x^2,"x3"=x^3, "xknot"=pos(x-knot)^3)
+
#X <- cbind("int"=1,"x"=x,"x2"=x^2,"x3"=x^3, "xknot"=pos(x-knot)^3)
fit <- lm(y~X-1)
+
#fit <- lm(y~X-1)
 +
 
 +
fit <- lm(y~1+x+x^2+x^3+pos(x-knot)^3)
 
xx <- seq(1,7, length=200)
 
xx <- seq(1,7, length=200)
XX <- cbind("Xint"=1,"Xx"=xx, "Xx2"=xx^2, "Xx3"=xx^3, "Xxknot"=pos(xx-knot)^3)
+
XX <- cbind("x"=xx)
 
lines(xx, predict(fit, as.data.frame(XX))) # Bug here
 
lines(xx, predict(fit, as.data.frame(XX))) # Bug here
 
abline(v=knot, lty=2)
 
abline(v=knot, lty=2)
 
X
 
X
 +
 
# p41
 
# p41
 
par(mfrow=c(2,1))
 
par(mfrow=c(2,1))

Revision as of 02:23, 17 October 2007

Code snipits and programs written in R, S or S-PLUS library(splines)

  1. p38

pos <- function(x) ifelse(x>0, x, 0) x <- 1:7 y <- c(8,3,8,5,9,14,11) knot <- 4 plot(x,y)

  1. X <- cbind("int"=1,"x"=x,"x2"=x^2,"x3"=x^3, "xknot"=pos(x-knot)^3)
  2. fit <- lm(y~X-1)

fit <- lm(y~1+x+x^2+x^3+pos(x-knot)^3) xx <- seq(1,7, length=200) XX <- cbind("x"=xx) lines(xx, predict(fit, as.data.frame(XX))) # Bug here abline(v=knot, lty=2) X

  1. p41

par(mfrow=c(2,1)) knots <- c(1:3,5,7,8,10) atx <- seq(0:11, by=0.1) for(ord in 1:4) { B <- bs(x=atx, degree=ord, knots=knots, intercept=TRUE, Boundary.knots=c(0,11)) matplot(atx, B[,1], type="l", ylim=0:1, lty=2) matlines(atx, B[,-1], col=1, lty=1) title(paste("B splines of order", ord)) abline(v=knots, lty=2) }