# LIVE CODING SCRIPT ------------------------------------------------------

# Website: https://sccwrp.github.io/CABW2018_R_training/index.html
# Live Coding Link: https://tinyurl.com/cabw-r-2018

# Data: https://sccwrp.github.io/CABW2018_R_training/data/datazip.zip


# here's a comment that is describing my code
print('hello world')

seq(1, 10)

rnorm(10)

# saving rnadom variables
a <- rnorm(10000, mean = 5, sd = 1)
mean(a)

# # install packages from CRAN
# install.packages("tidyverse")
# install.packages("sf")
# install.packages("mapview")
# install.packages("viridis")
# install.packages("USAboundaries")

library("tidyverse")
library("sf")
library("mapview")
library("viridis")
library("USAboundaries")

dbl_var <- c(1, 2.5, 4.5)
int_var <- c(1L, 6L, 10L)
log_var <- c(TRUE, FALSE, T, F)
chr_var <- c("a", "b", "c")

length(dbl_var)

ltrs <- c('a', 'b', 'c')
nums <- c(1, 2, 3)
logs <- c(T, F, T)
mydf <- data.frame(ltrs, nums, logs)
mydf

cscidat <- read.csv('data/cscidat.csv', stringsAsFactors = F)
ascidat <- read.csv('data/ascidat.csv', stringsAsFactors = F)

tmp <- select(cscidat, SampleID_old, New_Lat, New_Long, CSCI)

tmp2 <- select(cscidat, -E, -OE, -pMMI)

tmp3 <- select(cscidat, matches('c'))

tmp4 <- filter(cscidat, CSCI > 0.79)

tmp5 <- filter(cscidat, CSCI <= 0.79)

tmp6 <- filter(cscidat, !is.na(CSCI))

# check for NA's

summary(cscidat)

dplyr_mut1 <- mutate(cscidat, observed = OE * E)

newcolumn <- seq(1, nrow(cscidat))
dply_mut2 <- mutate(cscidat, newcolumn = newcolumn)

dply_mut3 <- mutate(cscidat, cscicat = ifelse(CSCI > 0.79, 'hi', 'lo'))

dplyr_rnm <- rename(cscidat, 
                    lat = New_Lat,
                    lon = New_Long
)

ggplot(data = ascidat) +
  geom_boxplot(mapping = aes(x = site_type, y = ASCI))


# Making Maps -------------------------------------------------------------

library(tidyverse)
library(sf)
library(mapview)
library(viridis)
library(USAboundaries)



# Load the Data -----------------------------------------------------------


# if reading locally (from your data folder in you RStudio project):
ascidat <- read_csv("data/ascidat.csv")
cscidat <- read_csv("data/cscidat.csv")
latlons <- read_csv("data/latlon.csv")

summary(cscidat)



# Make Data Spatial -------------------------------------------------------

# make data sf object:


df_sf <- st_as_sf(latlons,
                  # coords = c("Lon", "Lat"), # can use numbers here too
                  coords = c(3, 2), # can use numbers here too
                  remove = F, # don't remove these lat/lon cols from df
                  crs = 4326) # add projection (this is WGS84)

class(df_sf)
class(latlons)


st_crs(df_sf)
df_utms <- st_transform(df_sf, crs=3310)
st_crs(df_utms)



# Make a Map --------------------------------------------------------------



plot(df_sf)

plot(df_sf$geometry)

plot(df_sf$geometry, col = "dodgerblue")


# make an interactive map using the mapview package

mapview(df_sf)


# Spatial Stuff -----------------------------------------------------------


library(USAboundaries)

# Pick a State
state_names <- c("california") # notice we can add more states to this list if we want

# Download STATE data and add projection
CA<-us_states(resolution = "high", states = state_names) # what does resolution do?

st_crs(CA) # double check the CRS


plot(CA$geometry)
plot(df_sf$geometry, add=TRUE, col="orange")


# ggplot map
ggplot() + # can include data here, if only using one dataset
  geom_sf(data=CA, color="gray") +
  geom_sf(data=df_sf, aes(color=New_Long)) +
  labs(x="Longitude", y="Latitude", title= "Map of stuff")

ggsave(filename = "figures/map_of_cscsi_points.png", width = 8.5, height = 11, units = "in", dpi = 200)



# CROP --------------------------------------------------------------------

# Download outlines of all the counties for California
CA_counties<-us_counties(resolution = "high", states = "CA") # what does resolution do?
plot(CA_counties$geometry)

# Pick a county of interest
co_names <- c("Sacramento") # notice we can add more states to this list if we want

# filter to just that county:
sacto_co <- filter(CA_counties, name==co_names)

# quick map
ggplot() +
  geom_sf(data=CA_counties, color="gray")+
  geom_sf(data=sacto_co, color="purple")


sacto_pts <- st_intersection(df_sf, sacto_co)

# ggplot of cropped data
ggplot() +
  geom_sf(data=sacto_co, color="gray")+
  geom_sf(data=sacto_pts, fill="orange", pch=21) +
  theme_bw()



# JOIN WITH CSCI DATA -----------------------------------------------------

sacto_csci <- left_join(sacto_pts, cscidat, by=c("StationID"="StationCode"))

class(sacto_csci)


mapview(sacto_csci, zcol="CSCI", layer="CSCI") + mapview(sacto_co)




# SAVE OUT SPATIAL DATA ---------------------------------------------------

st_write(obj = sacto_csci, "data/sacramento_county_csci.shp", delete_layer = TRUE)


test_shp <- st_read("data/sacramento_county_csci.shp")