J. Snow Cholera Deaths (1854)
Data
mdsr(Modern Data Science with R) 패키지와 공간 데이터 분석 전용 sp 패키지 설치.
library(mdsr) ## Loading tidyverse: ggplot2
## Loading tidyverse: tibble
## Loading tidyverse: tidyr
## Loading tidyverse: readr
## Loading tidyverse: purrr
## Loading tidyverse: dplyr## Conflicts with tidy packages ----------------------------------------------## filter(): dplyr, stats
## lag(): dplyr, stats## Loading required package: lattice## Loading required package: ggformula##
## New to ggformula? Try the tutorials:
## learnr::run_tutorial("introduction", package = "ggformula")
## learnr::run_tutorial("refining", package = "ggformula")## Loading required package: mosaicData##
## The 'mosaic' package masks several functions from core packages in order to add
## additional features. The original behavior of these functions should not be affected by this.
##
## Note: If you use the Matrix package, be sure to load it BEFORE loading mosaic.library(sp)
head(CholeraDeaths)## coordinates Id Count
## 1 (529308.7, 181031.4) 0 3
## 2 (529312.2, 181025.2) 0 2
## 3 (529314.4, 181020.3) 0 1
## 4 (529317.4, 181014.3) 0 1
## 5 (529320.7, 181007.9) 0 4
## 6 (529336.7, 181006) 0 2str(CholeraDeaths)## Formal class 'SpatialPointsDataFrame' [package "sp"] with 5 slots
## ..@ data :'data.frame': 250 obs. of 2 variables:
## .. ..$ Id : int [1:250] 0 0 0 0 0 0 0 0 0 0 ...
## .. ..$ Count: int [1:250] 3 2 1 1 4 2 2 2 3 2 ...
## ..@ coords.nrs : num(0)
## ..@ coords : num [1:250, 1:2] 529309 529312 529314 529317 529321 ...
## .. ..- attr(*, "dimnames")=List of 2
## .. .. ..$ : NULL
## .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
## ..@ bbox : num [1:2, 1:2] 529160 180858 529656 181306
## .. ..- attr(*, "dimnames")=List of 2
## .. .. ..$ : chr [1:2] "coords.x1" "coords.x2"
## .. .. ..$ : chr [1:2] "min" "max"
## ..@ proj4string:Formal class 'CRS' [package "sp"] with 1 slot
## .. .. ..@ projargs: chr "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +units=m +no_defs"John Snow의 콜레라 데이터(사망자) 산점도 (1854년, 127명 사망)
knitr::include_graphics("../pics/Snow-cholera-map-1.jpg", dpi = 96)
plot(CholeraDeaths) ## mdsr 패키지에 있는 John Snow 의 콜레라 데이터 산점도 
공간 데이터 분석 전용 라이브러리를 담고 있는 rgdal 설치. 다운로드한 압축 파일을 data 폴더에 압축을 풀어 둠. CholeraDeaths 와 Pumps의 두 공간정보 데이터 읽어들임/
- GDAL : Geospatial Data Abstraction Library
- OGR : OpenGIS Simple Features Reference Implementation
- dsn : data source name
- SHP : shape file format
- KML : Keyhole Markup Language
- ESRI : Environmental Systems Research Institute, ArcView 개발
library(rgdal) ## rgdal: version: 1.2-8, (SVN revision 663)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 1.11.4, released 2016/01/25
## Path to GDAL shared files: /usr/share/gdal
## Loaded PROJ.4 runtime: Rel. 4.8.0, 6 March 2012, [PJ_VERSION: 480]
## Path to PROJ.4 shared files: (autodetected)
## Linking to sp version: 1.2-5##
## Attaching package: 'rgdal'## The following object is masked from 'package:mosaic':
##
## projectdownload.file("http://rtwilson.com/downloads/SnowGIS_SHP.zip", "../data/SnowGIS_SHP.zip")
dsn <- "../data/SnowGIS_SHP/"
list.files(dsn) ## dsn의 파일 리스트## [1] "Cholera_Deaths.dbf" "Cholera_Deaths.prj"
## [3] "Cholera_Deaths.sbn" "Cholera_Deaths.sbx"
## [5] "Cholera_Deaths.shp" "Cholera_Deaths.shx"
## [7] "OSMap_Grayscale.tfw" "OSMap_Grayscale.tif"
## [9] "OSMap_Grayscale.tif.aux.xml" "OSMap_Grayscale.tif.ovr"
## [11] "OSMap.tfw" "OSMap.tif"
## [13] "Pumps.dbf" "Pumps.prj"
## [15] "Pumps.sbx" "Pumps.shp"
## [17] "Pumps.shx" "README.txt"
## [19] "SnowMap.tfw" "SnowMap.tif"
## [21] "SnowMap.tif.aux.xml" "SnowMap.tif.ovr"ogrListLayers(dsn)## [1] "Cholera_Deaths" "Pumps"
## attr(,"driver")
## [1] "ESRI Shapefile"
## attr(,"nlayers")
## [1] 2ogrInfo(dsn, layer = "Cholera_Deaths")## Source: "../data/SnowGIS_SHP/", layer: "Cholera_Deaths"
## Driver: ESRI Shapefile; number of rows: 250
## Feature type: wkbPoint with 2 dimensions
## Extent: (529160.3 180857.9) - (529655.9 181306.2)
## CRS: +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +units=m +no_defs
## LDID: 87
## Number of fields: 2
## name type length typeName
## 1 Id 0 6 Integer
## 2 Count 0 4 IntegerreadOGR(dsn, "Pumps")## OGR data source with driver: ESRI Shapefile
## Source: "../data/SnowGIS_SHP/", layer: "Pumps"
## with 8 features
## It has 1 fields## coordinates Id
## 1 (529396.5, 181025.1) 0
## 2 (529192.5, 181079.4) 0
## 3 (529183.7, 181193.7) 0
## 4 (529748.9, 180924.2) 0
## 5 (529613.2, 180896.8) 0
## 6 (529453.6, 180826.4) 0
## 7 (529593.7, 180660.5) 0
## 8 (529296.1, 180794.8) 0CholeraDeaths <- readOGR(dsn, layer = "Cholera_Deaths")## OGR data source with driver: ESRI Shapefile
## Source: "../data/SnowGIS_SHP/", layer: "Cholera_Deaths"
## with 250 features
## It has 2 fieldsPumps <- readOGR(dsn, layer = "Pumps")## OGR data source with driver: ESRI Shapefile
## Source: "../data/SnowGIS_SHP/", layer: "Pumps"
## with 8 features
## It has 1 fieldssummary(CholeraDeaths)## Object of class SpatialPointsDataFrame
## Coordinates:
## min max
## coords.x1 529160.3 529655.9
## coords.x2 180857.9 181306.2
## Is projected: TRUE
## proj4string :
## [+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000
## +y_0=-100000 +ellps=airy +units=m +no_defs]
## Number of points: 250
## Data attributes:
## Id Count
## Min. :0 Min. : 1.000
## 1st Qu.:0 1st Qu.: 1.000
## Median :0 Median : 1.000
## Mean :0 Mean : 1.956
## 3rd Qu.:0 3rd Qu.: 2.000
## Max. :0 Max. :15.000summary(Pumps)## Object of class SpatialPointsDataFrame
## Coordinates:
## min max
## coords.x1 529183.7 529748.9
## coords.x2 180660.5 181193.7
## Is projected: TRUE
## proj4string :
## [+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000
## +y_0=-100000 +ellps=airy +units=m +no_defs]
## Number of points: 8
## Data attributes:
## Id
## Min. :0
## 1st Qu.:0
## Median :0
## Mean :0
## 3rd Qu.:0
## Max. :0str(CholeraDeaths@data)## 'data.frame': 250 obs. of 2 variables:
## $ Id : int 0 0 0 0 0 0 0 0 0 0 ...
## $ Count: int 3 2 1 1 4 2 2 2 3 2 ...str(Pumps@data)## 'data.frame': 8 obs. of 1 variable:
## $ Id: int 0 0 0 0 0 0 0 0ggplot
id와 Count 산점도
## coords.x1 coords.x2
## [1,] 529308.7 181031.4
## [2,] 529312.2 181025.2
## [3,] 529314.4 181020.3
## [4,] 529317.4 181014.3
## [5,] 529320.7 181007.9
## [6,] 529336.7 181006.0
## [7,] 529290.1 181024.4
## [8,] 529301.0 181021.2
## [9,] 529285.0 181020.2
## [10,] 529288.4 181031.8
## [11,] 529280.6 181026.6
## [12,] 529264.7 181035.2
## [13,] 529274.9 181008.0
## [14,] 529278.3 181003.0
## [15,] 529281.4 180997.1
## [16,] 529259.6 181010.3
## [17,] 529256.2 181001.8
## [18,] 529309.8 181004.7
## [19,] 529314.1 180997.4
## [20,] 529318.0 180991.3
## [21,] 529251.0 181032.5
## [22,] 529248.0 181049.2
## [23,] 529231.4 181027.7
## [24,] 529236.7 181025.3
## [25,] 529244.0 181028.9
## [26,] 529251.0 181054.8
## [27,] 529247.0 181060.5
## [28,] 529243.6 181066.1
## [29,] 529241.6 181058.4
## [30,] 529217.4 181048.2
## [31,] 529211.5 181056.2
## [32,] 529212.3 181038.3
## [33,] 529191.3 181046.7
## [34,] 529184.2 181041.9
## [35,] 529160.3 181014.2
## [36,] 529228.3 181015.2
## [37,] 529187.3 180983.3
## [38,] 529213.1 180978.5
## [39,] 529214.0 180970.6
## [40,] 529222.8 180978.7
## [41,] 529233.1 180985.7
## [42,] 529221.1 180985.6
## [43,] 529275.9 180965.2
## [44,] 529262.1 180970.8
## [45,] 529299.8 180931.9
## [46,] 529304.2 180925.2
## [47,] 529307.8 180921.2
## [48,] 529302.2 180919.6
## [49,] 529293.0 180929.9
## [50,] 529319.5 180915.9
## [51,] 529330.5 180923.9
## [52,] 529335.4 180928.9
## [53,] 529331.1 180915.7
## [54,] 529345.1 180912.3
## [55,] 529350.2 180915.2
## [56,] 529353.7 180910.2
## [57,] 529377.2 180935.5
## [58,] 529362.9 180947.1
## [59,] 529371.7 180952.2
## [60,] 529326.3 180980.5
## [61,] 529338.1 180959.4
## [62,] 529346.0 180961.4
## [63,] 529349.7 180955.4
## [64,] 529352.8 180950.5
## [65,] 529337.9 180976.0
## [66,] 529335.4 180990.7
## [67,] 529341.9 180993.9
## [68,] 529347.4 180995.6
## [69,] 529352.9 180999.5
## [70,] 529367.8 181007.2
## [71,] 529374.1 181010.4
## [72,] 529383.5 181014.5
## [73,] 529389.1 181017.1
## [74,] 529394.1 181020.4
## [75,] 529409.2 181005.2
## [76,] 529420.9 180989.5
## [77,] 529414.7 180978.8
## [78,] 529407.8 180975.2
## [79,] 529402.6 180971.6
## [80,] 529395.8 180968.2
## [81,] 529387.7 180962.6
## [82,] 529377.1 180956.5
## [83,] 529419.3 180973.1
## [84,] 529400.9 180956.4
## [85,] 529404.3 180948.8
## [86,] 529415.3 180955.6
## [87,] 529421.0 180944.9
## [88,] 529414.9 180930.9
## [89,] 529419.2 180926.4
## [90,] 529431.7 180928.4
## [91,] 529435.7 180922.8
## [92,] 529438.6 180916.7
## [93,] 529443.8 180907.2
## [94,] 529424.2 180916.0
## [95,] 529450.8 180876.2
## [96,] 529456.3 180886.5
## [97,] 529460.5 180880.9
## [98,] 529466.4 180874.2
## [99,] 529508.1 180865.5
## [100,] 529505.9 180884.2
## [101,] 529488.0 180906.0
## [102,] 529487.2 180896.1
## [103,] 529482.7 180901.6
## [104,] 529461.3 180941.2
## [105,] 529454.0 180941.3
## [106,] 529446.5 180951.9
## [107,] 529442.3 180957.8
## [108,] 529436.4 180965.2
## [109,] 529431.8 180970.9
## [110,] 529438.4 180974.8
## [111,] 529434.3 180981.3
## [112,] 529491.1 180963.4
## [113,] 529505.2 180972.5
## [114,] 529539.2 180938.7
## [115,] 529548.7 180901.7
## [116,] 529517.8 180902.8
## [117,] 529506.0 181007.0
## [118,] 529555.2 180997.4
## [119,] 529547.1 181009.3
## [120,] 529552.0 181002.0
## [121,] 529540.8 181017.3
## [122,] 529537.3 181023.3
## [123,] 529532.4 181032.8
## [124,] 529528.9 181038.4
## [125,] 529494.9 181042.4
## [126,] 529519.2 181037.1
## [127,] 529500.8 181048.5
## [128,] 529488.0 181043.0
## [129,] 529507.2 181056.0
## [130,] 529457.0 181013.9
## [131,] 529459.9 181008.7
## [132,] 529461.5 181003.5
## [133,] 529468.3 180994.3
## [134,] 529471.5 180989.7
## [135,] 529478.8 180976.1
## [136,] 529454.7 181018.6
## [137,] 529439.2 181039.2
## [138,] 529432.9 181036.0
## [139,] 529426.9 181032.6
## [140,] 529514.4 180944.8
## [141,] 529528.8 180944.3
## [142,] 529521.7 180940.7
## [143,] 529514.6 180936.8
## [144,] 529618.1 181008.0
## [145,] 529633.0 181018.5
## [146,] 529631.0 181010.0
## [147,] 529655.9 181022.4
## [148,] 529581.2 181052.4
## [149,] 529577.9 181061.4
## [150,] 529570.0 181062.4
## [151,] 529574.8 181050.0
## [152,] 529561.4 181082.5
## [153,] 529565.3 181071.9
## [154,] 529531.8 181068.9
## [155,] 529526.0 181082.9
## [156,] 529509.9 181091.9
## [157,] 529516.9 181056.6
## [158,] 529518.3 181105.1
## [159,] 529565.4 181090.6
## [160,] 529571.1 181075.9
## [161,] 529554.3 181109.6
## [162,] 529548.4 181118.5
## [163,] 529546.0 181124.5
## [164,] 529533.5 181127.5
## [165,] 529553.4 181137.8
## [166,] 529566.3 181144.7
## [167,] 529567.7 181142.5
## [168,] 529573.9 181145.4
## [169,] 529572.7 181147.6
## [170,] 529590.7 181161.4
## [171,] 529583.8 181158.5
## [172,] 529588.9 181155.9
## [173,] 529597.4 181159.6
## [174,] 529595.9 181164.2
## [175,] 529600.7 181166.0
## [176,] 529608.4 181169.9
## [177,] 529615.2 181168.6
## [178,] 529620.2 181171.5
## [179,] 529541.7 181306.2
## [180,] 529489.0 181233.7
## [181,] 529479.7 181228.3
## [182,] 529437.2 181190.5
## [183,] 529382.8 181191.9
## [184,] 529398.0 181200.1
## [185,] 529463.7 181155.2
## [186,] 529468.2 181182.9
## [187,] 529518.8 181153.5
## [188,] 529530.7 181152.2
## [189,] 529520.6 181195.5
## [190,] 529460.5 181143.3
## [191,] 529470.3 181144.7
## [192,] 529464.4 181137.4
## [193,] 529481.0 181125.5
## [194,] 529469.8 181124.7
## [195,] 529485.2 181116.8
## [196,] 529477.1 181112.5
## [197,] 529424.3 181139.0
## [198,] 529417.3 181135.2
## [199,] 529424.1 181159.6
## [200,] 529431.3 181162.6
## [201,] 529437.3 181164.6
## [202,] 529420.8 181167.1
## [203,] 529410.2 181162.7
## [204,] 529374.7 181154.2
## [205,] 529375.7 181128.3
## [206,] 529384.7 181130.6
## [207,] 529386.8 181104.6
## [208,] 529392.3 181114.7
## [209,] 529391.7 181094.0
## [210,] 529432.4 181102.3
## [211,] 529443.7 181108.0
## [212,] 529448.9 181110.3
## [213,] 529455.4 181114.8
## [214,] 529451.4 181103.5
## [215,] 529440.5 181099.1
## [216,] 529434.3 181094.9
## [217,] 529429.4 181093.1
## [218,] 529477.3 181080.5
## [219,] 529484.5 181100.5
## [220,] 529470.6 181077.4
## [221,] 529463.9 181074.1
## [222,] 529454.7 181069.3
## [223,] 529448.5 181065.4
## [224,] 529442.5 181063.1
## [225,] 529427.1 181054.9
## [226,] 529402.3 181040.9
## [227,] 529395.7 181038.1
## [228,] 529389.6 181034.7
## [229,] 529381.9 181031.7
## [230,] 529376.7 181026.6
## [231,] 529357.4 181028.2
## [232,] 529396.3 181082.2
## [233,] 529406.5 181084.7
## [234,] 529347.4 181048.7
## [235,] 529343.4 181054.3
## [236,] 529339.6 181068.2
## [237,] 529281.7 181086.4
## [238,] 529347.1 181071.8
## [239,] 529304.0 181040.6
## [240,] 529293.4 181063.7
## [241,] 529286.6 181055.9
## [242,] 529284.2 181063.4
## [243,] 529296.3 181102.6
## [244,] 529306.3 181179.3
## [245,] 529320.3 181184.8
## [246,] 529362.7 181156.1
## [247,] 529365.2 181176.1
## [248,] 529274.2 180907.3
## [249,] 529299.4 180873.2
## [250,] 529324.8 180857.9
## Google Maps API Terms of Service: http://developers.google.com/maps/terms.## Please cite ggmap if you use it: see citation("ggmap") for details.## Source : https://maps.googleapis.com/maps/api/staticmap?center=John+Snow,+London,+England&zoom=17&size=640x640&scale=2&maptype=roadmap&language=en-EN## Source : https://maps.googleapis.com/maps/api/geocode/json?address=John%20Snow%2C%20London%2C%20England
## Warning: Removed 250 rows containing missing values (geom_point).
## Id Count coords.x1 coords.x2
## 1 0 3 529308.7 181031.4
## 2 0 2 529312.2 181025.2
## 3 0 1 529314.4 181020.3
## 4 0 1 529317.4 181014.3
## 5 0 4 529320.7 181007.9
## 6 0 2 529336.7 181006.0## 'ggmap' chr [1:1280, 1:1280] "#EDD9B0" "#F8ECD4" "#F8ECD4" "#F8ECD4" ...
## - attr(*, "bb")='data.frame': 1 obs. of 4 variables:
## ..$ ll.lat: num 51.5
## ..$ ll.lon: num -0.14
## ..$ ur.lat: num 51.5
## ..$ ur.lon: num -0.133
## - attr(*, "source")= chr "google"
## - attr(*, "maptype")= chr "roadmap"
## - attr(*, "zoom")= num 17## ll.lat ll.lon ur.lat ur.lon
## 1 51.51113 -0.1400205 51.51541 -0.133154maps 패키지의 map()함수
library(maps)##
## Attaching package: 'maps'## The following object is masked from 'package:purrr':
##
## mapmap("world", projection = "mercator", wrap = TRUE)## Warning in map("world", projection = "mercator", wrap = TRUE): projection
## failed for some data
map("world", projection = "cylequalarea", param = 45, wrap = TRUE)
map("state", projection = "lambert",
parameters = c(lat0 = 20, lat1 = 50), wrap = TRUE)
map("state", projection = "albers",
parameters = c(lat0 = 20, lat1 = 50), wrap = TRUE)
- tmerc : Transverse Mercator projectiotn
- CRS : Coordinate Reference System
- EPSG : European Petroleum Survey Group
proj4string(CholeraDeaths) %>%
strwrap() ## [1] "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000"
## [2] "+y_0=-100000 +ellps=airy +units=m +no_defs"CRS("+init=epsg:4326") ## CRS arguments:
## +init=epsg:4326 +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84
## +towgs84=0,0,0CRS("+init=epsg:3857")## CRS arguments:
## +init=epsg:3857 +proj=merc +a=6378137 +b=6378137 +lat_ts=0.0
## +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null
## +no_defsCRS("+init=epsg:27700")## CRS arguments:
## +init=epsg:27700 +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717
## +x_0=400000 +y_0=-100000 +datum=OSGB36 +units=m +no_defs
## +ellps=airy
## +towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894cholera_latlong <- CholeraDeaths %>%
spTransform(CRS("+init=epsg:4326"))
help("spTransform-methods", package = "rgdal")
bbox(cholera_latlong)## min max
## coords.x1 -0.1384685 -0.1313274
## coords.x2 51.5113460 51.5153252ggmap(m) +
geom_point(data = as.data.frame(cholera_latlong),
aes(x = coords.x1,
y = coords.x2,
size = Count))## Warning: Removed 42 rows containing missing values (geom_point).
Snow data 를 지도 상에
CholeraDeaths %>%
proj4string() %>%
showEPSG()## [1] "OGRERR_UNSUPPORTED_SRS"proj4string(CholeraDeaths) <- CRS("+init=epsg:27700")## Warning in `proj4string<-`(`*tmp*`, value = <S4 object of class structure("CRS", package = "sp")>): A new CRS was assigned to an object with an existing CRS:
## +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +units=m +no_defs
## without reprojecting.
## For reprojection, use function spTransformcholera_latlong <- CholeraDeaths %>%
spTransform(CRS("+init=epsg:4326"))
snow <- ggmap(m) +
geom_point(data = as.data.frame(cholera_latlong),
aes(x = coords.x1,
y = coords.x2,
size = Count))
snow## Warning: Removed 3 rows containing missing values (geom_point).
Pumps
pumps <- readOGR(dsn,
layer = "Pumps")## OGR data source with driver: ESRI Shapefile
## Source: "../data/SnowGIS_SHP/", layer: "Pumps"
## with 8 features
## It has 1 fieldsproj4string(pumps) <- CRS("+init=epsg:27700")## Warning in `proj4string<-`(`*tmp*`, value = <S4 object of class structure("CRS", package = "sp")>): A new CRS was assigned to an object with an existing CRS:
## +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +units=m +no_defs
## without reprojecting.
## For reprojection, use function spTransformpumps_latlong <- pumps %>%
spTransform(CRS("+init=epsg:4326"))
snow +
geom_point(data = as.data.frame(pumps_latlong),
aes(x = coords.x1,
y = coords.x2,
size = 3,
colour = "red")) +
# scale_colour_manual(guide = NULL)
guides(colour = "none")## Warning: Removed 3 rows containing missing values (geom_point).## Warning: Removed 2 rows containing missing values (geom_point).