Install packages
install.packages("mapproj")
install.packages("ggmap")
install.packages("DeducerSpatial")
load packages
require(maps)
## Loading required package: maps
require(ggmap)
## Loading required package: ggmap
## Loading required package: ggplot2
par(mfrow = c(2, 1))
map("usa")
map("county")
map("world", "China")
map.cities(country = "China", capitals = 2)
map("state", "GEORGIA")
data(us.cities)
map.cities(us.cities, country = "GA")
Plot the unemployment in each county
data(unemp)
data(county.fips)
# Plot unemployment by country
colors = c("#F1EEF6", "#D4B9DA", "#C994C7", "#DF65B0", "#DD1C77",
"#980043")
head(unemp)
## fips pop unemp
## 1 1001 23288 9.7
## 2 1003 81706 9.1
## 3 1005 9703 13.4
## 4 1007 8475 12.1
## 5 1009 25306 9.9
## 6 1011 3527 16.4
head(county.fips)
## fips polyname
## 1 1001 alabama,autauga
## 2 1003 alabama,baldwin
## 3 1005 alabama,barbour
## 4 1007 alabama,bibb
## 5 1009 alabama,blount
## 6 1011 alabama,bullock
unemp$colorBuckets <- as.numeric(cut(unemp$unemp, c(0, 2, 4, 6, 8,
10, 100)))
colorsmatched <- unemp$colorBuckets[match(county.fips$fips, unemp$fips)]
map("county", col = colors[colorsmatched], fill = TRUE, resolution = 0,
lty = 0, projection = "polyconic")
## Loading required package: mapproj
# Add border around each State
map("state", col = "white", fill = FALSE, add = TRUE, lty = 1, lwd = 0.2,
projection = "polyconic")
title("unemployment by county, 2009")
leg.txt <- c("<2%", "2-4%", "4-6%", "6-8%", "8-10%", ">10%")
legend("topright", leg.txt, horiz = TRUE, fill = colors)
OpenStreetMap is a new package that accesses raster open street maps from Mapnik, and satellite imagery
from Bing.
Some features:
- Uses multiple map tiles stitched together to create high quality images.
- No files are created or stored on the hard drive.
- Tiles are cached, so downloads occur only when necessary.
- ggplot 0.9.0 integration
More information is available at: http://blog.fellstat.com/?p=130
DeducerSpatial is a package for spatial data analysis which includes the ability to plot and explore open street map and Bing satellite images.
library(UScensus2000)
lat <- c(43.834526782236814,30.334953881988564)
lon <- c(-131.0888671875 ,-107.8857421875)
southwest <- openmap(c(lat[1],lon[1]),c(lat[2],lon[2]),5,'bing')
data(california.tract)
california.tract <- spTransform(california.tract,osm())
plot(southwest,removeMargin=TRUE)
choro_plot(california.tract,dem = slot(california.tract,"data")[,'med.age'], legend.title = 'Median Age',alpha=1)
Maps are extracted from Google Maps, OpenStreetMap, or Stamen Maps server for a map. You can query the Google Maps, OpenStreetMap, or Stamen Maps server for a map at a certain location at a certain spatial zoom.
The geocode function will extract the position (latitude and longtitude) of a location using Google Maps
> geocode('CDC')
lon lat
1 -84.3258 33.7988
> geocode('Baylor University')
lon lat
1 -97.11332 31.54641
> geocode('the white house', messaging = TRUE)
contacting http://maps.googleapis.com/maps/api/geocode/json?address=the+white+house&sensor=false... done.
lon lat
1 -77.03282 38.89521
> geocode(c('baylor university', 'CDC'), output = 'latlona')
lon lat
1 -97.11332 31.54641
2 -84.32580 33.79880
address
1 baylor university, 1311 s 5th st, waco, tx 76706, usa
2 centers for disease control and prevention, 1600 clifton rd ne, atlanta, ga 30329, usa
> geocode(c('harvard university', 'the vatican'), output = 'more')
lon lat type loctype
1 -71.11847 42.37315 point_of_interest approximate
2 12.45813 41.90226 locality approximate
address
1 harvard university housing office, 1350 massachusetts ave, cambridge, ma 02138, usa
2 00120 vatican city
north south east west postal_code country
1 42.38139 42.36490 -71.10246 -71.13447 02138 united states
2 41.90747 41.73199 12.66542 12.44584 00120 vatican city
administrative_area_level_2 administrative_area_level_1 locality
1 middlesex massachusetts cambridge
2 <NA> <NA> vatican city
street streetNo point_of_interest
1 massachusetts ave 1350 harvard university housing office
2 <NA> NA <NA>
Get the geocode for the eiffel tower. Is there a unique map?
geocode('the eiffel tower', output = 'all')
mapdist(from, to,
mode = c("driving", "walking", "bicycling"),
output = c("simple", "all"), messaging = FALSE,
sensor = FALSE, language = "en-EN",
override_limit = FALSE)
Example, how far is it to walk from the CDC to the white house and map the route.
> mapdist('CDC', 'the white house', mode = 'walking')
from to m km miles seconds minutes
1 the white house CDC 1019454 1019.454 633.4887 731359 12189.32
hours
1 203.1553
Google Geocoding API is subject to a query limit of 2,500 geolocation requests per day
geocodeQueryCheck()
Extract location of crimes in houston
violent_crimes <- subset(crime, offense != "auto theft" & offense !=
"theft" & offense != "burglary")
# rank violent crimes
violent_crimes$offense <- factor(violent_crimes$offense, levels = c("robbery",
"aggravated assault", "rape", "murder"))
# restrict to downtown
violent_crimes <- subset(violent_crimes, -95.39681 <= lon & lon <=
-95.34188 & 29.73631 <= lat & lat <= 29.784)
Map these crimes on the map of the city
HoustonMap <- qmap('houston', zoom = 14,color = 'bw', legend = 'topleft')
HoustonMap +geom_point(aes(x = lon, y = lat,
size = offense,colour = offense), data = violent_crimes )
Results of qmap using ggmap of crimes in houston
Plot again but use stats$_$denisty layer
houston <- get_map('houston', zoom = 14)
HoustonMap <- ggmap(houston, extent = 'device', legend = 'topleft')
HoustonMap + stat_density2d(aes(x = lon, y = lat,
fill = ..level.. , alpha = ..level..),size = 2, bins = 4,
data = violent_crimes, geom = 'polygon')
scale_fill_gradient('Violent\nCrime\nDensity') +
scale_alpha(range = c(.4, .75), guide = FALSE) +
guides(fill = guide_colorbar(barwidth = 1.5, barheight = 10))
Results of qmap using ggmap of crimes in houston
Add the colorbar guide to the key
HoustonMap
+ stat_density2d(aes(x = lon, y = lat, fill = ..level.., alpha = ..level..),
size = 2, bins = 4, data = violent_crimes, geom = 'polygon')
+scale_fill_gradient('Violent\nCrime\nDensity')
+scale_alpha(range = c(.4, .75), guide = FALSE)
+guides(fill = guide_colorbar(barwidth = 1.5, barheight = 10))