bone = read.table("../data/bone.data",header=T)

library(splines)
d.pred = data.frame(age = sort(bone$age))
fit.m = lm(spnbmd~ns(age,df=12),data=bone,subset=(gender=="male"))
pred.m = predict(fit.m,d.pred)
fit.f = lm(spnbmd~ns(age,df=12),data=bone,subset=(gender=="female"))
pred.f = predict(fit.f,d.pred)

plot(bone$age,bone$spnbmd,col=as.numeric(bone$gender=="female")+1)
lines(d.pred$age,pred.m,type="l",col=1)
lines(d.pred$age,pred.f,type="l",col=2)
legend("topright",c("male","female"),lty=1,col=1:2)

N = model.matrix(~ns(age,df=12)-1,data=bone,subset=(gender=="male"))
lambda = 0.00022
S = N%*%solve(t(N)%*%N+lambda*diag(rep(1,12)))%*%t(N)
show(sum(diag(S)))

bone.m = bone[bone$gender=="male",]
fit.m.s = smooth.spline(bone.m$age,bone.m$spnbmd,df=12)
bone.f = bone[bone$gender=="female",]
fit.f.s = smooth.spline(bone.f$age,bone.f$spnbmd,df=12)
plot(bone$age,bone$spnbmd,col=as.numeric(bone$gender=="female")+1)
lines(fit.m.s,col=1,lwd=2)
lines(fit.f.s,col=2,lwd=2)
legend("topright",c("male","female"),lty=1,col=1:2)