In the HistData package for R (https://r-forge.r-project.org/R/?group_id=574) I have the data sets related to John Snow's map of the cholera outbreak in London, 1854. I believe they are authoritative, having been carefully digitized under supervision of Walter Tobler. Some details on these data sets are described by John Mackenzie, at http://www1.udel.edu/johnmack/frec480/cholera/cholera2.html.
Unfortunately, the coordinates of deaths, pumps, and streets use an arbitrary coordinate system, not map coordinates suitable for other GIS applications or mapping software in R (spatial packages, ggmap, etc.)
In http://freakonometrics.hypotheses.org/19213 Arthur Charpentier
uses ggmap with a version of the John Snow data from
Cholera_Deaths.shp file, however lists only 489 deaths, not the 578 I have recorded in
One idea is to find the relations between the means and standard deviations of the (x,y) coordinates and rescale linearly, but maybe there is a better way?
Here's what I've tried so far
> data(Snow.deaths, package="HistData") > D <- Snow.deaths[,2:3] > colMeans(D) x y 13.03312 11.69721 > var(D) x y x 3.8150987 0.3802654 y 0.3802654 2.7213828
Read the Cholera_deaths file
> folder <- "C:/Dropbox/R/data/Snow/SnowGIS_v2/SnowGIS" > library(maptools) > deaths <- readShapePoints(file.path(folder, "Cholera_Deaths")) > head(deaths@coords) coords.x1 coords.x2 0 529308.7 181031.4 1 529312.2 181025.2 2 529314.4 181020.3 3 529317.4 181014.3 4 529320.7 181007.9 5 529336.7 181006.0 > # deaths has only 250 observations; 489 deaths > sum(deaths@data$Count)  489 > # try to relate to Snow.deaths > X <- deaths@coords > colnames(X) <- c("x", "y") > > XX <- data.frame(X, Freq=deaths@data$Count) > XX <- vcdExtra::expand.dft(XX) > > colMeans(XX) x y 529414.8 181031.9 > var(XX) x y x 10813.816 1521.693 y 1521.693 6227.924 >
OK, then I try to re-scale
D to have the same means and standard
XX, but something doesn't work correctly here-- the column means of
Dscaled should have turned out to be equal to those of
> # scale D to have the same means and standard deviations as XX > Dscaled <- scale(D, center=TRUE, scale=TRUE) > Dscaled <- scale(Dscaled, center=colMeans(XX), scale=sqrt(diag(var(XX)))) > colMeans(Dscaled) x y -5091.040 -2293.947 >
EDIT: It might be helpful in this problem to see Snow's map as drawn by the new function,
SnowMap(axis.labels=TRUE) now in the development version of
HistData (rev 102) on R-Forge. The axis labels show the origin of the coordinate system in the lower left corner as they are in my data
Snow.* data sets.