<- read.table("../donnees/ozone_simple.txt",header=T,sep=";") ozone
16 Introduction à la régression spline
<- function(donnee,d=3){
polyreg <- sd(donnee[,"T12"])
sigmax <- seq(min(donnee[,"T12"])-sigmax,max(donnee[,"T12"])+sigmax,length=100)
grillex <- data.frame(T12=grillex)
aprevoir <- lm(O3~poly(T12,degree=d,raw=TRUE),data=donnee)
regpol <- predict(regpol,aprevoir)
prev return(list(grillex=grillex,grilley=prev))
}
plot(O3~T12,data=ozone,xlab="T12",ylab="O3",xlim=c(0,35), ylim=c(0,150))
<- 1
iter for(ii in c(1,2,3,9)){
<- polyreg(ozone,d=ii)
tmp lines(tmp$grillex,tmp$grilley,col=iter,lty=iter)
<- iter+1
iter
}legend("topleft",c("d=1","d=2","d=3","d=9"),col=1:4,lty=1:4)
<- which(ozone[,2]<23)
ind <- lm(O3~T12,data=ozone[ind,])
regd <- lm(O3~T12,data=ozone[-ind,])
regf <- seq(3,23,length=50)
gxd <- regd$coef[1]+gxd*regd$coef[2]
gyd <- seq(23,35,length=50)
gxf <- regf$coef[1]+gxf*regf$coef[2]
gyf plot(O3~T12,data=ozone)
lines(gxd,gyd,col=2,lty=1,lwd=2)
lines(gxf,gyf,col=2,lty=1,lwd=2)
abline(v=23)
library(splines)
<- bs(ozone[,2], knots=c(15,23), degree=2,Boundary.knots=c(5,32))
XB <- lm(ozone[,"O3"] ~ XB)
regs $coef regs
(Intercept) XB1 XB2 XB3 XB4
51.101947 61.543761 5.562286 70.459103 106.711539
<- seq(5,32,length=100)
grillex <- bs(grillex, knots=c(15,23), degree=2,Boundary.knots=c(5,32))
bgrillex <- bgrillex%*%as.matrix(regs$coeff[-1])+regs$coeff[1]
prev plot(O3~T12,data=ozone)
lines(grillex,prev,col=2)
abline(v=c(15,23))
<- smooth.spline(ozone[,2],ozone[,1],lambda =100)
regssplinel1 <- predict(regssplinel1,grillex)
prevl1 plot(O3~T12,data=ozone)
lines(prevl1$x,prevl1$y,col=2)
<- smooth.spline(ozone[,2],ozone[,1])
regsspline <- predict(regsspline,grillex)
prev plot(O3~T12,data=ozone)
lines(prev$x,prev$y,col=2)
regsspline
Call:
smooth.spline(x = ozone[, 2], y = ozone[, 1])
Smoothing Parameter spar= 0.9410342 lambda= 0.006833357 (15 iterations)
Equivalent Degrees of Freedom (Df): 4.156771
Penalized Criterion (RSS): 11036.88
GCV: 289.8012