This is my Assignment 3 RMD file! Name: Lyndsay Cuff

#load in libraries
library(terra)

## terra 1.7.71

library(tidyterra)

## 
## Attaching package: 'tidyterra'

## The following object is masked from 'package:stats':
## 
##     filter

library(sp)
library(ggplot2)
library(ggthemes)
library(urbnmapr)
library(sf)

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

library(leaflet)
library(dplyr)   

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:terra':
## 
##     intersect, union

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(spData) 

## To access larger datasets in this package, install the spDataLarge
## package with: `install.packages('spDataLarge',
## repos='https://nowosad.github.io/drat/', type='source')`

library(tidyr)

## 
## Attaching package: 'tidyr'

## The following object is masked from 'package:terra':
## 
##     extract

#load in the data first and foremost

setwd("/Users/lyndsaycuff/DIDA 370")
foster2014 <- read.csv("countyfoster2014.csv")

setwd("/Users/lyndsaycuff/DIDA 370/NE1_50M_SR_W")
rastuh <- rast("NE1_50M_SR_W.tif")
rastuh <- terra::project(rastuh, "EPSG:4326", method = "bilinear")

## 
|---------|---------|---------|---------|
=========================================
                                          

#grab states shape files for our map
states <- get_urbn_map("states", sf = TRUE)
#change projection while we're at it
states <- states %>% st_transform("EPSG:4326")

## old-style crs object detected; please recreate object with a recent sf::st_crs()

#load in some more data...I am not waiting all that time AGAIN

setwd("/Users/lyndsaycuff/DIDA 370")
foster2014 <- read.csv("nycountyfoster2014.csv")

#lets standardize the data and make some RATES
foster2014 <- foster2014 %>% mutate(percent_foster = round((Foster.Family.Homes/Children.in.Care)*100, 2),
                                    percent_institution = round((Institutions/Children.in.Care)*100, 2), 
                                    percent_grouphome = round((Group.Homes/Children.in.Care)*100, 2))

#lets make the raster map

#first have to crop the raster so that it is the same size as the entire US
#should be able to do this with one line SINCE they are projected the same way
cropped_raster <- crop(rastuh, states)

#now lets plot it all together
ggplot()+  
  geom_spatraster(data = cropped_raster)+
  geom_sf(states, mapping = aes(), fill = "transparent", color = "black")+
  theme_map()+
  scale_fill_whitebox_c(
    palette = "atlas",
    na.value = "white")+
  labs(title = "Landcover Across US States",
       fill = "Landcover")+
  theme(plot.title = element_text(hjust = 0.5, size = 16),
        legend.position="right")

## <SpatRaster> resampled to 500636 cells for plotting
## ! `tidyterra::geom_spatraster()`: Plotting 3 overlapping layers: NE1_50M_SR_W_1, NE1_50M_SR_W_2, and NE1_50M_SR_W_3. Either:
## 
##   Use `facet_wrap(~lyr)` for faceting or
## 
##   Use `aes(fill = <name_of_layer>)` for displaying single layers

#I am gonna make the other map

#grab counties BASE MAP
counties <- get_urbn_map("counties", sf = TRUE)

#filter the data to get just NYS
nycounties <- counties %>% 
  filter(state_abbv == "NY") %>% 
  st_transform("EPSG:4326")

## old-style crs object detected; please recreate object with a recent sf::st_crs()

#MERGE BABY MERGE
#remove the unnecessary words first... we get it you're a county
nycounties <- nycounties %>% 
  separate(county_name, c("county_name", "county"), sep = " County") %>% 
  select(-county)

#Now merge it!!
full_data2014 <- nycounties %>%
  left_join(foster2014, by = c("county_name" = "Social.Service.District"))

#lets make a plot tho
#start with the basics

#start with a popup
#I used this reddit forum to figure out how to put the popup label on two separate lines (using the <br>)
#https://www.reddit.com/r/Rlanguage/comments/6bsji1/add_line_break_to_leaflet_pop_up/ 

popupfoster <- paste0("<strong>Community: </strong>", full_data2014$county_name, "<br>", "<strong>Percentage of Children in Foster Family Homes: </strong>", full_data2014$percent_foster)

popupgrouphomes <- paste0("<strong>Community: </strong>", full_data2014$county_name, "<br>", "<strong>Percentage of Children in Group Homes: </strong>", full_data2014$percent_grouphome)

#and then a palette
palettefoster <- colorBin("Purples", full_data2014$percent_foster, bins = 5, pretty = FALSE, na.color = "black")

palettegrouphome <- colorBin("Blues", full_data2014$percent_grouphome, bins = 5, pretty = FALSE, na.color = "black")

#I guess we plot now
leaflet() %>%
  leaflet::addPolygons(data = full_data2014,
                       stroke = FALSE,
                       fillColor = ~palettefoster(full_data2014$percent_foster), 
                       fillOpacity = 1, 
                       smoothFactor = 0.5, 
                       popup = popupfoster, 
                       group = "Foster Homes") %>%
  leaflet::addPolygons(data = full_data2014, 
                       stroke = FALSE, 
                       fillColor = ~palettegrouphome(full_data2014$percent_grouphome), 
                       fillOpacity = 1, 
                       smoothFactor = 0.5, 
                       popup = popupgrouphomes, 
                       group = "Group Homes") %>%
  addLegend("bottomright",
            pal=palettefoster,
            values=full_data2014$percent_foster, 
            title = "Percentage of Children in Foster Care", 
            group = 'Foster Homes')  %>%
  addLegend("bottomleft", 
            pal = palettegrouphome, 
            values=full_data2014$percent_grouphome, 
            title = "Percentage of Children in Group Homes", 
            group = "Group Homes") %>%
    addProviderTiles(providers$OpenStreetMap) %>%
    addLayersControl(overlayGroups = c("Foster Homes", "Group Homes"),
                   options = layersControlOptions(collapsed = FALSE))