haiti_analysis.R

library(geoR)
library(geoRglm)

# Perfrom MCMC simulations
model<-list(cov.pars = c(1, 1), beta = 1, family = "poisson")
mcmc.test<-mcmc.control(S.scale = 0.45, thin = 1)
test.tune<-glsm.mcmc(haiti.geo, model = model, mcmc.input=mcmc.test)
haiti.mcmc<-prepare.likfit.glsm(test.tune)
prior<- prior.glm.control(phi.prior = "fixed", phi = .1)

# Create grid over the Port-au-Prince
min.long<-min(geo.coords[,1])
max.long<-max(geo.coords[,1])
min.lat<-min(geo.coords[,2])
max.lat<-max(geo.coords[,2])
grid.loc<-expand.grid(x = seq(min.long, max.long, l = 100),y = seq(min.lat, max.lat, l = 100))

# Generate probabiltiies for Category 1 reports in each grid
# location and plot choropleth
pkb<-pois.krige.bayes(haiti.geo, locations=grid.loc,prior = prior, mcmc.input = mcmc.test)
png("haiti_choropleth.png",width=1200,height=800,res=120)
plot(haiti.geo$coords,cex=0.4,col="black",main="Choropleth of Probability for Category 1\n(Emergncy) Ushahidi Post in Port-au-Prince",
    xlab="Longitude",ylab="Latitude")
image(pkb,col=rev(heat.colors(50,alpha=0.75)),add=TRUE)
dev.off()