# QUESTION 1
library(sp)
library(sf)
## Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 6.3.1
roads <- st_read("final/ExamData/SeminarR.gdb", layer="USA_roads")
## Reading layer `USA_roads' from data source `C:\Users\mwtro\OneDrive\Documents\r spatial\data\final\ExamData\SeminarR.gdb' using driver `OpenFileGDB'
## Simple feature collection with 679 features and 8 fields
## geometry type: MULTILINESTRING
## dimension: XY
## bbox: xmin: -6146140 ymin: 316215.5 xmax: 2242259 ymax: 5729382
## projected CRS: NAD83 / Conus Albers
counties <- st_read("final/ExamData/SeminarR.gdb", layer="USA_counties")
## Reading layer `USA_counties' from data source `C:\Users\mwtro\OneDrive\Documents\r spatial\data\final\ExamData\SeminarR.gdb' using driver `OpenFileGDB'
## Simple feature collection with 3140 features and 53 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -6293474 ymin: 311822.4 xmax: 2256319 ymax: 6198811
## projected CRS: NAD83 / Conus Albers
cities <- read.csv("final/ExamData/cities1.csv")
class(roads)
## [1] "sf" "data.frame"
class(counties)
## [1] "sf" "data.frame"
class(cities)
## [1] "data.frame"
colnames(cities)
## [1] "FID" "AREANAME" "CLASS" "ST" "STFIPS"
## [6] "HOUSEUNITS" "POPULATION" "lat" "long"
cities_sf <- st_as_sf(cities, coords = c("long", "lat")) # error in data? had to do x = long and y = lat
st_crs(cities_sf)
## Coordinate Reference System: NA
st_crs(counties)
## Coordinate Reference System:
## User input: NAD83 / Conus Albers
## wkt:
## PROJCRS["NAD83 / Conus Albers",
## BASEGEOGCRS["NAD83",
## DATUM["North American Datum 1983",
## ELLIPSOID["GRS 1980",6378137,298.257222101,
## LENGTHUNIT["metre",1]]],
## PRIMEM["Greenwich",0,
## ANGLEUNIT["degree",0.0174532925199433]],
## ID["EPSG",4269]],
## CONVERSION["Conus Albers",
## METHOD["Albers Equal Area",
## ID["EPSG",9822]],
## PARAMETER["Latitude of false origin",23,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8821]],
## PARAMETER["Longitude of false origin",-96,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8822]],
## PARAMETER["Latitude of 1st standard parallel",29.5,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8823]],
## PARAMETER["Latitude of 2nd standard parallel",45.5,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8824]],
## PARAMETER["Easting at false origin",0,
## LENGTHUNIT["metre",1],
## ID["EPSG",8826]],
## PARAMETER["Northing at false origin",0,
## LENGTHUNIT["metre",1],
## ID["EPSG",8827]]],
## CS[Cartesian,2],
## AXIS["easting (X)",east,
## ORDER[1],
## LENGTHUNIT["metre",1]],
## AXIS["northing (Y)",north,
## ORDER[2],
## LENGTHUNIT["metre",1]],
## USAGE[
## SCOPE["unknown"],
## AREA["USA - CONUS - onshore"],
## BBOX[24.41,-124.79,49.38,-66.91]],
## ID["EPSG",5070]]
st_crs(cities_sf)=5070
plot(cities_sf)
colnames(cities_sf)
## [1] "FID" "AREANAME" "CLASS" "ST" "STFIPS"
## [6] "HOUSEUNITS" "POPULATION" "geometry"
colnames(counties)
## [1] "NAME" "STATE_NAME" "STATE_FIPS" "CNTY_FIPS" "FIPS"
## [6] "POP1990" "POP1999" "POP90_SQMI" "HOUSEHOLDS" "MALES"
## [11] "FEMALES" "WHITE" "BLACK" "AMERI_ES" "ASIAN_PI"
## [16] "OTHER" "HISPANIC" "AGE_UNDER5" "AGE_5_17" "AGE_18_29"
## [21] "AGE_30_49" "AGE_50_64" "AGE_65_UP" "NEVERMARRY" "MARRIED"
## [26] "SEPARATED" "WIDOWED" "DIVORCED" "HSEHLD_1_M" "HSEHLD_1_F"
## [31] "MARHH_CHD" "MARHH_NO_C" "MHH_CHILD" "FHH_CHILD" "HSE_UNITS"
## [36] "VACANT" "OWNER_OCC" "RENTER_OCC" "MEDIAN_VAL" "MEDIANRENT"
## [41] "UNITS_1DET" "UNITS_1ATT" "UNITS2" "UNITS3_9" "UNITS10_49"
## [46] "UNITS50_UP" "MOBILEHOME" "NO_FARMS87" "AVG_SIZE87" "CROP_ACR87"
## [51] "AVG_SALE87" "Shape_Length" "Shape_Area" "Shape"
nrow(roads)
## [1] 679
roadbuffer <- st_buffer(roads, 2000) # 20,000 meters = 20 km, since projection is in m
head(roadbuffer)
## Simple feature collection with 6 features and 8 fields
## geometry type: POLYGON
## dimension: XY
## bbox: xmin: -2057948 ymin: 2430740 xmax: 230943.8 ymax: 2907057
## projected CRS: NAD83 / Conus Albers
## LENGTH TYPE ADMN_CLASS TOLL_RD RTE_NUM1 RTE_NUM2
## 1 140.931 Multi-Lane Divided Interstate N 82
## 2 186.060 Multi-Lane Divided Interstate N 84
## 3 8.642 Multi-Lane Divided Interstate N 94
## 4 9.439 Multi-Lane Divided Interstate N 94
## 5 8.208 Multi-Lane Divided Interstate N 35
## 6 22.106 Multi-Lane Divided Interstate N 494
## ROUTE Shape_Length Shape
## 1 Interstate 82 226804.48 POLYGON ((-1793423 2772752,...
## 2 Interstate (OR/ID/UT) 84 299427.61 POLYGON ((-1790088 2748427,...
## 3 Interstate 94 13908.32 POLYGON ((214800.5 2450524,...
## 4 Interstate 94 15190.89 POLYGON ((220305.8 2447169,...
## 5 Interstate 35 13208.87 POLYGON ((210561.9 2435146,...
## 6 Interstate 494 35576.04 POLYGON ((202050.2 2453856,...
st_crs(roadbuffer)
## Coordinate Reference System:
## User input: NAD83 / Conus Albers
## wkt:
## PROJCRS["NAD83 / Conus Albers",
## BASEGEOGCRS["NAD83",
## DATUM["North American Datum 1983",
## ELLIPSOID["GRS 1980",6378137,298.257222101,
## LENGTHUNIT["metre",1]]],
## PRIMEM["Greenwich",0,
## ANGLEUNIT["degree",0.0174532925199433]],
## ID["EPSG",4269]],
## CONVERSION["Conus Albers",
## METHOD["Albers Equal Area",
## ID["EPSG",9822]],
## PARAMETER["Latitude of false origin",23,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8821]],
## PARAMETER["Longitude of false origin",-96,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8822]],
## PARAMETER["Latitude of 1st standard parallel",29.5,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8823]],
## PARAMETER["Latitude of 2nd standard parallel",45.5,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8824]],
## PARAMETER["Easting at false origin",0,
## LENGTHUNIT["metre",1],
## ID["EPSG",8826]],
## PARAMETER["Northing at false origin",0,
## LENGTHUNIT["metre",1],
## ID["EPSG",8827]]],
## CS[Cartesian,2],
## AXIS["easting (X)",east,
## ORDER[1],
## LENGTHUNIT["metre",1]],
## AXIS["northing (Y)",north,
## ORDER[2],
## LENGTHUNIT["metre",1]],
## USAGE[
## SCOPE["unknown"],
## AREA["USA - CONUS - onshore"],
## BBOX[24.41,-124.79,49.38,-66.91]],
## ID["EPSG",5070]]
st_crs(cities_sf)
## Coordinate Reference System:
## User input: EPSG:5070
## wkt:
## PROJCRS["NAD83 / Conus Albers",
## BASEGEOGCRS["NAD83",
## DATUM["North American Datum 1983",
## ELLIPSOID["GRS 1980",6378137,298.257222101,
## LENGTHUNIT["metre",1]]],
## PRIMEM["Greenwich",0,
## ANGLEUNIT["degree",0.0174532925199433]],
## ID["EPSG",4269]],
## CONVERSION["Conus Albers",
## METHOD["Albers Equal Area",
## ID["EPSG",9822]],
## PARAMETER["Latitude of false origin",23,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8821]],
## PARAMETER["Longitude of false origin",-96,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8822]],
## PARAMETER["Latitude of 1st standard parallel",29.5,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8823]],
## PARAMETER["Latitude of 2nd standard parallel",45.5,
## ANGLEUNIT["degree",0.0174532925199433],
## ID["EPSG",8824]],
## PARAMETER["Easting at false origin",0,
## LENGTHUNIT["metre",1],
## ID["EPSG",8826]],
## PARAMETER["Northing at false origin",0,
## LENGTHUNIT["metre",1],
## ID["EPSG",8827]]],
## CS[Cartesian,2],
## AXIS["easting (X)",east,
## ORDER[1],
## LENGTHUNIT["metre",1]],
## AXIS["northing (Y)",north,
## ORDER[2],
## LENGTHUNIT["metre",1]],
## USAGE[
## SCOPE["unknown"],
## AREA["USA - CONUS - onshore"],
## BBOX[24.41,-124.79,49.38,-66.91]],
## ID["EPSG",5070]]
plot(cities_sf)
contains <- st_contains(cities_sf, roadbuffer)
nrow(contains)
## [1] 23434
intersects <- st_intersects(cities_sf, roadbuffer)
head(counties)
## Simple feature collection with 6 features and 53 fields
## geometry type: MULTIPOLYGON
## dimension: XY
## bbox: xmin: -1826165 ymin: 2819070 xmax: 117332.6 ymax: 3127678
## projected CRS: NAD83 / Conus Albers
## NAME STATE_NAME STATE_FIPS CNTY_FIPS FIPS POP1990 POP1999
## 1 Lake of the Woods Minnesota 27 077 27077 4076 4597
## 2 Ferry Washington 53 019 53019 6295 7150
## 3 Stevens Washington 53 065 53065 30948 39965
## 4 Okanogan Washington 53 047 53047 33350 38596
## 5 Pend Oreille Washington 53 051 53051 8915 11788
## 6 Boundary Idaho 16 021 16021 8332 9840
## POP90_SQMI HOUSEHOLDS MALES FEMALES WHITE BLACK AMERI_ES ASIAN_PI OTHER
## 1 2 1576 2037 2039 4042 1 19 10 4
## 2 3 2247 3280 3015 5084 20 1131 24 36
## 3 12 11241 15454 15494 28747 65 1807 179 150
## 4 6 12654 16828 16522 27615 52 3597 166 1920
## 5 6 3395 4426 4489 8640 12 206 25 32
## 6 7 2857 4252 4080 7950 3 150 26 203
## HISPANIC AGE_UNDER5 AGE_5_17 AGE_18_29 AGE_30_49 AGE_50_64 AGE_65_UP
## 1 25 337 791 549 1122 584 693
## 2 85 486 1499 887 1929 827 667
## 3 483 2271 7486 3586 9605 4145 3855
## 4 2779 2536 7051 4492 9749 4890 4632
## 5 120 660 1963 996 2670 1384 1242
## 6 310 635 2065 1138 2351 1121 1022
## NEVERMARRY MARRIED SEPARATED WIDOWED DIVORCED HSEHLD_1_M HSEHLD_1_F MARHH_CHD
## 1 538 2102 24 248 200 185 184 465
## 2 1118 2741 82 251 498 302 195 652
## 3 4037 14795 384 1486 2096 1146 1210 3618
## 4 5042 15320 561 1889 2406 1410 1633 3355
## 5 1124 4349 114 507 641 373 402 961
## 6 1247 3935 89 391 480 300 299 938
## MARHH_NO_C MHH_CHILD FHH_CHILD HSE_UNITS VACANT OWNER_OCC RENTER_OCC
## 1 565 26 57 3050 1474 1332 244
## 2 660 90 169 3239 992 1568 679
## 3 3591 292 731 14601 3360 8566 2675
## 4 3946 401 933 16629 3975 8439 4215
## 5 1157 75 229 5404 2009 2500 895
## 6 928 55 164 3242 385 2237 620
## MEDIAN_VAL MEDIANRENT UNITS_1DET UNITS_1ATT UNITS2 UNITS3_9 UNITS10_49
## 1 40900 185 1927 14 24 29 80
## 2 50100 197 2128 10 40 51 44
## 3 55900 231 10388 109 133 324 264
## 4 50300 222 11281 147 346 816 286
## 5 49500 237 3915 29 41 97 89
## 6 49500 217 2393 7 31 107 34
## UNITS50_UP MOBILEHOME NO_FARMS87 AVG_SIZE87 CROP_ACR87 AVG_SALE87
## 1 0 937 222 536 83787 27958
## 2 0 936 218 3489 29482 22155
## 3 0 3264 1073 490 131700 18138
## 4 0 3431 1476 907 144053 71970
## 5 0 1185 227 276 22923 10367
## 6 0 643 297 267 51806 29463
## Shape_Length Shape_Area Shape
## 1 370200.2 4620576797 MULTIPOLYGON (((49031.84 28...
## 2 361954.6 5905616645 MULTIPOLYGON (((-1704274 29...
## 3 453727.4 6552443738 MULTIPOLYGON (((-1598429 29...
## 4 582790.3 13742574263 MULTIPOLYGON (((-1713358 29...
## 5 298909.7 3742506837 MULTIPOLYGON (((-1574882 30...
## 6 254544.8 3313277935 MULTIPOLYGON (((-1549203 30...
# I've spent ~3 hours stuck on this one. I either get 679 (# of roads) or 23434 (# of cities). Whenever I try st_intersection I get 0 results. I also tried to merge the cities and county files via the state FIPS. I managed to do that and assign it geometry, but then whenever I tried st_intersection or intersects, I got a memory exceeded error and it crashed my computer each time and I had to restart it. Not sure what I'm doing wrong....
# Question 2
library(raster)
library(rgdal)
## rgdal: version: 1.5-18, (SVN revision 1082)
## Geospatial Data Abstraction Library extensions to R successfully loaded
## Loaded GDAL runtime: GDAL 3.0.4, released 2020/01/28
## Path to GDAL shared files: C:/Users/mwtro/OneDrive/Documents/R/win-library/4.0/rgdal/gdal
## GDAL binary built with GEOS: TRUE
## Loaded PROJ runtime: Rel. 6.3.1, February 10th, 2020, [PJ_VERSION: 631]
## Path to PROJ shared files: C:/Users/mwtro/OneDrive/Documents/R/win-library/4.0/rgdal/proj
## Linking to sp version:1.4-4
## To mute warnings of possible GDAL/OSR exportToProj4() degradation,
## use options("rgdal_show_exportToProj4_warnings"="none") before loading rgdal.
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:raster':
##
## intersect, select, union
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
elevation <- raster("final/ExamData/Terrain1.tif")
# cities_sf$elevation = raster::extract(elevation, cities_sf, buffer = 100)
# I ran this without the buffer and it yielded an elevation = 0 for everything, so I tried running it in with a buffer but processed for 30+ min so I decided to move on. What is next is what I would have tried to do
# statepop <- group_by(cities_sf, ID = state) %>%
# summarize_at(vars(elevation$ST), list(~mean(.), ~stdev(.)))
# Question 3
library(ggplot2)
library(RColorBrewer)
library(classInt)
states <- counties %>%
group_by(STATE_NAME)
statesyoung = states %>%
group_by(STATE_NAME) %>%
summarize(pop = sum(AGE_5_17))
## `summarise()` ungrouping output (override with `.groups` argument)
statesold = states %>%
group_by(STATE_NAME) %>%
summarize(pop = sum(AGE_65_UP))
## `summarise()` ungrouping output (override with `.groups` argument)
head(statesold)
## Simple feature collection with 6 features and 2 fields
## geometry type: GEOMETRY
## dimension: XY
## bbox: xmin: -5132790 ymin: 827627.6 xmax: 1044846 ymax: 6198811
## projected CRS: NAD83 / Conus Albers
## # A tibble: 6 x 3
## STATE_NAME pop Shape
## <chr> <int> <GEOMETRY [m]>
## 1 Alabama 522989 POLYGON ((720076.2 909241.3, 718504.1 922525.4, 714186 957~
## 2 Alaska 22369 MULTIPOLYGON (((-2422944 4041692, -2421227 4020023, -24279~
## 3 Arizona 478774 POLYGON ((-1365096 1022031, -1423154 1030895, -1449288 104~
## 4 Arkansas 350058 POLYGON ((462786.8 1241382, 466210.7 1238969, 469877.3 123~
## 5 California 3135552 MULTIPOLYGON (((-1709586 1320466, -1711344 1318502, -17187~
## 6 Colorado 329443 POLYGON ((-527954.4 1639437, -529558.1 1610882, -531964 15~
statesolddf <- st_drop_geometry(statesold)
statesage <- merge(statesyoung, statesolddf, by.x="STATE_NAME", by.y="STATE_NAME")
colnames(statesage) # popx = young people, popy = old people
## [1] "STATE_NAME" "pop.x" "pop.y" "geometry"
statesage$difference <- (statesage$pop.x-statesage$pop.y/statesage$pop.y)
colnames(statesage)
## [1] "STATE_NAME" "pop.x" "pop.y" "geometry" "difference"
ggplot(statesage) +
geom_sf(aes(fill=difference)) +
ggtitle("Difference Between Young (5-17) and Old (65)+ Populations")
