#lip<-read.table("c:\\slip.txt",header=TRUE,sep=",",row.names=1)

lip$observed<-c(9,39,11,9,15,8,26,7,6,20,13,5,3,8,17,9,2,7,9,7,16,31,11,7,19,15,7,10,16,11,
5,3,7,8,11,9,11,8,6,4,10,8,2,6,19,3,2,3,28,6,1,1,1,1,0,0)

lip$expected<-c(1.4, 8.7, 3.0, 2.5, 4.3, 2.4, 8.1, 2.3, 2.0, 6.6,4.4, 1.8, 1.1, 3.3, 7.8,
4.6, 1.1, 4.2, 5.5, 4.4,10.5,22.7, 8.8, 5.6,15.5,12.5, 6.0, 9.0,14.4,10.2,4.8,2.9,7.0,
8.5,12.3,10.1,12.7, 9.4, 7.2, 5.3,18.8,15.8, 4.3,14.6,50.7, 8.2, 5.6, 9.3,88.7,19.6,3.4,
3.6, 5.7, 7.0, 4.2, 1.8)

adj1<-c(5,9,11,19)
adj2<-c(7,10)
adj3<-c(6,12)
adj4<-c(18,20,28)
adj5<-c(1,11,12,13,19)
adj6<-c(3,8)
adj7<-c(2,10,13,16,17)
adj8<-c(6)
adj9<-c(1,11,17,19,23,29)
adj10<-c(2,7,16,22)
adj11<-c(1,5,9,12)
adj12<-c(3,5,11)
adj13<-c(5,7,17,19)
adj14<-c(31,32,35)
adj15<-c(25,29,50)
adj16<-c(7,10,17,21,22,29)
adj17<-c(7,9,13,16,19,29)
adj18<-c(4,20,28,33,55,56)
adj19<-c(1,5,9,13,17)
adj20<-c(4,18,55)
adj21<-c(16,29,50)
adj22<-c(10,16)
adj23<-c(9,29,34,36,37,39)
adj24<-c(27,30,31,44,47,48,55,56)
adj25<-c(15,26,29)
adj26<-c(25,29,42,43)
adj27<-c(24,31,32,55)
adj28<-c(4,18,33,45)
adj29<-c(9,15,16,17,21,23,25,26,34,43,50)
adj30<-c(24,38,42,44,45,56)
adj31<-c(14,24,27,32,35,46,47)
adj32<-c(14,27,31,35)
adj33<-c(18,28,45,56)
adj34<-c(23,29,39,40,42,43,51,52,54)
adj35<-c(14,31,32,37,46)
adj36<-c(23,37,39,41)
adj37<-c(23,35,36,41,46)
adj38<-c(30,42,44,49,51,54)
adj39<-c(23,34,36,40,41)
adj40<-c(34,39,41,49,52)
adj41<-c(36,37,39,40,46,49,53)
adj42<-c(26,30,34,38,43,51)
adj43<-c(26,29,34,42)
adj44<-c(24,30,38,48,49)
adj45<-c(28,30,33,56)
adj46<-c(31,35,37,41,47,53)
adj47<-c(24,31,46,48,49,53)
adj48<-c(24,44,47,49)
adj49<-c(38,40,41,44,47,48,52,53,54)
adj50<-c(15,21,29)
adj51<-c(34,38,42,54)
adj52<-c(34,40,49,54)
adj53<-c(41,46,47,49)
adj54<-c(34,38,49,51,52)
adj55<-c(18,20,24,27,56)
adj56<-c(18,24,30,33,45,55)

adj<-matrix(nrow=56,ncol=56)

adj[,1]<-makeadj(adj1,56)
adj[,2]<-makeadj(adj2,56)
adj[,3]<-makeadj(adj3,56)
adj[,4]<-makeadj(adj4,56)
adj[,5]<-makeadj(adj5,56)
adj[,6]<-makeadj(adj6,56)
adj[,7]<-makeadj(adj7,56)
adj[,8]<-makeadj(adj8,56)
adj[,9]<-makeadj(adj9,56)
adj[,10]<-makeadj(adj10,56)
adj[,11]<-makeadj(adj11,56)
adj[,12]<-makeadj(adj12,56)
adj[,13]<-makeadj(adj13,56)
adj[,14]<-makeadj(adj14,56)
adj[,15]<-makeadj(adj15,56)
adj[,16]<-makeadj(adj16,56)
adj[,17]<-makeadj(adj17,56)
adj[,18]<-makeadj(adj18,56)
adj[,19]<-makeadj(adj19,56)
adj[,20]<-makeadj(adj20,56)
adj[,21]<-makeadj(adj21,56)
adj[,22]<-makeadj(adj22,56)
adj[,23]<-makeadj(adj23,56)
adj[,24]<-makeadj(adj24,56)
adj[,25]<-makeadj(adj25,56)
adj[,26]<-makeadj(adj26,56)
adj[,27]<-makeadj(adj27,56)
adj[,28]<-makeadj(adj28,56)
adj[,29]<-makeadj(adj29,56)
adj[,30]<-makeadj(adj30,56)
adj[,31]<-makeadj(adj31,56)
adj[,32]<-makeadj(adj32,56)
adj[,33]<-makeadj(adj33,56)
adj[,34]<-makeadj(adj34,56)
adj[,35]<-makeadj(adj35,56)
adj[,36]<-makeadj(adj36,56)
adj[,37]<-makeadj(adj37,56)
adj[,38]<-makeadj(adj38,56)
adj[,39]<-makeadj(adj39,56)
adj[,40]<-makeadj(adj40,56)
adj[,41]<-makeadj(adj41,56)
adj[,42]<-makeadj(adj42,56)
adj[,43]<-makeadj(adj43,56)
adj[,44]<-makeadj(adj44,56)
adj[,45]<-makeadj(adj45,56)
adj[,46]<-makeadj(adj46,56)
adj[,47]<-makeadj(adj47,56)
adj[,48]<-makeadj(adj48,56)
adj[,49]<-makeadj(adj49,56)
adj[,50]<-makeadj(adj50,56)
adj[,51]<-makeadj(adj51,56)
adj[,52]<-makeadj(adj52,56)
adj[,53]<-makeadj(adj53,56)
adj[,54]<-makeadj(adj54,56)
adj[,55]<-makeadj(adj55,56)
adj[,56]<-makeadj(adj56,56)

rrgp<-gampoi(niter=20,obs=lip$observed,exd=lip$expected)$smrhat
rrlb<-global.lb(obs=lip$observed,exd=lip$expected)
rrllb<-local.lb(obs=lip$observed,exd=lip$expected,adj=adj)

par(mfrow=c(3,1))
plot(rrgp,rrlb,xlab="gamma-poisson",ylab="global bayes")
plot(rrgp,rrllb,xlab="gamma-poisson",ylab="local bayes")
plot(rrlb,rrllb,xlab="local bayes",ylab="global bayes")