HW6

This analysis is for bars in Morgantown, West Virginia

addr <- read.csv("BarMorgantown.csv")

addr <- addr%>%
  dplyr::select(names(addr)[c(4,6,14,23:24,41,42)])

results <- st_as_sf(addr, coords=c("Longitude", "Latitude"), crs=4269,agr="constant")

results.proj<-st_transform(results,
                           crs = 32150)

mv1 <- mapview(results.proj)
mapshot(mv1, file = paste0(getwd(), "/map1.png"))

## PhantomJS not found. You can install it with webshot::install_phantomjs(). If it is installed, please make sure the phantomjs executable can be found via the PATH variable.

mv1

library(tmap)
library(tmaptools)


mean_feature<-apply(st_coordinates(results.proj), MARGIN = 2, FUN = mean)
mean_feature<-data.frame(place="meanfeature", x=mean_feature[1], y= mean_feature[2])
mean_feature<-st_as_sf(mean_feature, coords = c("x", "y"), crs= 32150)


tmap_mode("view")

## tmap mode set to interactive viewing

tm_basemap("OpenStreetMap.Mapnik")+
  tm_shape(results.proj, size = .2)+
  tm_dots(col = "green")+
  tm_shape(mean_feature)+
  tm_dots(col = "red", size = .2)

median_feature<-apply(st_coordinates(results.proj),
                      MARGIN = 2,
                      FUN = median)

median_feature<-data.frame(place="medianfeature",
                           x=median_feature[1],
                           y= median_feature[2])

median_feature<-st_as_sf(median_feature,
                         coords = c("x", "y"),
                         crs= 32150)


tmap_mode("view")

## tmap mode set to interactive viewing

tm_basemap("OpenStreetMap.Mapnik")+
  tm_shape(results.proj)+
  tm_dots(col = "green", size = .2)+
  tm_shape(mean_feature)+
  tm_dots(col = "red", size=.2)+
  tm_shape(median_feature)+
  tm_dots(col = "blue", size = .2)

chull <- st_convex_hull(st_union(results))

tmap_mode("view")

## tmap mode set to interactive viewing

tm_basemap("OpenStreetMap.Mapnik")+
  tm_shape(results.proj)+
  tm_dots(col = "green")+
  tm_shape(chull)+
  tm_polygons(col = "grey", alpha = .5)

library(SpatialKDE)

grid_groc <- results.proj %>%
  create_grid_rectangular(cell_size = 1000, side_offset = 2000)

kde <- results.proj%>%
  kde(band_width = 2000, kernel= "quartic", grid = grid_groc)

tm_shape(kde)+
  tm_polygons(col = "kde_value", palette= "viridis", title = "Density Estimate")+
  tm_shape(results.proj)+
  tm_dots()

## Getting data from the 2015-2019 5-year ACS

## Downloading feature geometry from the Census website.  To cache shapefiles for use in future sessions, set `options(tigris_use_cache = TRUE)`.

## Using the ACS Data Profile

spjoin<-st_join(results.proj, mgt_acs2)
head(spjoin)

## Simple feature collection with 6 features and 7 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: 556581.7 ymin: 125256.8 xmax: 560676.8 ymax: 128394
## Projected CRS: NAD83 / West Virginia North
##                          Business.Name               Physical.Address
## 1 Boston Beanery Restaurant and Tavern                383 Patteson Dr
## 2                  LongHorn Steakhouse 1427 University Town Centre Dr
## 3                Morgantown Brewing Co            1291 University Ave
## 4            Mountain State Brewing Co                     54 Clay St
## 5                      Texas Roadhouse          3505 Monongahela Blvd
## 6         Mario's Fishbowl Bar & Grill            3117 University Ave
##   Physical.ZIP Corporate.Employee.Size Revenue...Yr       GEOID totpop
## 1        26505                           $5,018,000 54061012000   5363
## 2        26501                  23,724   $4,665,000 54061011200   4695
## 3        26505                           $2,718,000 54061010102   5417
## 4        26501                           $2,378,000 54061011000   4276
## 5        26505                  35,554   $2,262,000 54061010400   4421
## 6        26505                           $2,262,000 54061010400   4421
##                    geometry
## 1 POINT (559867.8 127985.2)
## 2 POINT (556581.7 128225.4)
## 3 POINT (560676.8 125605.3)
## 4 POINT (560373.9 125256.8)
## 5   POINT (558163.4 128394)
## 6 POINT (559606.8 128234.5)

tmap_mode("view")

## tmap mode set to interactive viewing

tm_basemap("OpenStreetMap.Mapnik")+
  tm_shape(spjoin, is.master = T)+
  tm_dots("totpop", size = .1)+
  tm_shape(mgt_acs2)+
  tm_polygons(palette = "red", alpha = .1)

mgt_acs2$nbar<- lengths(st_intersects(mgt_acs2, results.proj))

mgt_acs2$bar_pc <- 1000*(mgt_acs2$nbar/mgt_acs2$totpop)

tmap_mode("plot")

## tmap mode set to plotting

m<-tm_basemap("OpenStreetMap.Mapnik")+
  tm_shape(mgt_acs2)+
  tm_polygons("bar_pc")+
  tm_shape(spjoin, is.master = T)+
  tm_dots( size = .01)

m

library(spatstat)
bar.pp<-as.ppp(as(results.proj, "Spatial"))

plot(nearest.neighbour(bar.pp))