## library(igraph)
library(rpart)

######
## GAM example
######
set.seed(3)
n=50

x1<-runif(n,7*pi,13*pi)  ##7,13
x2<-runif(n,0,2*pi)
f<-function(x,y){
	
  f1<-sin(2*y)+pi/2
  f2<-exp(x/pi)
  f3<-1/(1+cos(y))
  eta=f1/f3+log(1+f2) 
  eta-mean(eta)
}
eta=f(x1,x2)+rnorm(n,0,1)
y=round( exp(eta)/(1+exp(eta)))


####
## Fit your GAM
####

## do the initialization
myX <- cbind(x1,x2)
myBeta = rep(0, 50)
myEta = myX*myBeta

## start iterating
for(i in 1:250){
	myP = exp(myEta)/(1+exp(myEta))
	myTmpW =(myP*(1-myP))
	myW1 = diag(myTmpW[,1])
	myW2 = diag(myTmpW[,2])
	myZ1 = myEta - solve(myW1)*(y-myP[1])
	myBeta = solve(t(myX)*myW*(myX))*t(myX)*myW*myZ
	myEta = myX*myBeta
}

####
## Plot the true border
####
plot(x1,x2,col=y+2, pch=16, cex=2, xlab="x1", ylab="x2", ylim=c(0,8))##,xlim=c(min(x1),45))
x11<-seq(7*pi,13*pi,length=100)
x22<-seq(0,8,length=100)
z=outer(x11,x22,FUN=f)
contour(x11,x22,z,levels=0,add=TRUE,lty=ltys[1],method="edge",drawlabels=FALSE,col="purple",lwd=2)

####
## Use your gam to predict the x11, x22 grid
####
lines(myEta)


###
## Plot your prediction Border
####

title("Binary Simulated Data")