library(ggplot2)
library(ggthemes)
library(socviz)
library(maps)
library(mapproj)
library(viridis)DV Lab 7.1
Arizona Hospitals
GET ARIZONA POPULATION DATA BY COUNTY DIRECTLY FROM THE US CENSUS VIA API
head(as_tibble(az_pop))# A tibble: 6 × 6
GEOID NAME variable estimate moe geometry
<chr> <chr> <chr> <dbl> <dbl> <MULTIPOLYGON [°]>
1 04012 La Paz County, Arizona B01003_… 16681 NA (((-114.7312 33.30404, -…
2 04001 Apache County, Arizona B01003_… 66054 NA (((-110.0007 36.99798, -…
3 04013 Maricopa County, Ariz… B01003_… 4430871 NA (((-113.3351 33.55195, -…
4 04015 Mohave County, Arizona B01003_… 214229 NA (((-114.7532 36.08951, -…
5 04005 Coconino County, Ariz… B01003_… 144705 NA (((-113.3541 36.04209, -…
6 04007 Gila County, Arizona B01003_… 53419 NA (((-111.7207 34.1627, -1…
Here’s a look at our Arizona County data.
CREATE A CHOLOPLETH USING LEAFLET FOR ARIZONA COUNTIES
MapPalette <- colorQuantile(palette = "inferno", domain = az_pop$estimate, n = 20)
library(sf)
az_pop %>%
st_transform(crs = "+proj=longlat +datum=WGS84") %>%
leaflet(width = "100%", height = 500) %>%
addProviderTiles(provider = "Esri.WorldPhysical") %>%
addPolygons(popup = ~NAME,
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.7,
color = ~ MapPalette(estimate)) %>%
addLegend("bottomright",
pal = MapPalette,
values = ~ estimate,
title = "Population Percentiles",
opacity = 1) Here are the population percentiles plotted by county using the ‘inferno’ palette from Viridis. Felt like this palette matched the Arizona desert.
library(readxl)
Hospitals <- read_excel("C:/Users/legol/OneDrive - John Brown University/Data Viz/Hospitals.xlsx")
HospitalsAZ <- Hospitals %>% filter(STATE=="AZ")
HospitalsAZ %>% leaflet(width = "100%") %>%
addTiles() %>%
setView(-111.0937, 34.0489, zoom = 6) %>%
addMarkers(lat = ~Y,
lng = ~X,
popup = HospitalsAZ$NAME)Here are Arizona hospitals with popups displaying the hospital name.
MapPalette <- colorQuantile(palette = "inferno", domain = az_pop$estimate, n = 20)
az_pop %>%
st_transform(crs = "+proj=longlat +datum=WGS84") %>%
leaflet(width = "100%", height = 500) %>%
addProviderTiles(provider = "Esri.WorldPhysical") %>%
addPolygons(popup = ~NAME,
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.7,
color = ~ MapPalette(estimate)) %>%
addLegend("bottomright",
pal = MapPalette,
values = ~ estimate,
title = "Population Percentiles",
opacity = 1) %>%
addCircleMarkers(data = HospitalsAZ,
lat = HospitalsAZ$Y,
lng = HospitalsAZ$X,
popup = HospitalsAZ$NAME,
weight = 1,
radius=4,
color = "blue",
opacity = 1)AZ hospitals overlaying the Choropleth. We see a huge concentration in the Phoneiz area which is to be expected.
az_pop %>%
st_transform(crs = "+proj=longlat +datum=WGS84") %>%
leaflet(width = "100%", height = 500) %>%
addProviderTiles(provider = "OpenStreetMap") %>%
addPolygons(popup = ~NAME,
stroke = FALSE,
smoothFactor = 0,
fillOpacity = 0.7,
color = ~MapPalette(estimate)) %>%
addLegend("bottomright",
pal = MapPalette,
values = ~estimate,
title = "Population Percentiles",
opacity = 1) %>%
addCircleMarkers(data = HospitalsAZ,
lat = ~Y,
lng = ~X,
popup = ~paste("Hospital:", NAME, "<br>",
"Beds:", BEDS, "<br>",
"Owner:", OWNER, "<br>",
"City:", CITY, "<br>",
"County:", COUNTY, "<br>",
"Trauma Level:", TRAUMA, "<br>",
"Helipad:", HELIPAD),
weight = 1,
radius = ~sqrt(BEDS) * 0.22,
color = ~OWNER, # Assign colors based on the OWNER variable
fillColor = ~colorFactor(palette = c("green", "blue","red", "cyan"),
domain = c("Government - District/Authority",
"Non-profit",
"Proprietary","Not Available"))(as.character(OWNER)), # Convert OWNER to character
fillOpacity = 0.8, # Increase fill opacity
opacity = 1, # Set opacity to 1
group = "hospital_markers") %>% # Add a group name for the main markers
addLegend("bottomleft",
pal = colorFactor(palette = c("green", "blue","red", "cyan"),
domain = c("Government - District/Authority",
"Non-profit",
"Proprietary",
"Not Available")),
values = HospitalsAZ$OWNER,
title = "Hospital Owner",
opacity = 1)#Esri.WorldStreetMap") %>%
#addProviderTiles(provider = "OpenStreetMap") %>%
#addProviderTiles(provider = "Esri.WorldPhysical") %>%
#addProviderTiles(provider = "Esri.WorldImagery") %>%
#addProviderTiles(provider = "Esri.WorldTopoMap") %>% Added additional dimensions with the bubble size indicating the hospital size by number of beds and the color indicating the type of owner that the hospital has. We see many more Non-profit and Propriety hospitals than government hospitals. We also see more and larger hospitals in high population areas.The popup also tells the user if the hospital has a helipad and what level of trauma it can treat. Chose legend colors that would be very easy to differentiate from the inferno palette. Chose the openstreetmap provider because it had city names and roads.