Surratt Homework 4

Setup

Geocoe addresses using the Census Bureau’s batch geocoding service. The work will be done by the censusxy package.

Read in the address data

Use the form the CSV created in the Lab 3 exercise.

library(censusxy)
library(data.table)
library(mapview)
library(readr)
library(sf)

## Linking to GEOS 3.11.0, GDAL 3.5.3, PROJ 9.1.0; sf_use_s2() is TRUE

addr<-read.csv("/Users/briansurratt/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/DEM 7093 GIS/Homework 4/wic_west_side.csv") # Note this is the RAW data*

# head(addr)

addr1<-addr[c(6, 12:14)]

names(addr1)<-c("street", "city", "st", "zip")

# (addr1)

wic_addr<-cxy_geocode(addr1,
                     street = "street",
                     city = "city",
                     state ="st",
                     zip = "zip",
                     class="sf",
                     output = "simple")

## 5 rows removed to create an sf object. These were addresses that the geocoder could not match.

Basic interactive map of the points

mapview(wic_addr, layer.name="WIC Services mapped by addresses", col.regions = "blue")

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

wic_latlong <- st_transform(wic_latlong, crs = 2278)

mapview(wic_latlong, layer.name="WIC Services mapped by lat/long", col.regions = "red") 

Save the results if you want

We can write the results out to a shapefile now

st_write(wic_addr, dsn="/Users/briansurratt/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/DEM 7093 GIS/Homework 4/", layer="wic_addr", driver = "ESRI Shapefile", delete_layer = T, append=T)

## Deleting layer `wic_addr' using driver `ESRI Shapefile'
## Updating layer `wic_addr' to data source `/Users/briansurratt/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/DEM 7093 GIS/Homework 4/' using driver `ESRI Shapefile'
## Writing 98 features with 4 fields and geometry type Point.

# I got an error when trying to write this file.
# st_write(wic_latlong, dsn="/Users/briansurratt/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/DEM 7093 GIS/Homework 4/", layer="wic_latlong", driver = "ESRI Shapefile", delete_layer = T, append=T)

g_addr<-read.csv("/Users/briansurratt/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/DEM 7093 GIS/Homework 4/grocery_west_side.csv")

# head(g_addr)

g_addr1 <- g_addr[c(6, 12:14)]

names(g_addr1)<-c("street", "city", "st", "zip")

# head(g_addr1)

gro_addr<-cxy_geocode(g_addr1,
                     street = "street",
                     city = "city",
                     state ="st",
                     zip = "zip",
                     class="sf",
                     output = "simple")

## 4 rows removed to create an sf object. These were addresses that the geocoder could not match.

Basic interactive map of the points

mapview(gro_addr, layer.name="Grocery Stores mapped by Address", col.regions = "green")

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

gro_latlong <- st_transform(gro_latlong, crs = 2278)

mapview(gro_latlong, layer.name="Grocery Stores mapped by lat/long", col.regions = "orange") 

st_write(gro_addr, dsn="/Users/briansurratt/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/DEM 7093 GIS/Homework 4/", layer="gro_addr", driver = "ESRI Shapefile", delete_layer = T, append=T)

## Deleting layer `gro_addr' using driver `ESRI Shapefile'
## Updating layer `gro_addr' to data source `/Users/briansurratt/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/DEM 7093 GIS/Homework 4/' using driver `ESRI Shapefile'
## Writing 54 features with 4 fields and geometry type Point.

# I got an error when trying to write this file.
# st_write(gro_latlong, dsn="/Users/briansurratt/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/DEM 7093 GIS/Homework 4/", layer="gro_latlong", driver = "ESRI Shapefile", delete_layer = T, append=T)

Using R, create one map that has all 4 point elements.

mapview(wic_addr, layer.name="WIC Services mapped by addresses", col.regions = "blue") +
  mapview(wic_latlong, layer.name="WIC Services mapped by lat/long", col.regions = "red") +
  mapview(gro_addr, layer.name="Grocery Stores mapped by Address", col.regions = "green") +
  mapview(gro_latlong, layer.name="Grocery Stores mapped by lat/long", col.regions = "orange")