This is a simple space-time analysis on hand-food-mouth desease case data using R packages such as GISTools, rgeos, spacetime, etc..
library(sp)
library(spacetime)
library(foreign)
library(grid)
library(lattice)
library(RColorBrewer)
library(MASS)
library(rgeos)
library(maptools)
library(maps)
library(zoo)
library(GISTools)
SDtown <- readShapePoints(file.choose())
SDcounty <- readShapePoly(file.choose())
proj4string(SDtown) = CRS(as.character(NA))
proj4string(SDcounty) = CRS(as.character(NA))
case20096 = read.csv("f:/2Work/spacetime/Data/2009_6_R.csv", header = TRUE,
sep = ",")
caseData = read.csv("f:/2Work/spacetime/Data/cases36mon.csv", header = TRUE,
sep = ",")
Display original disease case data at town level:
shades = auto.shading(case20096$CaseNumber)
choropleth(SDtown, case20096$CaseNumber, shading = shades, pch = 16)
plot(SDcounty, add = TRUE)
choro.legend(9520000, 3550000, shades)
Aggregate disease case data at county level:
points = gContains(SDcounty, SDtown, byid = TRUE)
pointVs = replace(points, points == FALSE, NA) * case20096$CaseNumber
caseSum = colSums(pointVs, na.rm = TRUE)
labs = 1:6
cuts = c(100, 200, 300, 400, 600, 800)
names(cuts) = paste("#", labs, sep = "")
cols = brewer.pal(7, "YlOrRd")
shades = list(breaks = cuts, cols = cols)
class(shades) = "shading"
This is the spatial aggragation result in June 2009
choropleth(SDcounty, caseSum, shading = shades)
choro.legend(9520000, 3550000, shades)
Spatial aggragation in 36 months.
tcaseData = t(caseData)
values = tcaseData[4:39, 1]
for (i in 2:1884) {
values = rbind(values, tcaseData[4:39, i])
}
IDs36 = paste("ID", 1:(36 * 1884), sep = "_")
time36 = seq(as.yearmon("2009-01"), as.yearmon("2011-12"), 1/12)
value36 = values[1:(36 * 1884)]
mydata36 = data.frame(value36, ID = IDs36)
stfdf36 = STFDF(SDtown, time36, mydata36)
county36 = aggregate(stfdf36[, , "value36"], SDcounty, sum, na.rm = TRUE)
Spatial aggragation result of the 36 months.
stplot(county36, col.regions = brewer.pal(7, "YlOrRd"), cuts = 7, at = c(0,
100, 200, 300, 400, 600, 800, 1200), lty = 1, lwd = 0.5, cex = 0.1, col = "lightgrey")