# do not forget to unzip neigh-parana.zip

require(aRT)
aRTsilent(FALSE)
con=openConn()
con

dbname="testeparana"

if(any(showDbs(con)==dbname))
  deleteDb(con, dbname, force=T)

db=createDb(con,db=dbname)
db

l=createLayer(db, l="cidades")
addShape(l, tab="tab", file="neigh-parana.shp", ID="ID_IBGE")
l

getPolygons(l)->pols

rel=getRelation(l,"touches",ID=c("412855","410530"))
rel2=getRelation(l,"touches",ID=rel)
rel3=getRelation(l,"touches",ID=rel2)

colors=terrain.colors(5)
plot(pols[rel3,],col=colors[4])
plot(pols[rel2,],col=colors[3],add=T)
plot(pols[rel,],col=colors[2],add=T)
plot(pols[c("412855","410530"),],col=colors[1],add=T)

ids = getID(pols)
ids

aRTsilent(TRUE)
cat("Calculating relations ... ")
li=lapply(ids, function(x) {cat(".");which(is.element(ids,getRelation(l,"touches",ID=x)))})
cat(" done\n")
aRTsilent(FALSE)
li

getOperation(l,"centroid")->centroid
pols=getPolygons(l)

par( mar=c(.0,.0,.0,.0))

plot(pols, col="gray")

# implement an algorithm to colour the map using four colours.
#lapply(1:399, function(x)
#aRTplot(SpatialPolygons(list(pols@polygons[[x]]), pO=1), col=terrain.colors(5)[x%%5], lty=0, add=T))

#aRTplot(centroid,add=T)

lapply(1:length(li), function(pos)
lapply(li[[pos]], function(i) plot(SpatialLines(list(Lines(Line(centroid[c(pos,i),])))), col="blue", add=T)))

li

t = openTable(l,"tab")
aRT::getData(t)[,c("ID_IBGE","NOME_MUNIC","NASCIDOS","OBITOS")]->data

head(data)
names(data) <- c("ID_IBGE","nome","nascidos","obitos")
head(data)
attach(data)
nome

t = createTable(l,"resultados",ID="ID_IBGE")
t

dadosnb=li

dadosnb

library(spdep)
library(RColorBrewer)

res=probmap(obitos,nascidos)
head(res)
df=cbind(ID_IBGE,res)
head(df)
updateColumns(t, df)
t
                            

##variável raw
brks=c(-Inf,4,13,20,28,Inf)
nclr=length(brks)-1
plotclr=brewer.pal(nclr,"Greens")
plotclr

th=createTheme(l, "raw")
th
v = visualRaster(col=plotclr)
v@color
setVisual(th,v,att="raw",mode="q")


th=createTheme(l, "relRisk")
th
v = visualRaster(col=plotclr)
v@color
setVisual(th,v,att="relRisk",mode="q")

#Bayes empírico
res=EBest(obitos,nascidos)

novodado=data.frame(ID_IBGE, estmm=res$estmm)
novodado
updateColumns(t,novodado)
aRT::getData(t)

th=createTheme(l,"estmm2")
th
v = visualRaster(col=plotclr)
v@color
setVisual(th,v,att="estmm",mode="q")

class(dadosnb)="nb"

#Bayes empírico local
res=EBlocal(obitos,nascidos,dadosnb,zero.policy=TRUE)

res


novodado=data.frame(ID_IBGE, est=res$est)
head(novodado)
updateColumns(t,novodado)
t
aRT::getData(t)

plotclr=brewer.pal(5,"Oranges")

th=createTheme(l, "est3")
v = visualRaster(col=plotclr)
v@color
setVisual(th,v,att="est",mode="s")