Wie übersetzt man John Snow-Datensätze in Kartenkoordinaten?

14

Im HistData-Paket für R ( https://r-forge.r-project.org/R/?group_id=574 ) sind die Datensätze enthalten, die sich auf John Schnees Karte des Cholera-Ausbruchs in London von 1854 beziehen sind maßgebend, sorgfältig digitalisiert unter Aufsicht von Walter Tobler. Einige Details zu diesen Datensätzen werden von John Mackenzie unter http://www1.udel.edu/johnmack/frec480/cholera/cholera2.html beschrieben .

Leider verwenden die Koordinaten von Todesfällen, Pumpen und Straßen ein willkürliches Koordinatensystem, keine Kartenkoordinaten, die für andere GIS-Anwendungen oder Kartierungssoftware in R geeignet sind (räumliche Pakete, ggmap usw.).

In http://freakonometrics.hypotheses.org/19213 verwendet Arthur Charpentier ggmap mit einer Version der John Snow-Daten von http://www.rtwilson.com/downloads/SnowGIS_v2.zip . Die Cholera_Deaths.shpDatei listet jedoch nur 489 Todesfälle auf, nicht die 578, in denen ich aufgezeichnet habe HistData::Snow.deaths.

Eine Idee ist, die Beziehungen zwischen den Mitteln und den Standardabweichungen der (x, y) -Koordinaten zu finden und linear neu zu skalieren, aber vielleicht gibt es einen besseren Weg?

Folgendes habe ich bisher versucht

> 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

Lesen Sie die Cholera_deaths-Datei

> 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)
[1] 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, dann versuche ich, neu zu skalieren D, um die gleichen XXMittelwerte und Standardabweichungen wie zu erhalten , aber hier funktioniert etwas nicht richtig - die Spaltenmittelwerte von Dscaledhätten gleich denen von sein müssen XX:

> # 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 
>

BEARBEITEN: Bei diesem Problem kann es hilfreich sein, die Karte von Snow so zu sehen, wie sie von der neuen Funktion SnowMap(axis.labels=TRUE)in der Entwicklungsversion von HistData(Rev. 102) auf R-Forge gezeichnet wurde . Die Achsenbeschriftungen zeigen den Ursprung des Koordinatensystems in der unteren linken Ecke wie in meinen Snow.*Datendatensätzen.

Schneekarte

user101089
quelle
Ich habe gerade versucht, die Pumpen aus jedem Datensatz so zu skalieren, dass sie übereinstimmen. Ich glaube nicht, dass die Hilfslinie (Snow.pumps) für die Koordinaten 100 Meter beträgt, da eine Skala von etwa 54 (und eine Übersetzung) diese am besten auf die britischen Gitterkoordinaten des Shapefiles abbildet (die definitiv in sind) meter). Selbst wenn sich die Punkte nicht genau überlappen, ist eine andere Drehung / Neigung deutlich vorhanden. Da es weniger Pumpen gibt, ist es möglich, entsprechende Pumpen in jedem Datensatz zu identifizieren und die Verschiebung / Übersetzung für diese zu berechnen.
Spacedman
Ich nehme an, Sie haben sich HistData / inst / doc / Snow_deaths-duplicates.html angesehen und fanden es nicht hilfreich?
Barrycarter
Mir fiel auch ein, dass ich die lineare Umwandlung der Koordinaten in meinen Snow.*Dateien in eine GIS-basierte Karte mit den Positionen von zwei oder drei Pumpen erhalten konnte, um die Genauigkeit zu überprüfen. Leider befinden sich in den SnowGISDateien keine Etiketten für die Pumpen , und ich habe kein Beispiel dafür gesehen, wie man sie plottet, damit ich sie visuell vergleichen kann.
user101089
1
Nachdem ich Ihren Titel gelesen hatte, dachte ich für eine Sekunde, Sie wollten die Koordinaten in Westeros kartieren .
user35594

Antworten:

4

Vielleicht bewerten Sie das Shapefile von http://donboyes.com/2011/10/14/john-snow-and-serendipity es hat 578 Punkte.

Ich glaube nicht, dass der Versuch, die HistData-Schneetodesfälle mit Robin Wilsons (@robintw) Version in Beziehung zu setzen, funktioniert, da das Shapefile eine einzelne Punktkoordinate für mehrere Todesfälle an einer einzelnen Adresse enthält, anstatt der mehreren Punkte, die von der Straße in der gestapelt wurden Karte .

Robins Version fehlt definitiv eine Menge Punkte. Ein kurzer Blick zeigt, dass nicht wenige einzelne Todesfälle übersehen wurden. Ein weiteres Problem liegt näher an der Kartenmitte, wo es beim Zusammenfügen nicht richtig kantenangepasst war (dies ist auch in der Wikipedia- Karte sichtbar ) und dies verdeckt eine Reihe von Punkten.

Ausschnitt der im Download mitgelieferten Karte :

Bildbeschreibung hier eingeben

Auszug aus der UCLA-Version :

Bildbeschreibung hier eingeben

user2856
quelle
Groß! Der Vollständigkeit.shp
halber
2

Um die Antwort auf diese Frage zu vervollständigen, findet der folgende Code die lineare Transformation der Koordinaten in den ursprünglichen Tobler-Dateien (in HistData) und den von Don Boyes bereitgestellten.

folder <- "C:/Dropbox/R/data/Snow/snow_shp"
library(maptools)
deaths <- readShapePoints(file.path(folder, "deaths_gcs"))
data(Snow.deaths, package="HistData")
X <- deaths@coords
D <- Snow.deaths[,2:3]

Korrelieren Sie dann D [, 1] auf X [, 1] und D [, 2] auf X [, 2] und bilden Sie eine Regression aus. Die lineare Transformation ist durch die Regressionskoeffizienten gegeben.

> cor(D[,1], X[,1])
[1] 0.9999664
> cor(D[,2], X[,2])
[1] 0.9995559
> 
> # linear transformations to GIS coords
> lm(D[,1] ~ X[,1])

Call:
lm(formula = D[, 1] ~ X[, 1])

Coefficients:
(Intercept)       X[, 1]  
      185.4       1264.7  

> 
> lm(D[,2] ~ X[,2])

Call:
lm(formula = D[, 2] ~ X[, 2])

Coefficients:
(Intercept)       X[, 2]  
    -105441         2047  
user101089
quelle