GIS HW 5
Jules Gonzalez
2023-04-03
library(mapview, quietly = TRUE)
library(sf,quietly = TRUE)## Linking to GEOS 3.11.0, GDAL 3.5.3, PROJ 9.1.0; sf_use_s2() is TRUElibrary(censusxy,quietly = TRUE)
library(dplyr,quietly = TRUE)##
## 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(tidycensus,quietly = TRUE)
library(maptools,quietly = TRUE)## Checking rgeos availability: TRUE
## Please note that 'maptools' will be retired during 2023,
## plan transition at your earliest convenience;
## some functionality will be moved to 'sp'.addr<-read.csv("/Users/juliegonzalez/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/Spring 2023 GIS/drinking.csv")
addr<-addr[c(6, 12:14)]
names(addr)<-c("street", "city", "st", "zip")
head(addr)## street city st zip
## 1 108 King William San Antonio TX 78204
## 2 905 Nogalitos St San Antonio TX 78204
## 3 507 Ruiz St San Antonio TX 78207
## 4 514 W Commerce St San Antonio TX 78207
## 5 5721 W Commerce St San Antonio TX 78237
## 6 254 Hobart St San Antonio TX 78237Using Lat/Long
addr<-read.csv("/Users/juliegonzalez/Library/CloudStorage/OneDrive-UniversityofTexasatSanAntonio/Spring 2023 GIS/drinking.csv")
results <- st_as_sf(addr, coords=c("Longitude", "Latitude"), crs=4269,agr="constant")
results.proj<-st_transform(results,
crs = 2278)#results<-cxy_geocode(addr,
#street = "street",
#city = "city",
#state ="st",
#zip = "zip",
#class="sf",
#output = "simple")results.proj<-st_transform(results,
crs = 2278)mapview(results.proj)Homework 5:
Perform the following analysis of these data: 1. Mean Center of the points - Make a map showing this. 2. Convex hull of the points - Make a map showing this. 3. 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. 4. Nearest Neighbor Analysis of the bar layer. Report the Z-score.
Mean Center of the points
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= 2278)
mapview(mean_feature, col.regions="red")+mapview( results)Convex hull of the points
chull <- st_convex_hull(st_union(results))
mapview(chull)+
mapview(results, col.regions = "green")Spatial Join
library(tidycensus, quietly = TRUE)
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")## Warning: • You have not set a Census API key. Users without a key are limited to 500
## queries per day and may experience performance limitations.
## ℹ For best results, get a Census API key at
## http://api.census.gov/data/key_signup.html and then supply the key to the
## `census_api_key()` function to use it throughout your tidycensus session.
## This warning is displayed once per session.##
|
| | 0%
|
| | 1%
|
|= | 1%
|
|= | 2%
|
|== | 2%
|
|== | 3%
|
|== | 4%
|
|=== | 4%
|
|=== | 5%
|
|==== | 5%
|
|==== | 6%
|
|===== | 7%
|
|===== | 8%
|
|====== | 8%
|
|====== | 9%
|
|======= | 10%
|
|======= | 11%
|
|======== | 11%
|
|======== | 12%
|
|========= | 13%
|
|========== | 14%
|
|========== | 15%
|
|=========== | 16%
|
|============ | 17%
|
|============ | 18%
|
|============= | 18%
|
|============= | 19%
|
|============== | 21%
|
|=============== | 22%
|
|================ | 23%
|
|================= | 25%
|
|================== | 26%
|
|=================== | 27%
|
|==================== | 28%
|
|===================== | 29%
|
|===================== | 31%
|
|======================= | 33%
|
|======================== | 35%
|
|========================= | 36%
|
|=========================== | 38%
|
|============================ | 40%
|
|============================= | 41%
|
|=============================== | 44%
|
|=============================== | 45%
|
|================================ | 46%
|
|=================================== | 50%
|
|==================================== | 51%
|
|==================================== | 52%
|
|===================================== | 53%
|
|====================================== | 54%
|
|======================================= | 56%
|
|=============================================== | 67%
|
|================================================= | 70%
|
|================================================== | 72%
|
|=================================================== | 73%
|
|==================================================== | 74%
|
|===================================================== | 76%
|
|====================================================== | 77%
|
|======================================================= | 78%
|
|======================================================= | 79%
|
|============================================================== | 88%
|
|============================================================== | 89%
|
|=============================================================== | 90%
|
|=============================================================== | 91%
|
|================================================================ | 92%
|
|================================================================= | 93%
|
|================================================================== | 95%
|
|=================================================================== | 96%
|
|===================================================================== | 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(results.proj, sa_acs2)
#head(spjoin)mapview(spjoin["ppov"])+mapview(sa_trol)#summary(results.proj)
#code gave me 23 results Nearest Neighbor Analysis
library(maptools)
library(spatstat)## Loading required package: spatstat.data## Loading required package: spatstat.geom## spatstat.geom 3.1-0## Loading required package: spatstat.random## spatstat.random 3.1-4## Loading required package: spatstat.explore## 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## Loading required package: rpart## spatstat.model 3.2-1## Loading required package: spatstat.linnet## spatstat.linnet 3.0-6##
## spatstat 3.0-3
## For an introduction to spatstat, type 'beginner'wic.pp<-as.ppp(as(results.proj, "Spatial"))
plot(nearest.neighbour(wic.pp))
#install.packages("remotes")
remotes::install_github("paleolimbot/qgisprocess") #Select 1 for ALL; next there will be a pop-up box, select Yes.## Skipping install of 'qgisprocess' from a github remote, the SHA1 (464d75b3) has not changed since last install.
## Use `force = TRUE` to force installationforce=TRUElibrary(qgisprocess) #load the package## Attempting to load the cache ... Success!
## QGIS version: 3.28.2-Firenze
## Having access to 1019 algorithms from 6 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':
## /Applications/QGIS.app/Contents/MacOS/bin/qgis_process
## Trying command '/Applications/QGIS.app/Contents/MacOS/bin/qgis_process'
## Success!
## Now using 'qgis_process' at '/Applications/QGIS.app/Contents/MacOS/bin/qgis_process'.
## >>> 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.
## 3 out of 4 available processing provider plugins are enabled.
## You now have access to 1019 algorithms from 6 QGIS processing providers.
##
## >>> Run `qgis_enable_plugins()` to enable 1 disabled plugin(s) and
## access their algorithms: otbprovider
##
## Saving configuration to '~/Library/Caches/R-qgisprocess/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())## # A tibble: 6 × 24
## provider provider_title algorithm algorithm_id algorithm_title
## <chr> <chr> <chr> <chr> <chr>
## 1 3d QGIS (3D) 3d:tessellate tessellate Tessellate
## 2 gdal GDAL gdal:aspect aspect Aspect
## 3 gdal GDAL gdal:assignprojection assignproject… Assign project…
## 4 gdal GDAL gdal:buffervectors buffervectors Buffer vectors
## 5 gdal GDAL gdal:buildvirtualraster buildvirtualr… Build virtual …
## 6 gdal GDAL gdal:buildvirtualvector buildvirtualv… Build virtual …
## # ℹ 19 more variables: provider_can_be_activated <lgl>,
## # provider_is_active <lgl>, provider_long_name <chr>, provider_version <chr>,
## # provider_warning <chr>, can_cancel <lgl>, deprecated <lgl>, group <chr>,
## # has_known_issues <lgl>, help_url <chr>, requires_matching_crs <lgl>,
## # short_description <chr>, tags <list>, default_raster_file_extension <chr>,
## # default_vector_file_extension <chr>,
## # supported_output_raster_extensions <list>, …qgis_algorithms()[grepl(pattern = "nearest", x = qgis_algorithms()$algorithm),]## # A tibble: 10 × 24
## provider provider_title algorithm algorithm_id algorithm_title
## <chr> <chr> <chr> <chr> <chr>
## 1 gdal GDAL gdal:gridinversedist… gridinverse… Grid (IDW with…
## 2 gdal GDAL gdal:gridnearestneig… gridnearest… Grid (Nearest …
## 3 native QGIS (native c++) native:angletonearest angletonear… Align points t…
## 4 native QGIS (native c++) native:joinbynearest joinbyneare… Join attribute…
## 5 native QGIS (native c++) native:nearestneighb… nearestneig… Nearest neighb…
## 6 qgis QGIS qgis:distancetoneare… distanceton… Distance to ne…
## 7 qgis QGIS qgis:distancetoneare… distanceton… Distance to ne…
## 8 qgis QGIS qgis:knearestconcave… knearestcon… Concave hull (…
## 9 saga SAGA saga:knearestneighbo… knearestnei… K-nearest neig…
## 10 saga SAGA saga:nearestneighbour nearestneig… Nearest neighb…
## # ℹ 19 more variables: provider_can_be_activated <lgl>,
## # provider_is_active <lgl>, provider_long_name <chr>, provider_version <chr>,
## # provider_warning <chr>, can_cancel <lgl>, deprecated <lgl>, group <chr>,
## # has_known_issues <lgl>, help_url <chr>, requires_matching_crs <lgl>,
## # short_description <chr>, tags <list>, default_raster_file_extension <chr>,
## # default_vector_file_extension <chr>,
## # supported_output_raster_extensions <list>, …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-score# Run the native:nearestneighbouranalysis algorithm
wic_nn <- qgis_run_algorithm(alg = "native:nearestneighbouranalysis",
INPUT = results.proj,
OUTPUT_HTML_FILE = file.path(tempdir(), "wicnn.html"))
# Read the output HTML file
output_html_file <- wic_nn$OUTPUT_HTML_FILE
output_html_content <- readLines(output_html_file)## Warning in readLines(output_html_file): incomplete final line found on
## '/var/folders/2n/7p4156kj3h50qq59bkrzr6000000gn/T//Rtmp9XrbhD/wicnn.html'# Print the output
cat(output_html_content, sep = "\n")## <html><head><meta http-equiv="Content-Type" content="text/html;charset=utf-8"/></head><body>
## <p>Observed mean distance: 1981.16496462799</p>
## <p>Expected mean distance: 1965.10622343587</p>
## <p>Nearest neighbour index: 1.00817194562</p>
## <p>Number of points: 23</p>
## <p>Z-Score: 0.07497565504</p>
## </body></html>Z-score:
#Z-Score: Z-Score: 0.07497565504
#Z - scores shows how far a data point is from the average value, relative to the standard deviation.
#The Z-score of 0.07497565504 is really close to zero, which indicates that the point pattern of drinking businesses is not significantly different from a random distribution. The drinking locations in the selected zip codes are spread out. The drinking locations are not tightly clustered nor evenly spaced apart, they are more randomly located.