aa
aa=PISCO.GBY %>% select(site,side) %>% unique()
aa
bb= subset(PISCO.GBY,site=='LOVERS')
View(bb)
View(aa)
aa=PISCO.GBY %>% select(site,sideaa) %>% unique()
aa
aa
View(bb)
bb = subset(PISCO.GBY, grep(side,"OLD."))
?grep
bb = subset(PISCO.GBY, side == grep("OLD."))
library(stringr)
#look at all the "old transect"
PISCO.old = PISCO.GBY %>%
filter(str_detect(side, '^OLD'))
View(PISCO.old)
ANDREW_MOLERA.sub = PISCO.GBY %>%
filter(site=='ANDREW_MOLERA')
with(ANDREW_MOLERA.sub, table(side,year))
LOPEZ.sub = PISCO.GBY %>%
filter(site=='LOPEZ')
with(LOPEZ.sub, table(side,year))
##Drop sites with "old" sides transect names
PISCO.GBY1 = PISCO.GBY %>%
filter(!str_detect(side, '^OLD'))
summary(PISCO.GBY1$side)
PISCO.GBY = droplevels(PISCO.GBY)
aa = subset(PISCO.GBY, site=='BUTTERFLY')
View(aa)
#read in pisco sites
PISCO.sites = read.csv('C:/Users/melissa.monk/Desktop/PISCO_sites.csv')
PISCO.site.transect = inner_join(PISCO.GBY,PISCO.sites)
#look at all the "old transect"
PISCO.old = PISCO.GBY %>%
filter(str_detect(side, '^OLD'))
#check "old" sides for two sites
ANDREW_MOLERA.sub = PISCO.GBY %>%
filter(site=='ANDREW_MOLERA')
with(ANDREW_MOLERA.sub, table(side,year))
LOPEZ.sub = PISCO.GBY %>%
filter(site=='LOPEZ')
with(LOPEZ.sub, table(side,year))
##Drop sites with "old" sides transect names
PISCO.GBY = PISCO.GBY %>%
filter(!str_detect(side, '^OLD'))
PISCO.GBY = droplevels(PISCO.GBY)
PISCO.GBY = inner_join(PISCO.GBY,PISCO.sites)
##Drop sites with "old" sides transect names
PISCO.GBY = PISCO.GBY %>%
filter(!str_detect(side, '^OLD'))
PISCO.GBY = droplevels(PISCO.GBY)
#Look at inside/outside for an mpa
PISCO.mpa = PISCO.GBY %>%
group_by(MPAGroup,RESERVE) %>%
tally()
View(PISCO.mpa)
aa = subset(PISCO, site=='CAYUCOS')
View(aa)
#Positives
GBY.pos = subset(PISCO.GBY, GBY>0)
dim(PISCO.GBY)
dim(GBY.pos)
3843/11382
#tally how many transects at a single site in a day
#grouping over side, zone and transect
PISCO.daily.transect = PISCO.GBY %>%
group_by(campus,method,year,month,dat,site) %>%
tally()
#tally how many transects at a single site in a day
#grouping over side, zone and transect
PISCO.daily.transect = PISCO.GBY %>%
group_by(campus,method,year,month,day,site) %>%
tally()
View(PISCO.daily.transect)
#tally how many transects at a single site in a day
#grouping over side, zone and transect
PISCO.daily.transect = PISCO.GBY %>%
group_by(campus,method,year,month,day,site,RESERVE) %>%
tally()
View(PISCO.GBY)
#tally how many transects at a single site in a day
#grouping over side, zone and transect
PISCO.daily.transect = PISCO.GBY %>%
group_by(campus,method,year,month,day,site,side,RESERVE) %>%
tally()
View(PISCO.daily.transect)
#See if gby never appear in some combos
PISCO.GBY.transect = PISCO.GBY %>%
group_by(campus,method,site,side,RESERVE) %>%
summarise(GBYsum = sum(GBY))
View(PISCO.GBY.transect)
aa = subset(PISCO.GBY, site=='ANACAPA_WEST_ISLE' & side=='E')
View(aa)
library(dplyr)
library(tidyr)
library(tibble)
library(ggplot2)
library(stringr)
PISCO = read.csv('Z:/Gopher_BlackandYellow_2019/Indices_of_Abundance/PISCO/PISCO_FISH.2018.csv')
#classcode is species column
#gopher SCAR
#b&y SCHR
###Initial sample sizes
#take a look at where the black and yellow and gophers are
GBY = subset(PISCO, classcode %in% c('SCAR','SCHR'))
summary(GBY$level)  #only really in the bottom level - can subset all data to this
GBY = subset(GBY, level=='BOT')
summary(as.factor(GBY$month)) #probably remove june
summary(as.factor(GBY$year)) #keep 1999?  96 positives, probably ok
summary(GBY$campus)
summary(GBY$method)  #all the same
#come back to this and find site sampled at least x number of times, with gophers
#this df doesn't comine multiple fish on a transect, need to do this first
summary(GBY$site)
with(GBY, table(site, campus))
#side...many blank...need look up table
summary(GBY$side)
#General sample, not yet filtered to 1 row per transect
#summaries to check for missing essential data
summary(as.factor(PISCO$month))
summary(as.factor(PISCO$year))  #1 NA
summary(as.factor(PISCO$day))
summary(as.factor(PISCO$zone))
summary(as.factor(PISCO$site))
summary(as.factor(PISCO$vis)) #lots of NA's
summary(as.factor(PISCO$surge))  #lots of blanks
summary(as.factor(PISCO$depth)) #lots of NA's
#Try and find unique "transects" and only keep bottom transects
PISCO = PISCO %>%
filter(!is.na(year), level=='BOT')
PISCO = droplevels(PISCO)
PISCO.transect = PISCO %>%
select(campus, method, year, month, day, site, side, zone,
level, transect) %>%
unique() %>%
mutate(id = row_number())   #adds a row number as a trip identifier
#Merge in the row_number to the original dataframe
PISCO.keep = right_join(PISCO, PISCO.transect)
#check join
summary(PISCO.transect$id)
summary(PISCO.keep$id)
#Get fish by transect - keep all just in case
PISCO.fish  = PISCO.keep %>%
select(id, classcode, count) %>%  #keep just transect id, specides and count
group_by(id, classcode) %>%       #group by the id and speice and summ counts
summarise(sum_cnt = sum(count)) %>%
spread(classcode, sum_cnt) %>%
replace(is.na(.), 0)
#Now get the gophers plus black and yellows and keep only relevent columns
PISCO.GBYcnts = PISCO.fish %>%
mutate(GBY= SCAR + SCHR) %>%
select(id, GBY)
#Merge back in with PISCO.transect
PISCO.GBY = inner_join(PISCO.transect, PISCO.GBYcnts)
GBY.pos.site = PISCO.GBY %>%
filter(GBY>0) %>%   #get positives and look at which transects ever saw them
group_by(site) %>%
tally()
#Filter to reefs that have seen GBY at least once
PISCO.GBY = subset(PISCO.GBY, site %in% GBY.pos.site$site)
#look at % GBY by zone
with(subset(PISCO.GBY, GBY>0), table(zone))/with(PISCO.GBY, table(zone))
#look at % GBY by zone
with(subset(PISCO.GBY, GBY>0), table(zone))/with(PISCO.GBY, table(zone))
#not really in the mid zone and not well campled
PISCO.GBY = filter(PISCO.GBY,zone !='MID')
PISCO.GBY = droplevels(PISCO.GBY)
#Look at year and month
GBY.pos = subset(PISCO.GBY, GBY>0)
summary(as.factor(PISCO.GBY$month))
summary(as.factor(GBY.pos$month))
summary(as.factor(GBY.pos$month))/summary(as.factor(PISCO.GBY$month))
summary(as.factor(PISCO.GBY$year))
summary(as.factor(GBY.pos$year))
summary(as.factor(GBY.pos$year))/summary(as.factor(PISCO.GBY$year))
#Drop June
PISCO.GBY = filter(PISCO.GBY, month !=6)
#look to see reefs seen in at least 1/2 of all years - too strict
#either need to combine sites or figure out another way
GBY.site.years = PISCO.GBY %>%
group_by(site) %>%
summarise(yearcount = n_distinct(year),
n = n()) %>%
filter(yearcount>9)   #Filter to sites sampled in at least 1/2 of the years
#Subset to sites sampled these years
PISCO.GBY = subset(PISCO.GBY, site %in% GBY.site.years$site)
#how many positives left
GBY.pos = subset(PISCO.GBY, GBY>0)
#look at how many positives by number of years sampled and site
GBY.pos.site.year = GBY.pos %>%
group_by(site) %>%
summarise(yearcount = n_distinct(year),
n = n())
#Look at what's left by site
GBY.site.year = PISCO.GBY %>%
group_by(site) %>%
summarise(yearcount = n_distinct(year),
n = n())
#a good number of transects left
#look campus to get idea of north and south
summary(PISCO.GBY$campus) #a lot remains from UCSB
summary(GBY.pos$campus)  #not a lot of gophers left in southern CA
#Need lat/long for each site - found and added in below
#Are all the UCSB sites south of Conception??
#Model separately?
PISCO.UCSB = subset(PISCO.GBY, campus=='UCSB')
PISCO.UCSB = droplevels(PISCO.UCSB)
GBY.pos.UCSB = subset(PISCO.UCSB, GBY>0)
summary(PISCO.UCSB$site)
summary(GBY.pos.UCSB$site)
summary(GBY.pos.UCSB$site)/summary(PISCO.UCSB$site)
#UCSC
PISCO.UCSC = subset(PISCO.GBY, campus=='UCSC')
PISCO.UCSC = droplevels(PISCO.UCSC)
GBY.pos.UCSC = subset(PISCO.UCSC, GBY>0)
summary(PISCO.UCSC$site)
summary(GBY.pos.UCSC$site)
summary(GBY.pos.UCSC$site)/summary(PISCO.UCSC$site)
summary(GBY.pos.UCSB$site)/summary(PISCO.UCSB$site)
###More overall filtering
GBY.pos.site.20 = subset(GBY.pos.site, n>19)
View(GBY.pos.site.20)
View(GBY.pos.site.year)
###More overall filtering
GBY.pos.site.20 = subset(GBY.pos.site, n>19)  #subset to site that saw at least 20 gophers over time
PISCO.GBY = subset(PISCO.GBY, site %in% GBY.pos.site.20$site)
GBY.pos = subset(PISCO.GBY, GBY>0)
PISCO.GBY = droplevels(PISCO.GBY)
summary(PISCO.GBY$site)
#read in pisco sites
PISCO.sites = read.csv('C:/Users/melissa.monk/Desktop/PISCO_sites.csv')
PISCO.GBY = inner_join(PISCO.GBY,PISCO.sites)
#look at all the "old transect"
PISCO.old = PISCO.GBY %>%
filter(str_detect(side, '^OLD'))
#check "old" sides for two sites
ANDREW_MOLERA.sub = PISCO.GBY %>%
filter(site=='ANDREW_MOLERA')
with(ANDREW_MOLERA.sub, table(side,year))
LOPEZ.sub = PISCO.GBY %>%
filter(site=='LOPEZ')
with(LOPEZ.sub, table(side,year))
##Drop sites with "old" sides transect names
PISCO.GBY = PISCO.GBY %>%
filter(!str_detect(side, '^OLD'))
PISCO.GBY = droplevels(PISCO.GBY)
#Look at inside/outside for an mpa
PISCO.mpa = PISCO.GBY %>%
group_by(MPAGroup,RESERVE) %>%
tally()
View(PISCO.mpa)
#Positives
GBY.pos = subset(PISCO.GBY, GBY>0)  #3842 out of 11382 transects
#tally how many transects at a single site in a day
#grouping over side, zone and transect
PISCO.daily.transect = PISCO.GBY %>%
group_by(campus,method,year,month,day,site,side,RESERVE) %>%
tally()
#See if gby never appear in some combos
PISCO.GBY.transect = PISCO.GBY %>%
group_by(campus,method,site,side,RESERVE) %>%
summarise(GBYsum = sum(GBY))
View(PISCO.GBY.transect)
#See if gby never appear in some combos
PISCO.GBY.transect = PISCO.GBY %>%
group_by(campus,method,site,RESERVE) %>%
summarise(GBYsum = sum(GBY))
#See if gby never appear in some combos
PISCO.GBY.transect1 = PISCO.GBY %>%
group_by(campus,method,site,level,RESERVE) %>%
summarise(GBYsum = sum(GBY))
View(PISCO.GBY.transect1)
View(PISCO.GBY)
#See if gby never appear in some combos
PISCO.GBY.transect1 = PISCO.GBY %>%
group_by(campus,method,site,zone,RESERVE) %>%
summarise(GBYsum = sum(GBY))
View(PISCO.GBY.transect1)
#See if gby never appear in some combos
PISCO.GBY.transect1 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYsum = sum(GBY))
#See if gby never appear in some combos of variables
PISCO.GBY.transect1 = PISCO.GBY %>%
group_by(campus,site,zone,RESERVE) %>%
summarise(GBYsum = sum(GBY))
#See if gby never appear in some combos
PISCO.GBY.transect2 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYsum = sum(GBY))
View(PISCO.GBY.transect2)
#See if gby never appear in some combos
PISCO.GBY.transect2 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYsum = sum(GBY),
n = n_distinct())
#See if gby never appear in some combos
PISCO.GBY.transect2 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYsum = sum(GBY),
total_transect = n())
summary(PISCO.GBY$GBY)
aa = subset(PISCO.GBY, GBY>10)
View(aa)
bb = subset(PISCO, site=='TERRACE' & zone=='OUTER' & year==2003 & month==8 & day==28)
View(bb)
bb = subset(PISCO, site=='TERRACE' & zone=='OUTER' & year==2003 & month==8 & day==28 & classcode=='SCAR')
#See if gby never appear in some combos
PISCO.GBY.transect2 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYtransect = n_distinct(GBY>0),
total_transect = n())
#See if gby never appear in some combos
PISCO.GBY.transect2 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYtransect = count(GBY>0),
total_transect = n())
#See if gby never appear in some combos
PISCO.GBY.transect2 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYtransect = tally(GBY>0),
total_transect = n())
#See if gby never appear in some combos
PISCO.GBY.transect2 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYtransect = length(GBY>0),
total_transect = n())
View(PISCO.GBY)
#See if gby never appear in some combos
PISCO.GBY.transect2 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYtransect = length(value[GBY>0]),
total_transect = n())
#See if gby never appear in some combos
PISCO.GBY.transect2 = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(GBYtransect = length(value[GBY>0]),
total_transect = n())
#See if gby never appear in some combos
PISCO.GBY.transect2 = GBY.pos %>%
group_by(campus,zone) %>%
summarise(total_transect = n())
#See if gby never appear in some combos
PISCO.GBY.transect2a = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(total_transect = n())
View(PISCO.GBY.transect2a)
GBY.zone = inner_join(PISCO.GBY.transect2, PISCO.GBY.transect2a)
GBY.zone = cbind(PISCO.GBY.transect2, PISCO.GBY.transect2a)
View(GBY.zone)
GBY.zone.percent=GBY.zone$totalpos_transect/GBY.zone$total_transect
#See if gby never appear in some combos
PISCO.GBY.transect2 = GBY.pos %>%
group_by(campus,zone) %>%
summarise(totalpos_transect = n())
#See if gby never appear in some combos
PISCO.GBY.transect2a = PISCO.GBY %>%
group_by(campus,zone) %>%
summarise(total_transect = n())
GBY.zone = cbind(PISCO.GBY.transect2, PISCO.GBY.transect2a)
GBY.zone.percent = GBY.zone$totalpos_transect/GBY.zone$total_transect
GBY.zone$percent = GBY.zone$totalpos_transect/GBY.zone$total_transect
#See if gby never appear in some combos
PISCO.GBY.transect = PISCO.GBY %>%
group_by(campus,site,RESERVE) %>%
summarise(GBYsum = sum(GBY))
View(PISCO.GBY.transect)
#See if gby never appear in some combos
PISCO.GBY.transect = PISCO.GBY %>%
group_by(campus,site,RESERVE) %>%
summarise(GBYsum = sum(GBY))
#See if gby never appear in some combos of variables
PISCO.GBY.transect1 = PISCO.GBY %>%
group_by(campus,site,zone,RESERVE) %>%
summarise(GBYsum = sum(GBY))
View(PISCO.GBY.transect1)
#See if gby never appear in some combos of variables
PISCO.GBY.transect1 = GBY.pos %>%
group_by(campus,site,zone,RESERVE) %>%
summarise(GBYsum = n())
PISCO.GBY.transect1a = PISCO.GBY %>%
group_by(campus,site,zone,RESERVE) %>%
summarise(total = n())
View(PISCO.GBY.transect1a)
View(PISCO.GBY.transect1)
View(PISCO.GBY.transect1a)
PISCO.GBY.zone = inner_join(PISCO.GBY.transect1,PISCO.GBY.transect1a)
View(PISCO.GBY.zone)
PISCO.GBY.zone = outer_join(PISCO.GBY.transect1,PISCO.GBY.transect1a)
PISCO.GBY.zone = smie_join(PISCO.GBY.transect1,PISCO.GBY.transect1a)
PISCO.GBY.zone = semi_join(PISCO.GBY.transect1,PISCO.GBY.transect1a)
?inner_join
PISCO.GBY.zone = full_join(PISCO.GBY.transect1,PISCO.GBY.transect1a)
View(PISCO.GBY.zone)
#See if gby never appear in some combos of variables
PISCO.GBY.transect1 = GBY.pos %>%
group_by(campus,site,zone,RESERVE) %>%
summarise(GBYtotal = n())
PISCO.GBY.transect1a = PISCO.GBY %>%
group_by(campus,site,zone,RESERVE) %>%
summarise(total = n())
PISCO.GBY.zone = full_join(PISCO.GBY.transect1,PISCO.GBY.transect1a)
PISCO.GBY.zone = PISCO.GBY.zone$GBYtotal/PISCO.GBY.zone$total
PISCO.GBY.zone = full_join(PISCO.GBY.transect1,PISCO.GBY.transect1a)
PISCO.GBY.zone$percent = PISCO.GBY.zone$GBYtotal/PISCO.GBY.zone$total
PISCO.GBY.transect1b = PISCO.GBY %>%
group_by(campus,site,zone,RESERVE,year) %>%
summarise(total = n()) %>%
spread(year,total)
View(PISCO.GBY.transect1b)
getwd()
write.csv(PISCO.GBY.transect1b,'Z:/Gopher_BlackandYellow_2019/Indices_of_Abundance/PISCO/transect_by_zone_year')
write.csv(PISCO.GBY.transect1b,'Z:/Gopher_BlackandYellow_2019/Indices_of_Abundance/PISCO/transect_by_zone_year.csv')
PISCO.GBY.transect1c = GBY.pos %>%
group_by(campus,site,zone,RESERVE,year) %>%
summarise(GBYtotal = n()) %>%
spread(year,GBYtotal=n())
PISCO.GBY.transect1c = GBY.pos %>%
PISCO.GBY.transect1c = GBY.pos %>%
group_by(campus,site,zone,RESERVE,year) %>%
summarise(GBYtotal = n()) %>%
spread(year,GBYtotal)
PISCO.GBY.transect1c = GBY.pos %>%
group_by(campus,site,zone,RESERVE,year) %>%
summarise(GBYtotal = n()) %>%
spread(year,GBYtotal)
View(PISCO.GBY.transect1c)
View(PISCO.GBY.transect1c)
write.csv(PISCO.GBY.transect1b,'Z:/Gopher_BlackandYellow_2019/Indices_of_Abundance/PISCO/GBYtransect_by_zone_year.csv')
write.csv(PISCO.GBY.transect1c,'Z:/Gopher_BlackandYellow_2019/Indices_of_Abundance/PISCO/GBYtransect_by_zone_year.csv')
View(PISCO.GBY)
View(GBY.pos)
60*24
60+102
160+160
229/12500
install.packages("rmarkdown")
install.packages("tinytex")
25*.7
3860+50
20009+101431+47652
20009+101431+47652+75948
20009+101431+47652+75948+1558
20009+101431+47652+75948+1558+1582
50907
50907-34752
3910+13305
193585-191335
224247-218000
1804+3243
1804+3243+1032
1804+3243+1032+784
34752-29060
(696+1444+252)*2
4784+413+413
5610+696+696
5692-4413
5692-1804-3243
3242+1804
3243+1804
1972+50
972141+50907
109276-1023048
1009276-1023048
991832-991781
1023048-1009276
13772-13305
13772-13398
3485-2373
3485-2373-700
974524+34752
50+13305
50+13398
setwd("C:/Users/melissa.monk/Desktop/STAR_decision_table_final/Model9c_Bmark2")
getwd()
#Remove all variables and turn off open graphics
rm(list=ls(all=TRUE))
graphics.off()
library(r4ss)
output.dir.plots = file.path(getwd(),'r4ss','plots')
#name and create output directory for r4ss
output.dir = file.path(getwd(),'r4ss')
dir.create(file.path(output.dir))
dir.create(file.path(output.dir.plots))
#run r4ss
with.covar = T		# T: without -nohess; F: with -nohess
#Remove all variables and turn off open graphics
rm(list=ls(all=TRUE))
graphics.off()
library(r4ss)
#name and create output directory for r4ss
output.dir = file.path(getwd(),'r4ss')
dir.create(file.path(output.dir))
output.dir.plots = file.path(getwd(),'r4ss','plots')
dir.create(file.path(output.dir.plots))
#run r4ss
with.covar = T		# T: without -nohess; F: with -nohess
mod1 = SS_output(dir = getwd(), forecast=F, covar=with.covar, ncol=1000)
# plots
SS_plots(mod1,
png = TRUE,
html = TRUE,
datplot = TRUE,
uncertainty = TRUE,
maxrows = 6,
maxcols = 6,
maxrows2 = 4,
maxcols2 = 4,
printfolder = '',
dir = output.dir.plots)
