OBJECTIVE

Create a sereies of 4 maps that highlight the spatial differences in hospital service coverage for th state of PA. This is done using leaflet package.

library(foreign)
## Warning: package 'foreign' was built under R version 3.4.4
library(maptools)
## Loading required package: sp
## Checking rgeos availability: FALSE
##      Note: when rgeos is not available, polygon geometry     computations in maptools depend on gpclib,
##      which has a restricted licence. It is disabled by default;
##      to enable gpclib, type gpclibPermit()
library(leaflet)
## Warning: package 'leaflet' was built under R version 3.4.4
library(rgdal)
## rgdal: version: 1.2-16, (SVN revision 701)
##  Geospatial Data Abstraction Library extensions to R successfully loaded
##  Loaded GDAL runtime: GDAL 2.1.3, released 2017/20/01
##  Path to GDAL shared files: /Library/Frameworks/R.framework/Versions/3.4/Resources/library/rgdal/gdal
##  GDAL binary built with GEOS: FALSE 
##  Loaded PROJ.4 runtime: Rel. 4.9.3, 15 August 2016, [PJ_VERSION: 493]
##  Path to PROJ.4 shared files: /Library/Frameworks/R.framework/Versions/3.4/Resources/library/rgdal/proj
##  Linking to sp version: 1.2-5
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, union
PA_dat <- read.dbf("/Users/lekhana/Google Drive/_Landing Zone/HU/2018_05to08_DataVisualization/Data Viz in R/GIS data/hospitals (1)/pennsylv.dbf")
names(PA_dat)
##   [1] "acc_trauma" "air_amb"    "als"        "arc_street" "arc_zone"  
##   [6] "bas_ls"     "bassinets"  "bb_id"      "bc_beds"    "bc_sus_bed"
##  [11] "beds_sus"   "birthing_r" "bone_marro" "burn_car"   "burn_care" 
##  [16] "card_beds"  "card_surge" "card_sus_b" "cardiac"    "cardiac_ca"
##  [21] "cardio_reh" "chemo"      "city"       "clin_lab"   "clin_psych"
##  [26] "county"     "countyname" "ct_scan"    "cty_key"    "cystoscopi"
##  [31] "deliv_rms"  "dental"     "detox_alc_" "diag_radio" "diag_xray" 
##  [36] "doh_hosp"   "doh_phone"  "emer_dept"  "endoscopie" "fac_id"    
##  [41] "facility"   "flu_old"    "fred_con_1" "fred_conta" "fred_email"
##  [46] "fred_fax"   "fred_hosp"  "fred_pager" "fred_phone" "gamma_knif"
##  [51] "gen_outpat" "gene_counc" "heart_tran" "helipad"    "hemodial_c"
##  [56] "hemodial_m" "hosp_id"    "hospice"    "hyper_cham" "icu"       
##  [61] "icu_beds"   "icu_sus_be" "inpat_flu_" "inpat_pneu" "kidney_tra"
##  [66] "labor_rms"  "lic_beds"   "lic_dent"   "lic_dos"    "lic_mds"   
##  [71] "lic_pod"    "linear_acc" "lithotrips" "liver_tran" "loc_method"
##  [76] "ltc"        "mcd"        "mcd_key"    "mcd_name"   "medical"   
##  [81] "mob_ccu"    "mob_icu"    "mri"        "ms1"        "neo2_beds" 
##  [86] "neo2_sus_b" "neo3_beds"  "neo3_sus_b" "neuro_surg" "neurology" 
##  [91] "obs_gyn"    "occ_ther"   "optometry"  "organ_bank" "ped_trauma"
##  [96] "pediatric"  "pet"        "pharmacy"   "phys_med"   "phys_ther" 
## [101] "podiatry"   "providerid" "psych"      "psych_inpa" "reg_trauma"
## [106] "resp_ther"  "so_flu_65u" "social_wor" "speech_pat" "street"    
## [111] "surgical"   "surgical_s" "thera_radi" "typ_org"    "typ_serv"  
## [116] "ultrasound" "x"          "y"          "zip"
PAdistricts<- readOGR(dsn = path.expand("/Users/lekhana/Google Drive/_Landing Zone/HU/2018_05to08_DataVisualization/Data Viz in R/GIS data/cb_2016_42_unsd_500k"),layer="cb_2016_42_unsd_500k")
## OGR data source with driver: ESRI Shapefile 
## Source: "/Users/lekhana/Google Drive/_Landing Zone/HU/2018_05to08_DataVisualization/Data Viz in R/GIS data/cb_2016_42_unsd_500k", layer: "cb_2016_42_unsd_500k"
## with 501 features
## It has 8 fields
## Integer64 fields read as strings:  ALAND AWATER
PAcounty <- readOGR(dsn = path.expand("/Users/lekhana/Google Drive/_Landing Zone/HU/2018_05to08_DataVisualization/Data Viz in R/GIS data/cb_2016_us_county_500k"),layer="cb_2016_us_county_500k")
## OGR data source with driver: ESRI Shapefile 
## Source: "/Users/lekhana/Google Drive/_Landing Zone/HU/2018_05to08_DataVisualization/Data Viz in R/GIS data/cb_2016_us_county_500k", layer: "cb_2016_us_county_500k"
## with 3233 features
## It has 9 fields
## Integer64 fields read as strings:  ALAND AWATER

CREATING A POLYGON LAYER

Using the district shapefile to create a polygon layer over a leaflet map.

Districts<-leaflet(PAdistricts) %>% addTiles() %>%
        addPolygons(weight=.75,color="green",fillOpacity = .2) 
Districts
leaflet(PA_dat) %>% 
  addTiles() %>%
  addCircles(lng = ~y, lat = ~x, weight = 1,
    radius = ~beds_sus*40, popup = ~county)

NUMBER OF CARDIAC ARRESTS IN PHILADELPHIA

Using the hospitals data set, and adding awesomeIcons feature to the leaflet maps creating to plot the number of hospitals in Philadelphia which treat for Cardiac problems.

icons <- awesomeIcons(icon = "blue",
                      iconColor = "black",
                      library = "glyphicon",
                      markerColor = "purple")

noofcardiac <- subset(PA_dat,cardiac == "Y")
PA_cardiac <- subset(noofcardiac, county == "Philadelphia")
leaflet(PA_cardiac, width = "100%") %>%
        addTiles() %>%addProviderTiles(providers$Esri.WorldPhysical) %>%
        addProviderTiles(providers$Esri.WorldPhysical,options = providerTileOptions(opacity = 0.50)) %>%
       addProviderTiles(providers$Esri.WorldPhysical) %>%
       addAwesomeMarkers(~x, ~y, popup = ~as.character(facility), icon = icons)

##NUMBER OF CHEMO HOSPITALS IN PITTSBURGH Here I am using the greenleafIcon typical of leaflet package to check th number os hospitals that treat for chemo in Pittsburgh.

greenLeafIcon <- makeIcon(
  iconUrl = "http://leafletjs.com/examples/custom-icons/leaf-green.png",
  iconWidth = 38, iconHeight = 70,
  iconAnchorX = 22, iconAnchorY = 94,
  shadowUrl = "http://leafletjs.com/examples/custom-icons/leaf-shadow.png",
  shadowWidth = 50, shadowHeight = 64,
  shadowAnchorX = 4, shadowAnchorY = 62
)
library(leaflet)
noofchemo <- subset(PA_dat,chemo == "Y")
PA_chemo <- subset(noofchemo, county == "Allegheny")
leaflet(PA_chemo, width = "100%") %>%
        addTiles() %>%addProviderTiles(providers$HERE.normalDay) %>%
        addProviderTiles(providers$HERE.normalDay,options = providerTileOptions(opacity = 0.50)) %>%
       addProviderTiles(providers$HERE.normalDay) %>%
       addMarkers(~x, ~y,popup = ~as.character(facility), icon = greenLeafIcon)

##NUMBER OF PEDIATRIC HOSPITALS IN PHILADELPHIA AND PITTSBURGH Here we can zoom in and zoom out to see which of the 2 major cities have more number of pediatric hospitals.

icons <- awesomeIcons(icon = "home",
                      iconColor = "glass",
                      library = "glyphicon",
                      markerColor = "red")
noofpediatric <- subset(PA_dat,pediatric == "Y")
PA_pediatric <- subset(noofpediatric, county == c("Allegheny","Philadelphia"))
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length
## Warning in `==.default`(county, c("Allegheny", "Philadelphia")): longer
## object length is not a multiple of shorter object length
leaflet(PA_pediatric, width = "100%") %>%
        addTiles() %>%addProviderTiles(providers$
 Thunderforest.OpenCycleMap) %>%
        addProviderTiles(providers$
 Thunderforest.OpenCycleMap,options = providerTileOptions(opacity = 0.50)) %>%
       addProviderTiles(providers$HERE.normalDay) %>%
       addMarkers(~x, ~y,popup = ~as.character(facility), icon = icons) %>%
       addCircles(40.4406, -79.9959, weight = 1, radius = 30)

##NUMBER OF SURGICAL HOSPITALS IN PITTSBURGH

icons <- awesomeIcons(icon = "blue",
                      iconColor = "black",
                      library = "glyphicon",
                      markerColor = "red")
noofsurgical <- subset(PA_dat,surgical == "Y")
PA_surgical<- subset(noofsurgical, county == "Allegheny")
leaflet(PA_surgical, width = "100%") %>%
        addTiles() %>%addProviderTiles(providers$
 Thunderforest.OpenCycleMap) %>%
        addProviderTiles(providers$
 Thunderforest.OpenCycleMap,options = providerTileOptions(opacity = 0.50)) %>%
       addProviderTiles(providers$HERE.normalDay) %>%
       addMarkers(~x, ~y,popup = ~as.character(facility), icon = icons)