Homework5
Lydia Okabe
2023-04-03
##Load libraries
library(mapview)
library(sf)## Linking to GEOS 3.9.1, GDAL 3.4.3, PROJ 7.2.1; sf_use_s2() is TRUElibrary(censusxy)
library(dplyr)## Warning: package 'dplyr' was built under R version 4.2.3##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(readr)## Warning: package 'readr' was built under R version 4.2.3library(data.table)##
## Attaching package: 'data.table'## The following objects are masked from 'package:dplyr':
##
## between, first, last##LOADING WIC STORE DATA
Drinking <- read_csv("C:/Users/okabe/OneDrive/Pictures/Stats 2/drinking_places.csv")## New names:
## Rows: 23 Columns: 228
## ── Column specification
## ──────────────────────────────────────────────────────── Delimiter: "," chr
## (76): Source...1, Date, Obsolescence Date, Business Name, Legal Name, P... dbl
## (40): Physical Address Number, Physical ZIP, Physical ZIP 4, Location E... lgl
## (112): Physical Post Direction, Corporate Employee Size, Mailing Post Di...
## ℹ Use `spec()` to retrieve the full column specification for this data. ℹ
## Specify the column types or set `show_col_types = FALSE` to quiet this message.
## • `Source` -> `Source...1`
## • `Source` -> `Source...227`
## • `` -> `...228`head(Drinking)##OR
#Drinking1<-read.csv("C:/Users/okabe/OneDrive/Pictures/Stats 2/drinking_places.csv")
#drinkingresults <- st_as_sf(Drinking1, coords=c("Longitude", "Latitude"), crs=2278,agr="constant")
#drinkingresults.proj<-st_transform(drinkingresults,
#crs = 2278)Drinkingplaces<- Drinking[c(6, 12:14)]
names(Drinkingplaces)<-c("street", "city", "st", "zip" )
head(Drinkingplaces)library(censusxy)
drinkingresults<-cxy_geocode(Drinkingplaces,
street = "street",
city = "city",
state ="st",
zip = "zip",
class="sf",
output = "simple")## 2 rows removed to create an sf object. These were addresses that the geocoder could not match.drinkingresults.proj<-st_transform(drinkingresults,
crs = 2278)mapview(drinkingresults.proj)Perform the following analysis of these data: 1. Mean Center of the points - Make a map showing this.
mean_feature<-apply(st_coordinates(drinkingresults.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= 2278)
mapview(mean_feature, col.regions="red")+mapview( drinkingresults)- Convex hull of the points - Make a map showing this.
chull <- st_convex_hull(st_union(drinkingresults))
mapview(chull)+
mapview(drinkingresults, col.regions = "green")- Spatial join of the points to the ACS layer (in sa_tracts_2019.gpkg). Set the Symbology to Graduated, Pretty Breaks, 7 Classes, for the variable ppov. - Make a map showing this.
library(tidycensus)
library(dplyr)
#load census tract data
sa_acs<-get_acs(geography = "tract",
state="TX",
county = "Bexar",
year = 2019,
variables=c( "DP05_0001E", "DP03_0009P", "DP03_0062E", "DP03_0119PE",
"DP05_0001E","DP02_0009PE","DP02_0008PE","DP02_0040E","DP02_0038E",
"DP02_0066PE","DP02_0067PE","DP02_0080PE","DP02_0092PE",
"DP03_0005PE","DP03_0028PE","DP03_0062E","DP03_0099PE","DP03_0101PE",
"DP03_0119PE","DP04_0046PE","DP05_0072PE","DP05_0073PE",
"DP05_0066PE", "DP05_0072PE", "DP02_0113PE") ,
geometry = T, output = "wide")## 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
## Using the ACS Data Profile##
|
| | 0%
|
|= | 1%
|
|= | 2%
|
|== | 2%
|
|== | 3%
|
|=== | 4%
|
|=== | 5%
|
|==== | 5%
|
|==== | 6%
|
|===== | 7%
|
|===== | 8%
|
|====== | 8%
|
|====== | 9%
|
|======= | 10%
|
|======= | 11%
|
|======== | 11%
|
|======== | 12%
|
|========= | 12%
|
|========= | 13%
|
|========= | 14%
|
|========== | 14%
|
|========== | 15%
|
|=========== | 15%
|
|=========== | 16%
|
|============ | 16%
|
|============ | 17%
|
|============ | 18%
|
|============= | 18%
|
|============= | 19%
|
|============== | 19%
|
|============== | 20%
|
|============== | 21%
|
|=============== | 21%
|
|=============== | 22%
|
|================ | 22%
|
|================ | 23%
|
|================= | 24%
|
|================= | 25%
|
|================== | 25%
|
|================== | 26%
|
|=================== | 27%
|
|==================== | 28%
|
|==================== | 29%
|
|===================== | 30%
|
|====================== | 31%
|
|====================== | 32%
|
|======================= | 32%
|
|======================= | 33%
|
|======================== | 34%
|
|======================== | 35%
|
|========================= | 35%
|
|========================= | 36%
|
|========================== | 37%
|
|========================== | 38%
|
|=========================== | 38%
|
|=========================== | 39%
|
|============================ | 39%
|
|============================ | 40%
|
|============================ | 41%
|
|============================= | 41%
|
|============================= | 42%
|
|============================== | 43%
|
|=============================== | 44%
|
|=============================== | 45%
|
|================================ | 45%
|
|================================ | 46%
|
|================================= | 47%
|
|================================= | 48%
|
|================================== | 48%
|
|================================== | 49%
|
|=================================== | 49%
|
|=================================== | 50%
|
|=================================== | 51%
|
|==================================== | 51%
|
|==================================== | 52%
|
|===================================== | 52%
|
|===================================== | 53%
|
|====================================== | 54%
|
|====================================== | 55%
|
|======================================= | 55%
|
|======================================= | 56%
|
|======================================== | 57%
|
|======================================== | 58%
|
|========================================= | 58%
|
|========================================= | 59%
|
|========================================== | 59%
|
|========================================== | 60%
|
|========================================== | 61%
|
|=========================================== | 61%
|
|=========================================== | 62%
|
|============================================ | 62%
|
|============================================ | 63%
|
|============================================= | 64%
|
|============================================= | 65%
|
|============================================== | 65%
|
|============================================== | 66%
|
|=============================================== | 66%
|
|=============================================== | 67%
|
|=============================================== | 68%
|
|================================================ | 68%
|
|================================================ | 69%
|
|================================================= | 70%
|
|================================================= | 71%
|
|================================================== | 71%
|
|================================================== | 72%
|
|=================================================== | 72%
|
|=================================================== | 73%
|
|==================================================== | 74%
|
|==================================================== | 75%
|
|===================================================== | 75%
|
|===================================================== | 76%
|
|====================================================== | 77%
|
|====================================================== | 78%
|
|======================================================= | 78%
|
|======================================================= | 79%
|
|======================================================== | 79%
|
|======================================================== | 80%
|
|========================================================= | 81%
|
|========================================================= | 82%
|
|========================================================== | 82%
|
|========================================================== | 83%
|
|=========================================================== | 84%
|
|=========================================================== | 85%
|
|============================================================ | 85%
|
|============================================================ | 86%
|
|============================================================= | 87%
|
|============================================================= | 88%
|
|============================================================== | 88%
|
|============================================================== | 89%
|
|=============================================================== | 90%
|
|=============================================================== | 91%
|
|================================================================ | 91%
|
|================================================================ | 92%
|
|================================================================= | 92%
|
|================================================================= | 93%
|
|================================================================= | 94%
|
|================================================================== | 94%
|
|================================================================== | 95%
|
|=================================================================== | 95%
|
|=================================================================== | 96%
|
|==================================================================== | 97%
|
|==================================================================== | 98%
|
|===================================================================== | 98%
|
|===================================================================== | 99%
|
|======================================================================| 100%#rename variables and filter missing cases
sa_acs2<-sa_acs%>%
mutate(totpop= DP05_0001E, pwhite=DP05_0072PE,
pblack=DP05_0073PE , phisp=DP05_0066PE,
phsormore=DP02_0066PE,punemp=DP03_0009PE, medhhinc=DP03_0062E,
ppov=DP03_0119PE)%>%
dplyr::select(GEOID, totpop, pblack, pwhite, phisp, punemp, medhhinc, ppov)
sa_acs2<-st_transform(sa_acs2, crs = 2278)
sa_trol<-st_cast(sa_acs2, "MULTILINESTRING")spjoin<-st_join(drinkingresults.proj, sa_acs2)
head(spjoin)summary(drinkingresults.proj)## street city st zip
## Length:21 Length:21 Length:21 Min. :78204
## Class :character Class :character Class :character 1st Qu.:78204
## Mode :character Mode :character Mode :character Median :78207
## Mean :78212
## 3rd Qu.:78207
## Max. :78237
## geometry
## POINT :21
## epsg:2278 : 0
## +proj=lcc ...: 0
##
##
## mapview(spjoin["ppov"])+mapview(sa_trol)- Nearest Neighbor Analysis of the bar layer. Report the Z-score
#NOTE: this chunk was set to eval=F; it will not run for me. See link above for what Dr. Sparks got.
library(spatstat)## Warning: package 'spatstat' was built under R version 4.2.3## Loading required package: spatstat.data## Warning: package 'spatstat.data' was built under R version 4.2.3## Loading required package: spatstat.geom## Warning: package 'spatstat.geom' was built under R version 4.2.3## spatstat.geom 3.1-0##
## Attaching package: 'spatstat.geom'## The following object is masked from 'package:data.table':
##
## shift## Loading required package: spatstat.random## Warning: package 'spatstat.random' was built under R version 4.2.3## spatstat.random 3.1-4## Loading required package: spatstat.explore## Warning: package 'spatstat.explore' was built under R version 4.2.3## Loading required package: nlme##
## Attaching package: 'nlme'## The following object is masked from 'package:dplyr':
##
## collapse## spatstat.explore 3.1-0## Loading required package: spatstat.model## Warning: package 'spatstat.model' was built under R version 4.2.3## Loading required package: rpart## spatstat.model 3.2-1## Loading required package: spatstat.linnet## Warning: package 'spatstat.linnet' was built under R version 4.2.3## spatstat.linnet 3.0-6##
## spatstat 3.0-3
## For an introduction to spatstat, type 'beginner'library(maptools)## Loading required package: sp## Checking rgeos availability: TRUE
## Please note that 'maptools' will be retired by the end of 2023,
## plan transition at your earliest convenience;
## some functionality will be moved to 'sp'.drinking.pp<-as.ppp(as(drinkingresults.proj, "Spatial"))
plot(nearest.neighbour(drinking.pp))
install.packages("remotes")
remotes::install_github("paleolimbot/qgisprocess") #Select 1 for ALL; next there will be a pop-up box, select Yes.library(qgisprocess) #load the package## Attempting to load the cache ... Success!
## QGIS version: 3.28.2-Firenze
## Having access to 1169 algorithms from 7 QGIS processing providers.
## Run `qgis_configure(use_cached_data = TRUE)` to reload cache and get more details.
## >>> Run `qgis_enable_plugins()` to enable 1 disabled plugin(s) and
## access their algorithms: otbproviderqgis_configure() #set up QGIS - find the executable## getOption('qgisprocess.path') was not found.
## Sys.getenv('R_QGISPROCESS_PATH') was not found.
## Trying 'qgis_process' on PATH...
## 'qgis_process' is not available on PATH.
## Found 1 QGIS installation containing 'qgis_process':
## C:/Program Files/QGIS 3.28.2/bin/qgis_process-qgis.bat
## Trying command 'C:/Program Files/QGIS 3.28.2/bin/qgis_process-qgis.bat'
## Success!
## Now using 'qgis_process' at 'C:/Program Files/QGIS 3.28.2/bin/qgis_process-qgis.bat'.
## >>> If you need another installed QGIS instance, run `qgis_configure()`;
## see `?qgis_configure` if you need to preset the path of 'qgis_process'.
##
## QGIS version is now set to: 3.28.2-Firenze
## Using JSON for output serialization.
## Using JSON for input serialization.
## 4 out of 5 available processing provider plugins are enabled.
## You now have access to 1169 algorithms from 7 QGIS processing providers.
##
## >>> Run `qgis_enable_plugins()` to enable 1 disabled plugin(s) and
## access their algorithms: otbprovider
##
## Saving configuration to 'C:\Users\okabe\AppData\Local/R-qgisprocess/R-qgisprocess/Cache/cache-0.0.0.9000.rds'
## Use qgis_algorithms(), qgis_providers(), qgis_plugins(), qgis_use_json_output(),
## qgis_path() and qgis_version() to inspect the cache environment.# qgis_algorithms() lists all the available routines in QGIS
head(qgis_algorithms())algs<-qgis_algorithms()
algs[grepl(pattern = "Spatial", x = algs$algorithm ),]## Warning in `[<-.data.frame`(`*tmp*`, is_list, value = list(`18` = "<>", :
## replacement element 1 has 1 row to replace 0 rows## Warning in `[<-.data.frame`(`*tmp*`, is_list, value = list(`18` = "<>", :
## replacement element 2 has 1 row to replace 0 rows## Warning in `[<-.data.frame`(`*tmp*`, is_list, value = list(`18` = "<>", :
## replacement element 3 has 1 row to replace 0 rows## Warning in `[<-.data.frame`(`*tmp*`, is_list, value = list(`18` = "<>", :
## replacement element 4 has 1 row to replace 0 rowsqgis_show_help("qgis:heatmapkerneldensityestimation")## Heatmap (Kernel Density Estimation) (qgis:heatmapkerneldensityestimation)
##
## ----------------
## Description
## ----------------
## Creates a density (heatmap) raster of an input point vector layer using kernel density estimation. Heatmaps allow easy identification of hotspots and clustering of points.
## The density is calculated based on the number of points in a location, with larger numbers of clustered points resulting in larger values.
##
##
## ----------------
## Arguments
## ----------------
##
## INPUT: Point layer
## Argument type: source
## Acceptable values:
## - Path to a vector layer
## RADIUS: Radius
## Default value: 100
## Argument type: distance
## Acceptable values:
## - A numeric value
## RADIUS_FIELD: Radius from field (optional)
## Argument type: field
## Acceptable values:
## - The name of an existing field
## - ; delimited list of existing field names
## PIXEL_SIZE: Output raster size
## Default value: 0.1
## Argument type: number
## Acceptable values:
## - A numeric value
## WEIGHT_FIELD: Weight from field (optional)
## Argument type: field
## Acceptable values:
## - The name of an existing field
## - ; delimited list of existing field names
## KERNEL: Kernel shape
## Default value: 0
## Argument type: enum
## Available values:
## - 0: Quartic
## - 1: Triangular
## - 2: Uniform
## - 3: Triweight
## - 4: Epanechnikov
## Acceptable values:
## - Number of selected option, e.g. '1'
## - Comma separated list of options, e.g. '1,3'
## DECAY: Decay ratio (Triangular kernels only) (optional)
## Default value: 0
## Argument type: number
## Acceptable values:
## - A numeric value
## OUTPUT_VALUE: Output value scaling
## Default value: 0
## Argument type: enum
## Available values:
## - 0: Raw
## - 1: Scaled
## Acceptable values:
## - Number of selected option, e.g. '1'
## - Comma separated list of options, e.g. '1,3'
## OUTPUT: Heatmap
## Argument type: rasterDestination
## Acceptable values:
## - Path for new raster layer
##
## ----------------
## Outputs
## ----------------
##
## OUTPUT: <outputRaster>
## Heatmapalgs[grepl(pattern = "nearest", x = algs$algorithm ),]qgis_show_help("native:nearestneighbouranalysis")## Nearest neighbour analysis (native:nearestneighbouranalysis)
##
## ----------------
## Description
## ----------------
## This algorithm performs nearest neighbor analysis for a point layer.
##
## The output describes how the data are distributed (clustered, randomly or distributed).
##
## Output is generated as an HTML file with the computed statistical values.
##
## ----------------
## Arguments
## ----------------
##
## INPUT: Input layer
## Argument type: source
## Acceptable values:
## - Path to a vector layer
## OUTPUT_HTML_FILE: Nearest neighbour (optional)
## Argument type: fileDestination
## Acceptable values:
## - Path for new file
##
## ----------------
## Outputs
## ----------------
##
## OUTPUT_HTML_FILE: <outputHtml>
## Nearest neighbour
## OBSERVED_MD: <outputNumber>
## Observed mean distance
## EXPECTED_MD: <outputNumber>
## Expected mean distance
## NN_INDEX: <outputNumber>
## Nearest neighbour index
## POINT_COUNT: <outputNumber>
## Number of points
## Z_SCORE: <outputNumber>
## Z-scoredrinking_nn<-qgis_run_algorithm(alg="native:nearestneighbouranalysis",
INPUT=drinkingresults.proj,
OUTPUT_HTML_FILE=file.path(tempdir(), "drinkingnn.html"),
load_output = TRUE)## Ignoring unknown input 'load_output'drinking_nn## <Result of `qgis_run_algorithm("native:nearestneighbouranalysis", ...)`>
## List of 6
## $ EXPECTED_MD : num 2045
## $ NN_INDEX : num 1.03
## $ OBSERVED_MD : num 2113
## $ OUTPUT_HTML_FILE: 'qgis_outputHtml' chr "C:\\Users\\okabe\\AppData\\Local\\Temp\\RtmpgbMtIx/drinkingnn.html"
## $ POINT_COUNT : num 21
## $ Z_SCORE : num 0.29##The z score is 0.29
##This means that in the zip codes 78207,78204 & 78237 the 0.29 z-score shows that there is a slight difference in the point patterns when compared to a random distribution. Also based on the maps above, we see that the businesses are spread out rather than clustered