# data wrangling & plotting
library(tidyverse) # dplyr and ggplot2
library(reshape2)
# spatial analyses - for later
library(sf) # working with vectors - polygons, lines, etc
library(raster) # raster
library(USAboundaries)
library(USAboundariesData)
# color scales
library(viridis)
library(ggrepel)
library(grid)
library(gridExtra)
theme_set(theme_bw())# set plotting theme for session
# get MD boundaries - USA boundaries apckage
md.state<-us_states(resolution = 'high', states='louisiana')
glimpse(md.state) # look at spatial object - basically a dataframe
## Rows: 1
## Columns: 13
## $ statefp <chr> "22"
## $ statens <chr> "01629543"
## $ affgeoid <chr> "0400000US22"
## $ geoid <chr> "22"
## $ stusps <chr> "LA"
## $ name <chr> "Louisiana"
## $ lsad <chr> "00"
## $ aland <dbl> 111904912452
## $ awater <dbl> 23746303848
## $ state_name <chr> "Louisiana"
## $ state_abbr <chr> "LA"
## $ jurisdiction_type <chr> "state"
## $ geometry <MULTIPOLYGON [°]> MULTIPOLYGON (((-88.8677 29...
md.cos<-us_counties(resolution = 'high', states='louisiana')
unique(md.cos$name)
## [1] "Ascension" "Beauregard" "Bienville"
## [4] "Bossier" "Catahoula" "Concordia"
## [7] "LaSalle" "Livingston" "Natchitoches"
## [10] "Ouachita" "Rapides" "St. Landry"
## [13] "Vermilion" "Webster" "West Carroll"
## [16] "Sabine" "St. Martin" "St. Mary"
## [19] "Tangipahoa" "Terrebonne" "Vernon"
## [22] "West Feliciana" "East Carroll" "East Feliciana"
## [25] "Franklin" "Iberia" "Lafourche"
## [28] "Madison" "Orleans" "Pointe Coupee"
## [31] "Red River" "Avoyelles" "Caldwell"
## [34] "De Soto" "Jefferson Davis" "Plaquemines"
## [37] "St. Tammany" "Grant" "East Baton Rouge"
## [40] "Acadia" "Caddo" "Evangeline"
## [43] "St. Bernard" "Winn" "St. Helena"
## [46] "Lafayette" "Tensas" "Union"
## [49] "Jefferson" "Morehouse" "Richland"
## [52] "Claiborne" "Jackson" "St. James"
## [55] "Washington" "Iberville" "Allen"
## [58] "Calcasieu" "Cameron" "St. John the Baptist"
## [61] "West Baton Rouge" "Lincoln" "St. Charles"
## [64] "Assumption"
# make simple plot of MD counties using geom_sf
md.map<-ggplot(md.cos)+
geom_sf(fill='blue', color='white')
md.map
md<-md.state%>%st_transform(2248)
md.cos<-md.cos%>%st_transform(2248)
st_crs(md) # units are in us ft
## Coordinate Reference System:
## User input: EPSG:2248
## wkt:
## PROJCS["NAD83 / Maryland (ftUS)",
## GEOGCS["NAD83",
## DATUM["North_American_Datum_1983",
## SPHEROID["GRS 1980",6378137,298.257222101,
## AUTHORITY["EPSG","7019"]],
## TOWGS84[0,0,0,0,0,0,0],
## AUTHORITY["EPSG","6269"]],
## PRIMEM["Greenwich",0,
## AUTHORITY["EPSG","8901"]],
## UNIT["degree",0.0174532925199433,
## AUTHORITY["EPSG","9122"]],
## AUTHORITY["EPSG","4269"]],
## PROJECTION["Lambert_Conformal_Conic_2SP"],
## PARAMETER["standard_parallel_1",39.45],
## PARAMETER["standard_parallel_2",38.3],
## PARAMETER["latitude_of_origin",37.66666666666666],
## PARAMETER["central_meridian",-77],
## PARAMETER["false_easting",1312333.333],
## PARAMETER["false_northing",0],
## UNIT["US survey foot",0.3048006096012192,
## AUTHORITY["EPSG","9003"]],
## AXIS["X",EAST],
## AXIS["Y",NORTH],
## AUTHORITY["EPSG","2248"]]
md.pt<-md%>%st_cast('POINT') # turn MD border into points
glimpse(md.cos)
## Rows: 64
## Columns: 13
## $ statefp <chr> "22", "22", "22", "22", "22", "22", "22", "22", "...
## $ countyfp <chr> "005", "011", "013", "015", "025", "029", "059", ...
## $ countyns <chr> "00558403", "00558436", "00558445", "00558453", "...
## $ affgeoid <chr> "0500000US22005", "0500000US22011", "0500000US220...
## $ geoid <chr> "22005", "22011", "22013", "22015", "22025", "220...
## $ name <chr> "Ascension", "Beauregard", "Bienville", "Bossier"...
## $ lsad <chr> "15", "15", "15", "15", "15", "15", "15", "15", "...
## $ aland <dbl> 751029431, 2997502171, 2101207370, 2176181488, 18...
## $ awater <dbl> 33194722, 21999001, 27945994, 70514387, 81308667,...
## $ state_name <chr> "Louisiana", "Louisiana", "Louisiana", "Louisiana...
## $ state_abbr <chr> "LA", "LA", "LA", "LA", "LA", "LA", "LA", "LA", "...
## $ jurisdiction_type <chr> "state", "state", "state", "state", "state", "sta...
## $ geometry <MULTIPOLYGON [US_survey_foot]> MULTIPOLYGON (((-318063...
# use lon & lat for each centroid to place labels
md.cos$lon<-st_coordinates(st_centroid(md.cos))[,1] # add longitude to sf
md.cos$lat<-st_coordinates(st_centroid(md.cos))[,2] # add latitude to sf
# create mape of MD counties colored by land area (aland)
md.co.map<-ggplot()+
geom_sf(data=md, fill=NA, color='black')+
geom_sf(data=md.cos, aes(fill=aland))+
theme_bw()+
scale_fill_viridis(alpha=.5)+
theme(legend.position='none')
md.co.map
### grid.newpage() clears whatever has been drawn on the device before
grid.newpage()
vp1 <- viewport(x = 0, y = 1,
height = 0.5, width = 0.25,
just = c("left", "top"),
name = "upper left")
pushViewport(vp1)
print(md.co.map, newpage = FALSE)
upViewport(1)
vp2 <- viewport(x = 0, y = 0,
height = 0.5, width = 0.25,
just = c("left", "bottom"),
name = "lower left")
pushViewport(vp2)
print(md.co.map, newpage = FALSE)
upViewport(1)
vp3 <- viewport(x = 0.25, y = 0,
height = 1, width = 0.5,
just = c("left", "bottom"),
name = "centre")
pushViewport(vp3)
print(md.co.map, newpage = FALSE)
upViewport(1)
vp4 <- viewport(x = 1, y = 1,
height = 0.5, width = 0.25,
just = c("right", "top"),
name = "upper right")
pushViewport(vp4)
print(md.co.map, newpage = FALSE)
upViewport(1)
vp5 <- viewport(x = 1, y = 0,
height = 0.5, width = 0.25,
just = c("right", "bottom"),
name = "lower right")
pushViewport(vp5)
print(md.co.map, newpage = FALSE)
upViewport(1)
