Making Maps in R
Linda Angulo Lopez
15/10/2020
Reproducible Research Course by Eric C. Anderson for (NOAA/SWFSC)
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.
Using base maps from R’s maps package and also using the ggmap package, with the ggplot2 package.
Read in traditional GIS shapefiles
# Packages you will need:
#################################################################
### Don't bother installing if you already have these:
# install.packages(c("ggplot2", "devtools", "dplyr", "stringr"))
### some standard map packages:
#install.packages(c("maps", "mapdata"))
### the github version of ggmap:
#devtools::install_github("dkahle/ggmap")
#################################################################
# Load up a few of the libraries:
library(ggplot2)
library(ggmap)## Google's Terms of Service: https://cloud.google.com/maps-platform/terms/.## Please cite ggmap if you use it! See citation("ggmap") for details.library(maps)
library(mapdata)
#################################################################Use ggplot2 to plot the maps in the maps package.
# Translate the maps data into a data frame format the ggplot
# can use.
#################################################################
# map_data() function
#################################################################
usa <- map_data("usa")
dim(usa)## [1] 7243 6str(usa)## 'data.frame': 7243 obs. of 6 variables:
## $ long : num -101 -101 -101 -101 -101 ...
## $ lat : num 29.7 29.7 29.7 29.6 29.6 ...
## $ group : num 1 1 1 1 1 1 1 1 1 1 ...
## $ order : int 1 2 3 4 5 6 7 8 9 10 ...
## $ region : chr "main" "main" "main" "main" ...
## $ subregion: chr NA NA NA NA ...w2hr <- map_data("world2Hires")
str(w2hr)## 'data.frame': 2274539 obs. of 6 variables:
## $ long : num 227 227 227 227 227 ...
## $ lat : num 58.4 58.4 58.4 58.4 58.4 ...
## $ group : num 1 1 1 1 1 1 1 1 1 1 ...
## $ order : int 1 2 3 4 5 6 7 8 9 10 ...
## $ region : chr "Canada" "Canada" "Canada" "Canada" ...
## $ subregion: chr NA NA NA NA ...# ggplot2's functions can take a group argument which controls
# (amongst other things) whether adjacent points should be
# connected by lines. If they are in the same group, then they
# get connected, but if they are in different groups then they
# don't.
#################################################################- Plot a USA Map
#################################################################
# geom_polygon()
# coord_fixed()
#################################################################
#map the group aesthetic to the group column:
usa <- map_data("usa") # we already did this, but we can do it again
ggplot() + geom_polygon(data = usa,
aes(x=long, y = lat, group = group)) +
coord_fixed(1.3)
#################################################################- Style the USA Map
#################################################################
# aes((x=long, y = lat, group = group)
#################################################################
ggplot() +
geom_polygon(data = usa, aes(x=long, y = lat, group = group), fill = NA, color = "red") +
coord_fixed(1.3)
#################################################################- Adding points to the USA Map
#################################################################
# data.frame( c(long(x1,x2),lat(y1,y2)), c(name1,name2),
# stringsAsFactors = FALSE)
#################################################################
gg1 <- ggplot() +
geom_polygon(data = usa, aes(x=long, y = lat, group = group),
fill = "violet", color = "blue") +
coord_fixed(1.3)
gg1
labs <- data.frame(
long = c(-122.064873, -122.306417),
lat = c(36.951968, 47.644855),
names = c("SWFSC-FED", "NWFSC"),
stringsAsFactors = FALSE
)
gg1 +
geom_point(data = labs, aes(x = long, y = lat), color = "black", size = 5) +
geom_point(data = labs, aes(x = long, y = lat), color = "yellow", size = 4)
#################################################################- Plot USA States on USA Map
#################################################################
# map_data("state")
# fill = region
#################################################################
states <- map_data("state")
dim(states)## [1] 15537 6str(states)## 'data.frame': 15537 obs. of 6 variables:
## $ long : num -87.5 -87.5 -87.5 -87.5 -87.6 ...
## $ lat : num 30.4 30.4 30.4 30.3 30.3 ...
## $ group : num 1 1 1 1 1 1 1 1 1 1 ...
## $ order : int 1 2 3 4 5 6 7 8 9 10 ...
## $ region : chr "alabama" "alabama" "alabama" "alabama" ...
## $ subregion: chr NA NA NA NA ...ggplot(data = states) +
geom_polygon(aes(x = long, y = lat, fill = region, group = group),
color = "white") +
coord_fixed(1.3) +
guides(fill=FALSE) # do this to leave off the color legend
#################################################################- Plot a subset of USA States Map
#################################################################
# subset()
# aes(group=group)
# coord_fixed()
#################################################################
west_coast <- subset(states, region %in% c("california", "oregon", "washington"))
ggplot(data = west_coast) +
geom_polygon(aes(x = long, y = lat, group=group),
fill = "palegreen", color = "black") +
coord_fixed(1.3)
#################################################################- Plot the Counties of the State of California
#################################################################
# subset(states, region == "california")
# map_data("county")
# theme_nothing()
#################################################################
ca_df <- subset(states, region == "california")
head(ca_df)## long lat group order region subregion
## 667 -120.0060 42.00927 4 667 california <NA>
## 668 -120.0060 41.20139 4 668 california <NA>
## 669 -120.0060 39.70024 4 669 california <NA>
## 670 -119.9946 39.44241 4 670 california <NA>
## 671 -120.0060 39.31636 4 671 california <NA>
## 672 -120.0060 39.16166 4 672 california <NA>counties <- map_data("county")
ca_county <- subset(counties, region == "california")
head(ca_county)## long lat group order region subregion
## 6965 -121.4785 37.48290 157 6965 california alameda
## 6966 -121.5129 37.48290 157 6966 california alameda
## 6967 -121.8853 37.48290 157 6967 california alameda
## 6968 -121.8968 37.46571 157 6968 california alameda
## 6969 -121.9254 37.45998 157 6969 california alameda
## 6970 -121.9483 37.47717 157 6970 california alamedaca_base <- ggplot(data = ca_df, mapping = aes(x = long, y = lat,
group = group)) +
coord_fixed(1.3) +
geom_polygon(color = "black", fill = "gray")
# grey map of California
ca_base + theme_nothing()
# white county border back on top
ca_base + theme_nothing() +
geom_polygon(data = ca_county, fill = NA, color = "white") +
geom_polygon(color = "black", fill = NA) 
#################################################################- Add data to the Map
#################################################################
## installed.packages("readr")
## devtools::install_github("tidyverse/readr")
### install.packages("tidyverse")
### devtools::install_github("tidyverse")
#################################################################
# library(stringr)
# library(dplyr)
# library(tidyverse)
# readlines()
# as.data.frame(stringsAsFactors = FALSE)
#################################################################
library(stringr)
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(tidyverse)## -- Attaching packages --------------------------------------- tidyverse 1.3.0 --## v tibble 3.0.4 v purrr 0.3.4
## v tidyr 1.1.2 v forcats 0.5.0
## v readr 1.4.0.9000## Warning: package 'tibble' was built under R version 4.0.3## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
## x purrr::map() masks maps::map()#import data set and assign column name
datF <- read.csv("C:/Users/angul/OneDrive/R/ReproducibleResearch/cc-est2019-agesex-06_subregion.csv", stringsAsFactors = FALSE)
names(datF)## [1] "SUMLEV" "STATE" "COUNTY" "subregion"
## [5] "STNAME" "CTYNAME" "YEAR" "POPESTIMATE"
## [9] "POPEST_MALE" "POPEST_FEM" "UNDER5_TOT" "UNDER5_MALE"
## [13] "UNDER5_FEM" "AGE513_TOT" "AGE513_MALE" "AGE513_FEM"
## [17] "AGE1417_TOT" "AGE1417_MALE" "AGE1417_FEM" "AGE1824_TOT"
## [21] "AGE1824_MALE" "AGE1824_FEM" "AGE16PLUS_TOT" "AGE16PLUS_MALE"
## [25] "AGE16PLUS_FEM" "AGE18PLUS_TOT" "AGE18PLUS_MALE" "AGE18PLUS_FEM"
## [29] "AGE1544_TOT" "AGE1544_MALE" "AGE1544_FEM" "AGE2544_TOT"
## [33] "AGE2544_MALE" "AGE2544_FEM" "AGE4564_TOT" "AGE4564_MALE"
## [37] "AGE4564_FEM" "AGE65PLUS_TOT" "AGE65PLUS_MALE" "AGE65PLUS_FEM"
## [41] "AGE04_TOT" "AGE04_MALE" "AGE04_FEM" "AGE59_TOT"
## [45] "AGE59_MALE" "AGE59_FEM" "AGE1014_TOT" "AGE1014_MALE"
## [49] "AGE1014_FEM" "AGE1519_TOT" "AGE1519_MALE" "AGE1519_FEM"
## [53] "AGE2024_TOT" "AGE2024_MALE" "AGE2024_FEM" "AGE2529_TOT"
## [57] "AGE2529_MALE" "AGE2529_FEM" "AGE3034_TOT" "AGE3034_MALE"
## [61] "AGE3034_FEM" "AGE3539_TOT" "AGE3539_MALE" "AGE3539_FEM"
## [65] "AGE4044_TOT" "AGE4044_MALE" "AGE4044_FEM" "AGE4549_TOT"
## [69] "AGE4549_MALE" "AGE4549_FEM" "AGE5054_TOT" "AGE5054_MALE"
## [73] "AGE5054_FEM" "AGE5559_TOT" "AGE5559_MALE" "AGE5559_FEM"
## [77] "AGE6064_TOT" "AGE6064_MALE" "AGE6064_FEM" "AGE6569_TOT"
## [81] "AGE6569_MALE" "AGE6569_FEM" "AGE7074_TOT" "AGE7074_MALE"
## [85] "AGE7074_FEM" "AGE7579_TOT" "AGE7579_MALE" "AGE7579_FEM"
## [89] "AGE8084_TOT" "AGE8084_MALE" "AGE8084_FEM" "AGE85PLUS_TOT"
## [93] "AGE85PLUS_MALE" "AGE85PLUS_FEM" "MEDIAN_AGE_TOT" "MEDIAN_AGE_MALE"
## [97] "MEDIAN_AGE_FEM"CAgender <-subset(datF, select=c("subregion", "POPESTIMATE",
"POPEST_MALE", "POPEST_FEM", "YEAR", "STNAME" ))
table(CAgender$subregion)##
## fresno mariposa sierra yolo alameda
## 12 12 12 12 12
## alpine amador butte calaveras colusa
## 12 12 12 12 12
## contra costa del norte el dorado glenn humboldt
## 12 12 12 12 12
## imperial inyo kern kings lake
## 12 12 12 12 12
## lassen los angeles madera marin mendocino
## 12 12 12 12 12
## merced modoc mono monterey napa
## 12 12 12 12 12
## nevada orange placer plumas riverside
## 12 12 12 12 12
## sacramento san benito san bernardino san diego san francisco
## 12 12 12 12 12
## san joaquin san luis obispo san mateo santa barbara santa clara
## 12 12 12 12 12
## santa cruz shasta siskiyou solano sonoma
## 12 12 12 12 12
## stanislaus sutter tehama trinity tulare
## 12 12 12 12 12
## tuolumne ventura yuba
## 12 12 12table(ca_county$subregion)##
## alameda alpine amador butte calaveras
## 36 36 33 72 35
## colusa contra costa del norte el dorado fresno
## 41 50 45 52 101
## glenn humboldt imperial inyo kern
## 35 66 25 63 37
## kings lake lassen los angeles madera
## 25 67 25 51 68
## marin mariposa mendocino merced modoc
## 62 45 73 43 7
## mono monterey napa nevada orange
## 54 90 49 30 37
## placer plumas riverside sacramento san benito
## 46 66 56 60 54
## san bernardino san diego san francisco san joaquin san luis obispo
## 58 56 18 48 62
## san mateo santa barbara santa clara santa cruz shasta
## 51 50 71 42 65
## sierra siskiyou solano sonoma stanislaus
## 34 72 43 60 43
## sutter tehama trinity tulare tuolumne
## 43 68 94 48 70
## ventura yolo yuba
## 34 62 50CA_popG <- inner_join(ca_county, CAgender, by = "subregion")
head(CA_popG)## long lat group order region subregion POPESTIMATE POPEST_MALE
## 1 -121.4785 37.4829 157 6965 california alameda 1510271 740573
## 2 -121.4785 37.4829 157 6965 california alameda 1510258 740568
## 3 -121.4785 37.4829 157 6965 california alameda 1512986 741756
## 4 -121.4785 37.4829 157 6965 california alameda 1530915 750107
## 5 -121.4785 37.4829 157 6965 california alameda 1553764 760739
## 6 -121.4785 37.4829 157 6965 california alameda 1579593 773696
## POPEST_FEM YEAR STNAME
## 1 769698 1 California
## 2 769690 2 California
## 3 771230 3 California
## 4 780808 4 California
## 5 793025 5 California
## 6 805897 6 California#################################################################- Plot a the percent of women in the State of California, by county
#################################################################
# geom_polygon()
#################################################################
#add a column of PCwomen:
CA_popG$PCwomen <- (CA_popG$POPEST_FEM/CA_popG$POPESTIMATE)*100
head(CA_popG)## long lat group order region subregion POPESTIMATE POPEST_MALE
## 1 -121.4785 37.4829 157 6965 california alameda 1510271 740573
## 2 -121.4785 37.4829 157 6965 california alameda 1510258 740568
## 3 -121.4785 37.4829 157 6965 california alameda 1512986 741756
## 4 -121.4785 37.4829 157 6965 california alameda 1530915 750107
## 5 -121.4785 37.4829 157 6965 california alameda 1553764 760739
## 6 -121.4785 37.4829 157 6965 california alameda 1579593 773696
## POPEST_FEM YEAR STNAME PCwomen
## 1 769698 1 California 50.96423
## 2 769690 2 California 50.96414
## 3 771230 3 California 50.97403
## 4 780808 4 California 51.00270
## 5 793025 5 California 51.03896
## 6 805897 6 California 51.01928summary(CA_popG)## long lat group order
## Min. :-124.4 Min. :32.54 Min. :157.0 Min. :6965
## 1st Qu.:-122.5 1st Qu.:36.40 1st Qu.:173.0 1st Qu.:7793
## Median :-121.4 Median :38.03 Median :188.0 Median :8541
## Mean :-120.9 Mean :37.79 Mean :186.8 Mean :8495
## 3rd Qu.:-119.9 3rd Qu.:39.28 3rd Qu.:199.0 3rd Qu.:9123
## Max. :-114.1 Max. :42.02 Max. :211.0 Max. :9849
## region subregion POPESTIMATE POPEST_MALE
## Length:27396 Length:27396 Min. : 8800 Min. : 4396
## Class :character Class :character 1st Qu.: 54221 1st Qu.: 28002
## Mode :character Mode :character Median : 256721 Median : 126617
## Mean : 787657 Mean : 391346
## 3rd Qu.: 743296 3rd Qu.: 368198
## Max. :10105708 Max. :4980981
## POPEST_FEM YEAR STNAME PCwomen
## Min. : 4393 Min. : 1.00 Length:27396 Min. :35.52
## 1st Qu.: 25974 1st Qu.: 3.75 Class :character 1st Qu.:49.28
## Median : 129849 Median : 6.50 Mode :character Median :50.02
## Mean : 396312 Mean : 6.50 Mean :49.59
## 3rd Qu.: 376048 3rd Qu.: 9.25 3rd Qu.:50.52
## Max. :5125335 Max. :12.00 Max. :51.25#PCwomen
#Min. :35.52
#1st Qu.:49.28
#Median :50.02
#Mean :49.59
#3rd Qu.:50.52
#Max. :51.25
# prepare to drop the axes and ticks but leave the guides and legends
# We can't just throw down a theme_nothing()!
ditch_the_axes <- theme(
axis.text = element_blank(),
axis.line = element_blank(),
axis.ticks = element_blank(),
panel.border = element_blank(),
panel.grid = element_blank(),
axis.title = element_blank()
)
CA_women <- ca_base +
geom_polygon(data = CA_popG, aes(fill = PCwomen), color = "white") +
geom_polygon(color = "black", fill = NA) +
theme_bw() +
ditch_the_axes
CA_women
################################################################# CAnW <-(CA_women + coord_fixed(xlim = c(-130.0, -110.0), ylim = c(42.5, 38.5), ratio = 1.3))## Coordinate system already present. Adding new coordinate system, which will replace the existing one. CAnW 
#################################################################ggmap
#################################################################
# install.packages("viridis")
# library("viridis")
# viridis(n), magma(n), inferno(n) , plasma(n)
#################################################################
library("viridis")## Loading required package: viridisLiteCAwc <- CA_women +
scale_fill_gradientn(colours = viridis(10),
breaks = c(34, 36,38, 40,42, 44, 46, 48, 50, 52) ,
trans = "log10")
CAwc 
##################################################################################################################################
#install.packages("RColorBrewer")
#library(RColorBrewer)
#display.brewer.all()
# Base R / rainbow(n), heat.colors(n), terrain.colors(n),
#topo.colors(n), cm.colors(n).
#################################################################
library(RColorBrewer)
CAwc2 <- CA_women +
scale_fill_gradientn(colours = rev(topo.colors(10)),
breaks = c(34, 36,38, 40,42, 44, 46, 48, 50, 52) ,
trans = "log10")
CAwc2 
#################################################################to use ggmap you will need to get a Google API
#################################################################
#install.packages("leaflet")
#library(leaflet)
#################################################################
library(leaflet)
leaflet() %>%
addTiles() %>%
addMarkers(lat = 48.8583701, lng = 2.2922926,
popup = "Eiffel Tower")register_google(key = “[your key]”, write = TRUE)
#################################################################
#register_google(key = "[your key]", write = TRUE)
#get_googlemap()
#geocode()
#################################################################
library("ggmap")
us <- c(left = -125, bottom = 25.75, right = -67, top = 49)
get_stamenmap(us, zoom = 5, maptype = "toner-lite") %>% ggmap() ## Map tiles by Stamen Design, under CC BY 3.0. Data by OpenStreetMap, under ODbL.
Three ggmap examples:
- Named “sampling” points on the Sisquoc River from the “Sisquoctober Adventure”
- A GPS track from a short bike ride in Wilder Ranch.
- Fish sampling locations from the coded wire tag data base.
Sisquoctober Adventure
#################################################################
#ggmap()
#make_bbox()
#################################################################
sisquoc <- read.table("C:/Users/angul/OneDrive/R/ReproducibleResearch/sisquoc-points.txt", sep = "\t", header = TRUE)
sbbox <- make_bbox(lon = sisquoc$lon, lat = sisquoc$lat, f = .1)
sbbox## left bottom right top
## -119.76198 34.75111 -119.74201 34.75507# left bottom right top
# -119.76198 34.75111 -119.74201 34.75507
# compute the mean lat and lon
ll_means <- sapply(sisquoc[2:3], mean)
sq_map2 <- get_map(location = ll_means, maptype = "satellite", source = "google", zoom = 15)## Source : https://maps.googleapis.com/maps/api/staticmap?center=34.753117,-119.751324&zoom=15&size=640x640&scale=2&maptype=satellite&language=en-EN&key=xxx#> Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=34.753117,-119.751324&zoom=15&size=640x640&scale=2&maptype=satellite&language=en-EN&sensor=false
ggmap(sq_map2) +
geom_point(data = sisquoc, color = "red", size = 4) +
geom_text(data = sisquoc, aes(label = paste(" ", as.character(name), sep="")), angle = 60, hjust = 0, color = "yellow")
##################################################################################################################################
sq_map3 <- get_map(location = ll_means, maptype = "terrain", source = "google", zoom = 15)## Source : https://maps.googleapis.com/maps/api/staticmap?center=34.753117,-119.751324&zoom=15&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx#> Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=34.753117,-119.751324&zoom=15&size=640x640&scale=2&maptype=terrain&language=en-EN&sensor=false
ggmap(sq_map3) +
geom_point(data = sisquoc, color = "red", size = 4) +
geom_text(data = sisquoc, aes(label = paste(" ", as.character(name), sep="")), angle = 60, hjust = 0, color = "yellow")
#################################################################Plot a bike route
#################################################################
bike <- read.csv("C:/Users/angul/OneDrive/R/ReproducibleResearch/bike-ride.csv")
head(bike)## lon lat elevation time
## 1 -122.0646 36.95144 15.8 2011-12-08T19:37:56Z
## 2 -122.0646 36.95191 15.5 2011-12-08T19:37:59Z
## 3 -122.0645 36.95201 15.4 2011-12-08T19:38:04Z
## 4 -122.0645 36.95218 15.5 2011-12-08T19:38:07Z
## 5 -122.0643 36.95224 15.7 2011-12-08T19:38:10Z
## 6 -122.0642 36.95233 15.8 2011-12-08T19:38:13Zbikemap1 <- get_map(location = c(-122.080954, 36.971709), maptype = "terrain", source = "google", zoom = 14)## Source : https://maps.googleapis.com/maps/api/staticmap?center=36.971709,-122.080954&zoom=14&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx#> Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=36.971709,-122.080954&zoom=14&size=640x640&scale=2&maptype=terrain&language=en-EN&sensor=false
ggmap(bikemap1) +
geom_path(data = bike, aes(color = elevation), size = 3, lineend = "round") +
scale_color_gradientn(colours = rainbow(7), breaks = seq(25, 200, by = 25))
#################################################################Exploring Recovery and Catch-Sample Data
#################################################################
bc <- readRDS("C:/Users/angul/OneDrive/R/ReproducibleResearch/bc_sites.rds")
# look at some of it:
bc %>% select(state_or_province:sub_location, longitude, latitude)## # A tibble: 1,113 x 9
## state_or_provin~ water_type sector region area location sub_location
## <chr> <chr> <chr> <chr> <chr> <chr> <chr>
## 1 2 M S 22 "016~ "THOR I~ " 01"
## 2 2 M N 26 "012~ "MITC B~ " 18"
## 3 2 M S 22 "015~ "HARW I~ " 02"
## 4 2 M N 26 "006~ "HOPK P~ " 01"
## 5 2 M S 23 "017~ "TENT I~ " 06"
## 6 2 M S 28 "23A~ "NAHM B~ " 02"
## 7 2 M N 26 "006~ "GIL IS~ " 06"
## 8 2 M S 27 "024~ "CLEL I~ " 06"
## 9 2 M S 27 "23B~ "SAND I~ " 04"
## 10 2 M N 26 "012~ "DUVA I~ " 16"
## # ... with 1,103 more rows, and 2 more variables: longitude <dbl>,
## # latitude <dbl>bc %>% group_by(sector, region, area) %>% tally()## # A tibble: 42 x 4
## # Groups: sector, region [11]
## sector region area n
## <chr> <chr> <chr> <int>
## 1 " " 48 "008 " 1
## 2 " " 48 "028 " 1
## 3 " " 48 "311 " 1
## 4 "N" 25 "001 " 33
## 5 "N" 25 "003 " 15
## 6 "N" 25 "004 " 44
## 7 "N" 25 "02E " 2
## 8 "N" 25 "02W " 34
## 9 "N" 26 "006 " 28
## 10 "N" 26 "007 " 23
## # ... with 32 more rows# compute the bounding box
bc_bbox <- make_bbox(lat = latitude, lon = longitude, data = bc)
bc_bbox## left bottom right top
## -133.63297 47.92497 -122.33652 55.80833# grab the maps from google
bc_big <- get_map(location = bc_bbox, source = "google", maptype = "terrain")## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=51.86665,-127.98475&zoom=6&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx#> Warning: bounding box given to google - spatial extent only approximate.
#> converting bounding box to center/zoom specification. (experimental)
#> Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=51.86665,-127.98475&zoom=6&size=640x640&scale=2&maptype=terrain&language=en-EN&sensor=false
# plot the points and color them by sector
ggmap(bc_big) +
geom_point(data = bc, mapping = aes(x = longitude, y = latitude, color = sector))
#Coloring it by region
ggmap(bc_big) +
geom_point(data = bc, mapping = aes(x = longitude, y = latitude, color = region))
#################################################################Zooming in on each region and coloring by area
#################################################################
region_plot <- function(MyRegion) {
tmp <- bc %>% filter(region == MyRegion)
bbox <- make_bbox(lon = longitude, lat = latitude, data = tmp)
mymap <- get_map(location = bbox, source = "google", maptype = "terrain")
# now we want to count up how many areas there are
NumAreas <- tmp %>% summarise(n_distinct(area))
NumPoints <- nrow(tmp)
the_map <- ggmap(mymap) +
geom_point(data = tmp, mapping = aes(x = longitude, y = latitude), size = 4, color = "black") +
geom_point(data = tmp, mapping = aes(x = longitude, y = latitude, color = area), size = 3) +
ggtitle(
paste("BC Region: ", MyRegion, " with ", NumPoints, " locations in ", NumAreas, " area(s)", sep = "")
)
ggsave(paste("bc_region", MyRegion, ".pdf", sep = ""), the_map, width = 9, height = 9)
}
dump <- lapply(unique(bc$region), region_plot)## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=49.925,-124.85835&zoom=9&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 18 rows containing missing values (geom_point).
## Warning: Removed 18 rows containing missing values (geom_point).## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=52.225,-127.76665&zoom=8&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 93 rows containing missing values (geom_point).## Warning: Removed 93 rows containing missing values (geom_point).## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=48.96665,-123.51665&zoom=9&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 2 rows containing missing values (geom_point).## Warning: Removed 2 rows containing missing values (geom_point).## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=49.1153,-124.69525&zoom=11&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=49.65445,-127.11665&zoom=8&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 89 rows containing missing values (geom_point).## Warning: Removed 89 rows containing missing values (geom_point).## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=48.475,-124&zoom=10&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 62 rows containing missing values (geom_point).## Warning: Removed 62 rows containing missing values (geom_point).## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=54.10835,-131.3514&zoom=8&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 11 rows containing missing values (geom_point).## Warning: Removed 11 rows containing missing values (geom_point).## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=50.3,-125.43335&zoom=10&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 16 rows containing missing values (geom_point).## Warning: Removed 16 rows containing missing values (geom_point).## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=49.78335,-124.5583&zoom=10&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 31 rows containing missing values (geom_point).## Warning: Removed 31 rows containing missing values (geom_point).## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=50.56665,-124.98605&zoom=8&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 3 rows containing missing values (geom_point).## Warning: Removed 3 rows containing missing values (geom_point).## Bounding box given to Google - spatial extent only approximate.## Source : https://maps.googleapis.com/maps/api/staticmap?center=49.25,-124.925&zoom=14&size=640x640&scale=2&maptype=terrain&language=en-EN&key=xxx## Warning: Removed 3 rows containing missing values (geom_point).
## Warning: Removed 3 rows containing missing values (geom_point).#################################################################