**Pirate Science » R**, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

**We want to evaluate the seriousness of the threat** posed by the recent ebola outbreak in western Uganda. The outbreak appeared in Kigadi, a small village in the Kibaale district. The disease was first confirmed by the government on 28 July. By 11 Aug, Medecins Sans Frontieres epidemiologist Dr Paul Roddy said the disease was under control, with no new deaths in August. It will be classified as contained if there are no confirmed cases by middle September. In total, 19 people died from the disease in the outbreak.

The disease spreads between people by contact with contaminated body fluids. Transmission, then, is slow (in comparison to disease with aerosol transmission, i.e. influenza). We reasoned that if the disease was to spread to other communities, the most likely route would be along roadways.

We obtained a shapefile describing the Uganda road network from Humanitarian Response. We plot it as follows:

### use the maptools library

require(maptools)## load the data:

ugnat <- readShapePoly('Uganda_national2006')

ugroads <- readShapeLines('Uganda_Roads_Feb2009')

ugtowns <- readShapePoints('Uganda_Towns')

ugvil <- readShapePoints('Uganda_villages_27Jan09')## plot it

png(“ugroads.png”,width=1000,height=1000)

plot(ugnat,lwd=3)

plot(ugroads[ugroads$RDLNTYPETX==”Unpaved Trails/Footpaths”,],add=TRUE,lwd=1,col=’lightgray’)

plot(ugroads[ugroads$RDLNTYPETX==”Tracks trails or footpaths”,],add=TRUE,lwd=1,col=’lightgray’)

plot(ugroads[ugroads$RDLNTYPETX==”Secondary road”,],add=TRUE,lwd=2,col=’burlywood3′)

plot(ugroads[ugroads$RDLNTYPETX==”Primary and secondary roads”,],add=TRUE,lwd=3,col=’burlywood3′)

plot(ugroads[is.na(ugroads$RDLNTYPETX),],add=TRUE,lwd=2,col=’burlywood3′)

plot(ugroads[ugroads$RDLNTYPETX==”Connector”,],add=TRUE,lwd=2,col=’burlywood3′)

plot(ugtowns,pch=16,col=colors()[131],cex=1,add=TRUE)

kagadi <- grep('Kagadi',ugvil$FULL_NAME,ignore.case=TRUE)

points(ugvil[kagadi,],pch=1,col=’red’,cex=1.5)

points(ugvil[kagadi,],pch=1,col=’red’,cex=2.5)

points(ugvil[kagadi,],pch=1,col=’red’,cex=3.5)

points(ugvil[kagadi,],pch=1,col=’red’,cex=4.5)

points(ugvil[kagadi,],pch=16,col=’red’,cex=0.75)

dev.off()

Converting it to a graph requires matching endpoints of the roads. Simple enough, except that the road coordinates are buried three classes deep in the data structure. We write two utility functions:

gthead <- function(indx){

return([email protected][indx][[1]]@Lines[[1]]@coords[1,])}gttail <- function(indx){

return(tail([email protected][indx][[1]]@Lines[[1]]@coords,n=1))}mkhash <- function(coord){

return(paste(coord[1],coord[2],sep=’_’))}

Then, to make the graph

nroads <- dim(ugroads)[1]

edges <- matrix('',nrow=nroads,ncol=2)for(i in 1:nroads){

edges[i,1] <- mkhash(gthead(ugroads,i))

edges[i,2] <- mkhash(gttail(ugroads,i))

}gr <- graph.edgelist(edges,directed=FALSE)

where the hash function pastes together the x and y coordinates of an endpoint. The resulting graph has 2011 nodes and 1992 edges. The maximum degree node has degree 7 (if one allows a maximum of one edge between any two nodes). Almost half of the nodes have degree one — they are terminal. Analysis of the graph shows that none of the nodes have strong spreading power, whether measured by *expected reach*, path counting, degree, or k-shell.

The model has obvious limitations. Nodes are road junctures, not population centers. No allowance is made for road length or size. our conclusion should be taken with a grain of salt.

**Conclusion: The disease is a local tragedy, but is unlikely to become a national one.**

**leave a comment**for the author, please follow the link and comment on their blog:

**Pirate Science » R**.

R-bloggers.com offers

**daily e-mail updates**about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.