R Session Info
sessionInfo()
## R version 3.4.3 (2017-11-30)
## Platform: x86_64-apple-darwin15.6.0 (64-bit)
## Running under: macOS High Sierra 10.13.3
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/3.4/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/3.4/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## loaded via a namespace (and not attached):
## [1] compiler_3.4.3 backports_1.1.2 magrittr_1.5 rprojroot_1.3-2
## [5] tools_3.4.3 htmltools_0.3.6 yaml_2.1.17 Rcpp_0.12.15
## [9] stringi_1.1.6 rmarkdown_1.9 knitr_1.20 stringr_1.3.0
## [13] digest_0.6.15 evaluate_0.10.1
Draft in progress 2/7/18. Quick graphic exploration of 2015 diet data.
#require packages
library(googlesheets)
library(tidyverse)
library(forcats)
suppressMessages(library(dplyr))
# detach("package:MASS", unload=TRUE) #pkg masks some dplyr functions
# Clear environment
rm(list=ls())
###################################################################
# Create Stacked Bar Charts of Diet Proportions by Sampling Event #
#
# a.) Finer-level catgories
# b.) coarser-level categories
###################################################################
# identify gs sheet object
SCTC2015 <- gs_title("2015 EPSCoR SCTC.xlsx")
#read in worksheet containing diet data
diet <- SCTC2015 %>%
gs_read(ws = "2015 Diet Contents Data",col_names=TRUE) %>%
#rename needed columns
rename(sample.id=Sample.ID,
spp=Fish_Species,
river=River,
reach=Reach,
site=Site,
date=Date,
sampledate=Sample_Date,
taxkey=Taxon_Keyed,
prey=Prey_Type_Used,
preycat=PreyCategory,
dm_tot=Total_Prey_Dry_Mass_mg,
sample.event=Sample.Event,
event.num=sample.event.num,
TA=TA_Final) %>%
select(sample.id,sample.event,river,
reach,site,sampledate,preycat,spp,dm_tot,TA,event.num)
#filter out non-diet samples and blank data
diet <- diet %>% filter(site!="DR1",
!is.na(spp),
!is.na(sample.event),
!is.na(preycat),
!is.na(dm_tot),
!is.na(event.num)
)
#transform columns to desired classes
diet <- transform(diet,
dm_tot = as.numeric(dm_tot),
sample.event = as.factor(sample.event),
event.num = as.factor(event.num),
spp = as.factor(spp)) %>%
mutate(season = fct_recode(event.num, #code sampling events as seasons
"Early Summer"="1",
"Mid-Summer" = "2",
"Late Summer" ="3"))
#create sum of pooled prey DM for each prey category by fish species and sample event
dm_sum_by_event <- diet %>%
group_by(river,reach,sample.event,spp,preycat,TA,season) %>%
summarise(sum(dm_tot)) %>%
rename(dm_sum=`sum(dm_tot)`)
#create table of summed dry mass for each fish species by sampling event
dm_sum_spp <- dm_sum_by_event %>%
group_by(sample.event,spp) %>%
summarise(pool_sum = sum(dm_sum))
#join pooled prey DM values to table grouped by fish species and prey category
dm_sum_final <- dm_sum_spp %>%
left_join(dm_sum_by_event,by=c("sample.event","spp")) %>%
mutate(pct = (dm_sum/pool_sum)*100) %>%
filter(sample.event!="#N/A")
#create table of quantity of Chinook and coho non-empty diet samples from each SITE VISIT
diet_sample_count <- SCTC2015 %>%
gs_read(ws = "C 2015 Diet & LW Data",col_names=TRUE) %>%
filter(DS=="Y") %>%
mutate(season = fct_recode(as.factor(sample.event.num),
"Early Summer"="1",
"Mid-Summer" = "2",
"Late Summer" ="3")) %>%
group_by(sample.event,season,spp) %>%
summarise(diet_count=n())
#create table of quantity of Chinook and coho non-empty diet samples from each SEASON
diet_sample_season <- SCTC2015 %>%
gs_read(ws = "C 2015 Diet & LW Data",col_names=TRUE) %>%
filter(DS=="Y") %>%
mutate(season = fct_recode(as.factor(sample.event.num),
"Early Summer"="1",
"Mid-Summer" = "2",
"Late Summer" ="3")) %>%
group_by(river,season,spp) %>%
summarise(diet_count=n()) %>%
mutate(event_label = paste(as.character(spp), "( n=", as.character(diet_count),")"))
#join diet sample dry mass calculations by sampling event with fish count table by event
#e.g., how many fish of each species were diet-sampled at each sampling event?
#data is lumped by sampling event
dm_sum_final_by_event <- dm_sum_final %>%
left_join(diet_sample_count,
by=c("season","spp", "sample.event")) %>%
mutate(event_label = paste("n =", as.character(diet_count)))
#join diet sample dry mass calculations by sampling event with fish count table by season
#e.g., how many fish of each species were diet-sampled in each season?
#data is lumped for each drainage by `season`
dm_sum_final_by_season <- dm_sum_final_by_event %>%
group_by(season, river, spp, preycat, pct) %>%
summarise(dry_mass_pooled = sum(dm_sum))
# create objects that contain only one species per drainage, summarizing 2015
# diet data
### Coho ###
# results for Beaver Creek 2015 only, Coho only
dm_sum_final_beaver_coho_2015 <- dm_sum_final_by_event %>%
filter(river=="Beaver Creek", spp=="Coho")
# results for Russian River 2016 only, Coho only
dm_sum_final_russian_coho_2015 <- dm_sum_final_by_event %>%
filter(river=="Russian River", spp=="Coho")
# results for Ptarmigan Creek 2015 only, Coho only
dm_sum_final_ptarmigan_coho_2015 <- dm_sum_final_by_event %>%
filter(river=="Ptarmigan Creek", spp=="Coho")
# results for Kenai River 2015 only, Coho only
dm_sum_final_kenai_coho_2015 <- dm_sum_final_by_event %>%
filter(river=="Kenai River", spp=="Coho")
### Chinook ###
# results for Beaver Creek 2015 only, Chinook only
dm_sum_final_beaver_chinook_2015 <- dm_sum_final_by_event %>%
filter(river=="Beaver Creek", spp=="Chinook")
# results for Russian River 2016 only, Chinook only
dm_sum_final_russian_chinook_2015 <- dm_sum_final_by_event %>%
filter(river=="Russian River", spp=="Chinook")
# results for Ptarmigan Creek 2015 only, Chinook only
dm_sum_final_ptarmigan_chinook_2015 <- dm_sum_final_by_event %>%
filter(river=="Ptarmigan Creek", spp=="Chinook")
# results for Kenai River 2015 only, Chinook only
dm_sum_final_kenai_chinook_2015 <- dm_sum_final_by_event %>%
filter(river=="Kenai River", spp=="Chinook")
##################
# #
# Visualizations #
# #
##################
#######################################
# Chinook/Coho diet proportions 2015
#
# a.) FINE-level diet data (~5 prey categories)
#
# pooled by SITE and SEASON
#
########################################
# Create format for displaying fine-level diet data by site/season
finelevel15 <- ggplot() +
labs(x="",y="Percent Dry Mass (%)") +
scale_fill_discrete(name="Prey Category",
breaks=c("Fish",
"InvertAquatic",
"InvertTerrestrial_AqOrigin",
"InvertTerrestrial_TerOrigin",
"InvertTerrestrial_UnkOrigin",
"InvertUnknown",
"FishEggs"),
labels=c("Fish",
"Aquatic Invertebrate",
"Terrestrial Invertebrate,\nAquatic Origin",
"Terrestrial Invertebrate,\nTerrestrial Origin ",
"Terrestrial Invertebrate,\nUnknown Origin",
"Unknown Invertebrate",
"Fish Eggs")) +
theme(legend.key.size = unit(1.7, 'lines')) +
theme(plot.title = element_text(size = 24, face = "bold"),
axis.text.x = element_text(size=14),
axis.title.y = element_text(size = 14, face = "bold")) +
facet_grid(season~reach) +
theme(strip.text.x = element_text(size=10, face="bold"),
strip.text.y = element_text(size=10, face="bold"))
################
# Beaver Creek #
################
#### Coho ####
#Create faceted bar charts of dry mass proportions by species and sampling event
finelevel15 + geom_bar(aes(y=pct, x=spp, fill=preycat),
data = dm_sum_final_beaver_coho_2015,
stat="identity") +
geom_text(data=dm_sum_final_beaver_coho_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Beaver Creek 2015\n Juvenile Coho Diet (% Dry Mass)")
#### Chinook ####
#Create faceted bar charts of dry mass proportions by species and sampling event
finelevel15 + geom_bar(aes(y=pct, x=spp, fill=preycat),
data = dm_sum_final_beaver_chinook_2015,
stat="identity") +
geom_text(data=dm_sum_final_beaver_chinook_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Beaver Creek 2015\n Juvenile Chinook Diet (% Dry Mass)")
#################
# Russian River #
#################
#### Coho ####
#Create faceted bar charts of dry mass proportions by species and sampling event
finelevel15 + geom_bar(aes(y=pct, x=spp, fill=preycat),
data = dm_sum_final_russian_coho_2015,
stat="identity") +
geom_text(data=dm_sum_final_russian_coho_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Russian River 2015\n Juvenile Coho Diet (% Dry Mass)")
#### Chinook ####
#Create faceted bar charts of dry mass proportions by species and sampling event
finelevel15 + geom_bar(aes(y=pct, x=spp, fill=preycat),
data = dm_sum_final_russian_chinook_2015,
stat="identity") +
geom_text(data=dm_sum_final_russian_chinook_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Russian River 2015\n Juvenile Chinook Diet (% Dry Mass)")
###################
# Ptarmigan Creek #
###################
#### Coho ####
#Create faceted bar charts of dry mass proportions by species and sampling event
finelevel15 + geom_bar(aes(y=pct, x=spp, fill=preycat),
data = dm_sum_final_ptarmigan_coho_2015,
stat="identity") +
geom_text(data=dm_sum_final_ptarmigan_coho_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Ptarmigan Creek 2015\n Juvenile Coho Diet (% Dry Mass)")
#### Chinook ####
#Create faceted bar charts of dry mass proportions by species and sampling event
finelevel15 + geom_bar(aes(y=pct, x=spp, fill=preycat),
data = dm_sum_final_ptarmigan_chinook_2015,
stat="identity") +
geom_text(data=dm_sum_final_ptarmigan_chinook_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Ptarmigan Creek 2015\n Juvenile Chinook Diet (% Dry Mass)")
###############
# Kenai River #
###############
#### Coho ####
#Create faceted bar charts of dry mass proportions by species and sampling event
finelevel15 + geom_bar(aes(y=pct, x=spp, fill=preycat),
data = dm_sum_final_kenai_coho_2015,
stat="identity") +
geom_text(data=dm_sum_final_kenai_coho_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Kenai River 2015\n Juvenile Coho Diet (% Dry Mass)")
#### Chinook ####
#Create faceted bar charts of dry mass proportions by species and sampling event
finelevel15 + geom_bar(aes(y=pct, x=spp, fill=preycat),
data = dm_sum_final_kenai_chinook_2015,
stat="identity") +
geom_text(data=dm_sum_final_kenai_chinook_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Kenai River 2015\n Juvenile Chinook Diet (% Dry Mass)")
# ------------------------------------------------------------------------ #
#######################################
#
# b.) COARSE-level prey category data
#
# e.g., Terrestrial vs. Aquatic source
#
#######################################
# Create format for displaying diet data by site/season
coarselevel15 <- ggplot() +
labs(x="",y="Percent Dry Mass (%)") +
scale_fill_discrete(name="Prey Category",
breaks=c("T","A","M","UNK"),
labels=c("Terrestrial",
"Aquatic",
"Marine",
"Unknown")) +
theme(legend.key.size = unit(1.7, 'lines')) +
theme(plot.title = element_text(size = 24, face = "bold"),
axis.text.x = element_text(size=14),
axis.title.y = element_text(size = 14, face = "bold")) +
facet_grid(season~reach) +
theme(strip.text.x = element_text(size=10, face="bold"),
strip.text.y = element_text(size=10, face="bold"))
################
# Beaver Creek #
################
#### Coho ####
#Create faceted bar charts of dry mass proportions by species and sampling event
coarselevel15 + geom_bar(aes(y=pct, x=spp, fill=TA, order=TA),
data = dm_sum_final_beaver_coho_2015,
stat="identity") +
geom_text(data=dm_sum_final_beaver_coho_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Beaver Creek 2015\n Juvenile Coho Diet (% Dry Mass)")
#### Chinook ####
#Create faceted bar charts of dry mass proportions by species and sampling event
coarselevel15 + geom_bar(aes(y=pct, x=spp, fill=TA, order=TA),
data = dm_sum_final_beaver_chinook_2015,
stat="identity") +
geom_text(data=dm_sum_final_beaver_chinook_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Beaver Creek 2015\n Juvenile Chinook Diet (% Dry Mass)")
#################
# Russian River #
#################
#### Coho ####
#Create faceted bar charts of dry mass proportions by species and sampling event
coarselevel15 + geom_bar(aes(y=pct, x=spp, fill=TA, order=TA),
data = dm_sum_final_russian_coho_2015,
stat="identity") +
geom_text(data=dm_sum_final_russian_coho_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Russian River 2015\n Juvenile Coho Diet (% Dry Mass)")
#### Chinook ####
#Create faceted bar charts of dry mass proportions by species and sampling event
coarselevel15 + geom_bar(aes(y=pct, x=spp, fill=TA, order=TA),
data = dm_sum_final_russian_chinook_2015,
stat="identity") +
geom_text(data=dm_sum_final_russian_chinook_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Russian River 2015\n Juvenile Chinook Diet (% Dry Mass)")
###################
# Ptarmigan Creek #
###################
#### Coho ####
#Create faceted bar charts of dry mass proportions by species and sampling event
coarselevel15 + geom_bar(aes(y=pct, x=spp, fill=TA, order=TA),
data = dm_sum_final_ptarmigan_coho_2015,
stat="identity") +
geom_text(data=dm_sum_final_ptarmigan_coho_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Ptarmigan Creek 2015\n Juvenile Coho Diet (% Dry Mass)")
#### Chinook ####
#Create faceted bar charts of dry mass proportions by species and sampling event
coarselevel15 + geom_bar(aes(y=pct, x=spp, fill=TA, order=TA),
data = dm_sum_final_ptarmigan_chinook_2015,
stat="identity") +
geom_text(data=dm_sum_final_ptarmigan_chinook_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Ptarmigan Creek 2015\n Juvenile Chinook Diet (% Dry Mass)")
###############
# Kenai River #
###############
#### Coho ####
#Create faceted bar charts of dry mass proportions by species and sampling event
coarselevel15 + geom_bar(aes(y=pct, x=spp, fill=TA, order=TA),
data = dm_sum_final_kenai_coho_2015,
stat="identity") +
geom_text(data=dm_sum_final_kenai_coho_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Kenai River 2015\n Juvenile Coho Diet (% Dry Mass)")
#### Chinook ####
#Create faceted bar charts of dry mass proportions by species and sampling event
coarselevel15 + geom_bar(aes(y=pct, x=spp, fill=TA, order=TA),
data = dm_sum_final_kenai_chinook_2015,
stat="identity") +
geom_text(data=dm_sum_final_kenai_chinook_2015,
aes(x=1, y=4, label=event_label),
vjust = 0.5) +
ggtitle("Kenai River 2015\n Juvenile Chinook Diet (% Dry Mass)")
# ------------------------------------------------------------------------ #
#######################################
#
# c.) Diet mass per fish mass by season (very coarse gut "fullness" visualization)
#
# e.g., Terrestrial vs. Aquatic source
#
#######################################
# join T/A/M diet dry mass diet proportions of each fish to individual fish mass data.
# read in desired worksheet of diet data
diet_ta <- SCTC2015 %>%
gs_read(ws = "diet_T_A_M",col_names=TRUE) %>%
#rename needed columns
rename(sample.id=X1,
UNK=X2,
total_mass=`Grand Total`) %>%
# select needed columns
select(sample.event,sample.id,UNK,A,M,`T`,total_mass) %>%
#transform columns to desired classes
transform(sample.id=as.factor(sample.id)) %>%
#filter out NA site.id's
filter(!is.na(sample.id),
!is.na(sample.event),
!is.na(total_mass))
# read in fish length and weight data
lw <- SCTC2015 %>%
gs_read(ws = "C 2015 Diet & LW Data",col_names=TRUE) %>%
rename(wt=Weight_g,len=Length_mm) %>%
select(river,sample.event.num,sample.event,sample.id,wt,len,spp) %>%
filter(!is.na(sample.id),
!is.na(sample.event),
!is.na(wt))
# join diet data with fish length/weight data.
fnl <- left_join(lw,diet_ta) %>%
# Mutate new column for grams of food per grams of fish.
mutate(dm_per_fish = as.numeric(total_mass)/wt) %>%
filter(total_mass != 0) %>% #considering NON-EMPTY fish only
filter(sample.id != "LBC-F1-20150605-D11") %>% # this sample is a HUGE outlier; piscivory event
mutate(Season = fct_recode(as.factor(sample.event.num), #code sampling events as seasons
"Early Summer"="1",
"Mid-Summer" = "2",
"Late Summer" ="3"))
# STUCK ON ERROR 8/17/17
# Spent a long time here attempting to create a theme object that could easily be
# propograted throughout many plots for diet "fullness." When I try to create this
# theme object, something about the xlab, ylab, and ggtitle commands don't jive.
# Something about applying a non-binary command to non-binary object error repeatedly
# appears despite many rearrangements.
# for now, repeat code in each iteration... argh.
# create object for each drainage
fnl_beaver <- fnl %>%
filter(river == "Beaver Creek")
fnl_russian <- fnl %>%
filter(river == "Russian River")
fnl_ptarmigan <- fnl %>%
filter(river == "Ptarmigan Creek")
fnl_kenai <- fnl %>%
filter(river == "Kenai River")
# Plot total dry mass in diet per grams of fish for each event.
# 2015 beaver creek g diet mass per g fish
bc15 <- ggplot(fnl_beaver,aes(x=factor(Season),y=dm_per_fish,color=Season)) +
geom_jitter(size=4) +
theme_bw() +
theme(plot.title = element_text(size = 26, face = "bold"),
axis.text.x = element_text(size=14, face="bold"),
axis.text.y = element_text(size=16, face="bold"),
axis.title.x = element_text(size=22, face="bold"),
axis.title.y = element_text(size = 22, face = "bold")) +
theme(legend.title = element_text(size=30, face="bold")) +
theme(legend.text = element_text(size = 16, face = "bold")) +
ggtitle ("Beaver Creek 2015\nDiet Mass per Fish Mass ") +
xlab(" ") +
ylab ("mg Diet Mass / g Fish Mass")
# plot
bc15
# plot 2015 russian river g diet mass per g fish
rr15 <- ggplot(fnl_russian,aes(x=factor(Season),y=dm_per_fish,color=Season)) +
geom_jitter(size=4) +
theme_bw() +
theme(plot.title = element_text(size = 26, face = "bold"),
axis.text.x = element_text(size=14, face="bold"),
axis.text.y = element_text(size=16, face="bold"),
axis.title.x = element_text(size=22, face="bold"),
axis.title.y = element_text(size = 22, face = "bold")) +
theme(legend.title = element_text(size=30, face="bold")) +
theme(legend.text = element_text(size = 16, face = "bold")) +
ggtitle ("Russian River 2015\nDiet Mass per Fish Mass ") +
xlab(" ") +
ylab ("mg Diet Mass / g Fish Mass")
# plot
rr15
# plot 2015 ptarmigan creek g diet mass per g fish
pc15 <- ggplot(fnl_ptarmigan,aes(x=factor(Season),y=dm_per_fish,color=Season)) +
geom_jitter(size=4) +
theme_bw() +
theme(plot.title = element_text(size = 26, face = "bold"),
axis.text.x = element_text(size=14, face="bold"),
axis.text.y = element_text(size=16, face="bold"),
axis.title.x = element_text(size=22, face="bold"),
axis.title.y = element_text(size = 22, face = "bold")) +
theme(legend.title = element_text(size=30, face="bold")) +
theme(legend.text = element_text(size = 16, face = "bold")) +
ggtitle ("Ptarmigan Creek 2015\nDiet Mass per Fish Mass ") +
xlab(" ") +
ylab ("mg Diet Mass / g Fish Mass")
# plot
pc15
# plot 2015 kenai river g diet mass per g fish
kr15 <- ggplot(fnl_kenai,aes(x=factor(Season),y=dm_per_fish,color=Season)) +
geom_jitter(size=4) +
theme_bw() +
theme(plot.title = element_text(size = 26, face = "bold"),
axis.text.x = element_text(size=14, face="bold"),
axis.text.y = element_text(size=16, face="bold"),
axis.title.x = element_text(size=22, face="bold"),
axis.title.y = element_text(size = 22, face = "bold")) +
theme(legend.title = element_text(size=30, face="bold")) +
theme(legend.text = element_text(size = 16, face = "bold")) +
ggtitle ("Kenai River 2015\nDiet Mass per Fish Mass ") +
xlab(" ") +
ylab ("mg Diet Mass / g Fish Mass")
# plot
kr15
# To do:
#- what % of diet items are calculted from etsimated lengths?
#-are there particular paces/times where UNID Parts dominate sample composition?