Population Assessment

Packages Required

Data Formatting

walleye<-read.csv("E:/Grad School/Spring 2020/Quantitative Fisheries Analysis/Population Assessment/walleye.csv")
sampling<-read.csv("E:/Grad School/Spring 2020/Quantitative Fisheries Analysis/Population Assessment/merritsamples.csv")
left_join(sampling,walleye,by=c("Linker","WaterbodyCode","Area","MethodCode","Station","Effort"))->combined
combined %<>% 
  tidyr::separate(surveydate.x, c("month","day","year"),remove=F) %>% 
  mutate(year=as.numeric(year),
    lencat10=lencat(FishLength,w=10),
         psd=lencat(FishLength,breaks=psdVal("Walleye",units = "mm"),use.names=T))

Expected 3 pieces. Additional pieces discarded in 1928 rows [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, ...].

Takeaways: 1. All sampling events were successfull in capturing walleye, therefore there is no need to add CPUE of 0 to any sampling event if all sizes of walleye are being used. I will refer to this in teh Relative abundance section.

General Summary Statistics

#Otoliths collected summary
combined %>% filter(is.na(LengthGroup),!is.na(Age)) %>% group_by(year) %>% dplyr::summarise(n=n()) %>% summarise(mean=mean(n),sd=sd(n),n=sum(n)) #n=1409, mean=108 sd=24.65

#number of nets set = 51
length(unique(combined$Linker))

[1] 51

#number of nets set per year
combined %>% group_by(year) %>% dplyr::summarise(fish=n(),length(unique(Linker))) 

General Length Frequency

#Length Freq for all years to determine length at recruit and visualize data
combined %>% dplyr::select(Linker,year,LengthGroup,FishCount,FishLength,FishWeight,lencat10)->lengthdata
#Count fish by length category 
lengthdata %>%filter(is.na(LengthGroup)) %>% group_by( lencat10) %>% 
  dplyr::summarise(n=n())->no_LengthGroup1 
#Count fish by length category in Length Group column
lengthdata %>% group_by(LengthGroup,FishCount,year,Linker) %>% 
  dplyr::summarise() %>%  filter(!is.na(LengthGroup)) ->LengthGroup1 
#Combing for full count (use for plotting & analysis)
full_join(no_LengthGroup1,LengthGroup1,by=c("lencat10"="LengthGroup","n"="FishCount"))->fish_counts1
#Total fish Captured
fish_counts1 %>% summarise(n=sum(n)) #n = 1949

#Giving each fish unique identity (LengthGroup)
uncount(LengthGroup1, FishCount)->uncount_lengthgroup # n=32
full_join(subset(lengthdata,is.na(LengthGroup)),uncount_lengthgroup,by=c("FishLength"="LengthGroup","year","Linker"))->all_fish_lengths#n=1949
#Mean # of fish captured per year
all_fish_lengths %>% group_by(year) %>% dplyr::summarise(n=n()) %>% summarise(mean(n),sd(n),sum(n))

#Mean length (no Counted fish, minimal effect) n=1949, mean=443.7019, sd=97.06711
all_fish_lengths %<>% dplyr::select(-LengthGroup,-FishCount,-lencat10) %>% mutate(lencat10=lencat(FishLength,w=10))
all_fish_lengths %>% group_by() %>% dplyr::summarise(n=n(),mean=mean(FishLength),sd(FishLength))#1,949

#Length Freq plot
fish_counts1 %>% dplyr::select(-year,-Linker) %>% group_by(lencat10,n) %>% dplyr::summarise() %>% 
ggplot(aes(x=lencat10,y=n))+theme_classic()+
  geom_vline(aes(xintercept=635),linetype="dashed")+
  geom_vline(aes(xintercept=457),linetype="dashed")+
    geom_rect(aes(xmin=457,xmax=635,ymin=0,ymax=110),fill="grey50",alpha=.01)+
    geom_rect(aes(xmin=457,xmax=750,ymin=0,ymax=110),fill="grey50",alpha=.01)+
  geom_col(fill="grey60",col="black",position = "dodge")+
scale_fill_manual(values = alpha("white", .8))  +
geom_vline(aes(xintercept=444),cex=1.2)+
  labs(x="Total Length (mm)",y="Count")+
  scale_x_continuous(breaks = seq(0,750,50))+
  theme(axis.line = element_blank(),
        axis.text = element_text(family = "serif",size = 12))->all_lengths_histo
ggsave(all_lengths_histo,file="E:/Grad School/Spring 2020/Quantitative Fisheries Analysis/Population Assessment/tiffs/all_lengths_histo.tiff",dpi=500,width = 9,height = 6)

Takeaways: 1. First, there were a number of fish and LINKERS that used counts instead of measuring each fish. This means that to plot and analyze length properly, I must combined the two types of data. 2. Total fish captured was 1,949. Using this number, I made sure the data that was analyzed and plot always had a n=1,949. 3. Mean total lengths of combined sampling years=444 sd=97.06711.

Stat Analysis of Length Here I analyzed mean lengths across years and plotted means with 95% CI’s

#test of normality
all_fish_lengths$year<-as.factor(all_fish_lengths$year)#making year a factor to allow stat testing
shapiro.test(all_fish_lengths$FishLength)

    Shapiro-Wilk normality test

data:  all_fish_lengths$FishLength
W = 0.99016, p-value = 3.14e-10

(ggdensity(all_fish_lengths$FishLength))#normality will be assumed
Summarize(all_fish_lengths$FishLength)

         n       mean         sd        min         Q1     median         Q3        max 
1949.00000  443.70190   97.06711  179.00000  380.00000  454.00000  508.00000  734.00000 

lm(FishLength~year,data=subset(all_fish_lengths,FishLength>=250 ))->lm_lengths_year
summary(lm_lengths_year)

Call:
lm(formula = FishLength ~ year, data = subset(all_fish_lengths, 
    FishLength >= 250))

Residuals:
     Min       1Q   Median       3Q      Max 
-205.124  -60.518   -4.191   49.550  280.541 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  439.459      6.682  65.766  < 2e-16 ***
year1995      11.798      9.593   1.230  0.21888    
year1996      41.666      9.085   4.586 4.82e-06 ***
year1997     -53.227      9.030  -5.894 4.45e-09 ***
year1998       3.868      8.500   0.455  0.64911    
year1999     -17.762     10.356  -1.715  0.08648 .  
year2000      25.624      9.626   2.662  0.00784 ** 
year2001     -12.035     11.275  -1.067  0.28592    
year2002       6.590     11.363   0.580  0.56203    
year2003       9.541      9.950   0.959  0.33770    
year2004      53.392     10.303   5.182 2.43e-07 ***
year2005      67.606      9.751   6.933 5.64e-12 ***
year2006      18.340      9.496   1.931  0.05359 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 83.73 on 1883 degrees of freedom
Multiple R-squared:  0.1211,    Adjusted R-squared:  0.1155 
F-statistic: 21.62 on 12 and 1883 DF,  p-value: < 2.2e-16

leng_tuk <- multcomp::glht(lm_lengths_year,mcp(year="Tukey"))#no difference in age
summary(leng_tuk)

Completion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > abseps

     Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts


Fit: lm(formula = FishLength ~ year, data = subset(all_fish_lengths, 
    FishLength >= 250))

Linear Hypotheses:
                 Estimate Std. Error t value Pr(>|t|)    
1995 - 1994 == 0   11.798      9.593   1.230   0.9911    
1996 - 1994 == 0   41.666      9.085   4.586    <0.01 ***
1997 - 1994 == 0  -53.227      9.030  -5.894    <0.01 ***
1998 - 1994 == 0    3.868      8.500   0.455   1.0000    
1999 - 1994 == 0  -17.762     10.356  -1.715   0.8859    
2000 - 1994 == 0   25.624      9.626   2.662   0.2698    
2001 - 1994 == 0  -12.035     11.275  -1.067   0.9976    
2002 - 1994 == 0    6.590     11.363   0.580   1.0000    
2003 - 1994 == 0    9.541      9.950   0.959   0.9991    
2004 - 1994 == 0   53.392     10.303   5.182    <0.01 ***
2005 - 1994 == 0   67.606      9.751   6.933    <0.01 ***
2006 - 1994 == 0   18.340      9.496   1.931   0.7720    
1996 - 1995 == 0   29.868      9.234   3.235   0.0630 .  
1997 - 1995 == 0  -65.025      9.180  -7.084    <0.01 ***
1998 - 1995 == 0   -7.930      8.658  -0.916   0.9995    
1999 - 1995 == 0  -29.560     10.486  -2.819   0.1895    
2000 - 1995 == 0   13.825      9.766   1.416   0.9711    
2001 - 1995 == 0  -23.833     11.395  -2.092   0.6646    
2002 - 1995 == 0   -5.209     11.482  -0.454   1.0000    
2003 - 1995 == 0   -2.257     10.085  -0.224   1.0000    
2004 - 1995 == 0   41.594     10.434   3.987    <0.01 ** 
2005 - 1995 == 0   55.808      9.889   5.643    <0.01 ***
2006 - 1995 == 0    6.542      9.638   0.679   1.0000    
1997 - 1996 == 0  -94.893      8.648 -10.973    <0.01 ***
1998 - 1996 == 0  -37.798      8.093  -4.671    <0.01 ***
1999 - 1996 == 0  -59.428     10.024  -5.928    <0.01 ***
2000 - 1996 == 0  -16.042      9.269  -1.731   0.8786    
2001 - 1996 == 0  -53.701     10.971  -4.895    <0.01 ***
2002 - 1996 == 0  -35.076     11.061  -3.171   0.0748 .  
2003 - 1996 == 0  -32.124      9.604  -3.345   0.0439 *  
2004 - 1996 == 0   11.727      9.969   1.176   0.9940    
2005 - 1996 == 0   25.940      9.398   2.760   0.2171    
2006 - 1996 == 0  -23.326      9.133  -2.554   0.3340    
1998 - 1997 == 0   57.095      8.031   7.109    <0.01 ***
1999 - 1997 == 0   35.465      9.974   3.556   0.0229 *  
2000 - 1997 == 0   78.851      9.215   8.557    <0.01 ***
2001 - 1997 == 0   41.192     10.926   3.770   0.0106 *  
2002 - 1997 == 0   59.817     11.016   5.430    <0.01 ***
2003 - 1997 == 0   62.768      9.552   6.571    <0.01 ***
2004 - 1997 == 0  106.619      9.919  10.749    <0.01 ***
2005 - 1997 == 0  120.833      9.345  12.930    <0.01 ***
2006 - 1997 == 0   71.567      9.078   7.883    <0.01 ***
1999 - 1998 == 0  -21.630      9.497  -2.278   0.5273    
2000 - 1998 == 0   21.755      8.696   2.502   0.3670    
2001 - 1998 == 0  -15.903     10.492  -1.516   0.9511    
2002 - 1998 == 0    2.721     10.586   0.257   1.0000    
2003 - 1998 == 0    5.673      9.052   0.627   1.0000    
2004 - 1998 == 0   49.524      9.439   5.247    <0.01 ***
2005 - 1998 == 0   63.738      8.834   7.215    <0.01 ***
2006 - 1998 == 0   14.472      8.551   1.692   0.8953    
2000 - 1999 == 0   43.386     10.517   4.125    <0.01 ** 
2001 - 1999 == 0    5.727     12.044   0.476   1.0000    
2002 - 1999 == 0   24.352     12.127   2.008   0.7227    
2003 - 1999 == 0   27.304     10.814   2.525   0.3510    
2004 - 1999 == 0   71.154     11.139   6.388    <0.01 ***
2005 - 1999 == 0   85.368     10.631   8.030    <0.01 ***
2006 - 1999 == 0   36.102     10.398   3.472   0.0293 *  
2001 - 2000 == 0  -37.659     11.423  -3.297   0.0521 .  
2002 - 2000 == 0  -19.034     11.510  -1.654   0.9098    
2003 - 2000 == 0  -16.082     10.117  -1.590   0.9308    
2004 - 2000 == 0   27.769     10.465   2.654   0.2735    
2005 - 2000 == 0   41.983      9.922   4.231    <0.01 ** 
2006 - 2000 == 0   -7.283      9.671  -0.753   0.9999    
2002 - 2001 == 0   18.625     12.920   1.441   0.9666    
2003 - 2001 == 0   21.576     11.697   1.845   0.8223    
2004 - 2001 == 0   65.427     11.999   5.453    <0.01 ***
2005 - 2001 == 0   79.641     11.529   6.908    <0.01 ***
2006 - 2001 == 0   30.375     11.314   2.685   0.2556    
2003 - 2002 == 0    2.952     11.782   0.251   1.0000    
2004 - 2002 == 0   46.803     12.081   3.874    <0.01 ** 
2005 - 2002 == 0   61.017     11.614   5.253    <0.01 ***
2006 - 2002 == 0   11.751     11.401   1.031   0.9983    
2004 - 2003 == 0   43.851     10.763   4.074    <0.01 ** 
2005 - 2003 == 0   58.065     10.236   5.673    <0.01 ***
2006 - 2003 == 0    8.799      9.993   0.880   0.9996    
2005 - 2004 == 0   14.214     10.580   1.344   0.9810    
2006 - 2004 == 0  -35.052     10.345  -3.388   0.0396 *  
2006 - 2005 == 0  -49.266      9.796  -5.029    <0.01 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)

(cld(leng_tuk))

Completion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > abseps

 1994  1995  1996  1997  1998  1999  2000  2001  2002  2003  2004  2005  2006 
 "bc"  "bd"  "df"   "a"  "bc"   "b" "cde"  "bc"  "bd"  "bc"  "ef"   "f"  "cd" 

results_len<-data.frame(year=as.factor(c(1994,  1995,  1996  ,1997 , 1998  ,1999  ,2000 , 2001 ,2002 , 2003,  2004 , 2005  ,2006 )),label=c(  "bc" , "bd",  "df"   ,"a" , "bc"  , "b" ,"cde" ,"bc" , "bd" , "bc" , "ef" ,  "f"  ,"cd" ))
all_fish_lengths  %>% filter(FishLength>=250) %>% 
group_by(year) %>% dplyr::summarise(mean=mean(FishLength),
                   difference=508-mean,#difference between mean and preferred length (managment goal)
                   n=n(),
                   sd=sd(FishLength),
                   error= qnorm(0.975)*sd/sqrt(n),
                   LCI = mean - error,
                   UCI = mean + error)->all_fish_mean_lengths
ggplot()+theme_classic()+
 geom_point(data=all_fish_mean_lengths,aes(x=year,y=mean),cex=2)+
  geom_errorbar(data=all_fish_mean_lengths,aes(x=year,y=mean,ymin=LCI,ymax=UCI))+
  labs(x="Year",y="Mean Total Length (mm)")+
  geom_hline(aes(yintercept=457),linetype="dotted")+
  geom_smooth(data=subset(all_fish_lengths,FishLength>=250),
              aes(x=as.numeric(year),y=lencat10),method="lm",cex=1,linetype="dashed",col="grey50",se=F)+
  geom_text(data=results_len,aes(x=as.numeric(year), y=all_fish_mean_lengths$UCI, label=label),fontface="bold",vjust=-.5)+
  theme(text=element_text(family = "serif",size=12),
        axis.text.x = element_text(angle = 45,hjust=.8,vjust=.8))->len_means
ggsave(len_means,file="E:/Grad School/Spring 2020/Quantitative Fisheries Analysis/Population Assessment/tiffs/len_means.tiff",dpi=500,width=8,height=6)

summary(lm(FishLength~as.numeric(year),data=subset(all_fish_lengths,FishLength>=250 )))#Multiple R-squared:  0.02115,   Adjusted R-squared:  0.02064 F-statistic: 40.93 on 1 and 1894 DF,  p-value: 1.983e-10

Call:
lm(formula = FishLength ~ as.numeric(year), data = subset(all_fish_lengths, 
    FishLength >= 250))

Residuals:
     Min       1Q   Median       3Q      Max 
-219.757  -61.795    6.823   57.613  288.297 

Coefficients:
                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)      428.2817     4.0339 106.171  < 2e-16 ***
as.numeric(year)   3.4212     0.5348   6.398 1.98e-10 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 88.1 on 1894 degrees of freedom
Multiple R-squared:  0.02115,   Adjusted R-squared:  0.02064 
F-statistic: 40.93 on 1 and 1894 DF,  p-value: 1.983e-10

#creating PSD categories
all_fish_lengths %<>% mutate(psd=lencat(FishLength,breaks=psdVal("Walleye",units = "mm"),use.names=T))

Takeaways: 1. Shapiro-wilks indicated abnormal, however visual appears to fit assumption of Normality well. 2. Assuming Normality, I preceded with an one-way Analysis of Variance (ANOVA) to see if there was differences in the mean total length between years. 3. Anova indicated significant difference between years(DF 12-1936, Fvalue 17.55 ,<2e-16 ***) 4. Next, I conducted a Tukey’s pos hoc test to find which years were significantly different from one another. There were 35 combinations of years that were different, signifying high variability in the mean length between years.

Discussion: High variability in mean length led to many significant differences between years. One possible reason for high variability is low sample size. ALl years had four nets sets, except 1999 which had three. The mean number of fish caught per years was 149.9231 sd=47.87042 n=1949. Although the mean number per year was high, there are many reasons why 4 nets is not adequate to sample Walleye. This warrents discussion in the report.

Relative Abundance Calculations Before starting: All sampling events in the sampling data file (n=51) were successful in catching at least one walleye. Therefore subsequent analysis which look at a particular size category must have at least 51 sampling events. If NA’s are present, replace with a 0.

#Each subset should be 51 nets set, regardless of success
combined %>% group_by(year,Linker) %>% dplyr::summarise() ->netset_success#allnets set
combined %>% group_by(year,Linker) %>% dplyr::summarise(fish=n())%>% mutate(model="All Fish")->fishcaptured
left_join(fishcaptured,netset_success,by=c("year","Linker"))->effort_all #n = 51? YE
effort_all %<>% mutate(fish=ifelse(is.na(fish),0,fish),model=ifelse(is.na(model),"All Fish",model)) 
effort_all %>% 
  group_by(year) %>% dplyr::summarise(n=n(),
                                      fish=sum(fish),
                                      CPUE=fish/n,model="All Fish")->CPUE_all
combined %>% filter(FishLength>=250) %>%  group_by(year,Linker) %>% dplyr::summarise( fish=n()) %>% mutate(model="Stock")->fishcaptured_stock
left_join(fishcaptured_stock,netset_success,by=c("year","Linker"))->effort_stock#n = 51? YES
effort_stock %<>% mutate(fish=ifelse(is.na(fish),0,fish),model=ifelse(is.na(model),"Stock",model)) 
effort_stock %>% 
  group_by(year) %>% dplyr::summarise(n=n(),
                                      fish=sum(fish),
                                      CPUE=fish/n,model="Stock")->CPUE_stock
combined %>% filter(FishLength>=510) %>%  group_by(year,Linker) %>% dplyr::summarise( fish=n())%>% mutate(model="Preferred")->fishcaptured_pref
left_join(fishcaptured_pref,netset_success,by=c("year","Linker"))->effort_pref#n = 51? YES
effort_pref %<>% mutate(fish=ifelse(is.na(fish),0,fish),model=ifelse(is.na(model),"Preferred",model))
effort_pref %>% 
  group_by(year) %>% dplyr::summarise(n=n(),
                                      fish=sum(fish),
                                      CPUE=fish/n,model="Preferred")->CPUE_pref
combined %>% filter(FishLength<=380) %>%  group_by(year,Linker) %>% dplyr::summarise( fish=n())%>% mutate(model="Quality")->fishcaptured_qual
left_join(netset_success,fishcaptured_qual,by=c("year","Linker"))->effort_qual#n = 51? YES
effort_qual %<>% mutate(fish=ifelse(is.na(fish),0,fish),model=ifelse(is.na(model),"Quality",model))
effort_qual %>% 
  group_by(year) %>% dplyr::summarise(n=n(),
                                      fish=sum(fish),
                                      CPUE=fish/n,model="Quality")->CPUE_qual
rbind(CPUE_stock,CPUE_pref,CPUE_qual)->effort
(ggplot(data=effort)+theme_classic()+
  geom_point(aes(x=year,y=CPUE,pch=model))+
  scale_shape_manual(values = c(15,1,2,0))+
  geom_path(aes(x=year,y=CPUE,linetype=model,group=model))+
    labs(x="Year",y="CPUE (fish per net night)")+
  scale_linetype_manual(values=c("solid","dashed","dotted"))+
  scale_x_continuous(breaks = seq(1994,2006,2))->CPUE_plot)

 # ggsave(CPUE_plot,file="E:/Grad School/Spring 2020/Quantitative Fisheries Analysis/Population Assessment/tiffs/CPUE.tiff",dpi=500,width = 9,height = 6)

Takeaways:

Discussion: Abundance of Preferred-length individuals is variable during the early years of the survey, a trough occured around 10 years in to the restriction and then as been very steadily on the rise, with little variation. Something that must be looked into is the number of years which is required for individuals to reach preferred-length in this system. Additionally, poor reservoir conditions accompanied by poor recruit years may have influenced the variation seen in the early years of the regulation.

Relative Abundance stats

#visual inspection of CPUE for different size classes
ggdensity(effort_all$fish)

#negative bionomial CPUE  
#Log10 transformation
effort_all$log_fish<-(log10(effort_all$fish)+1)
effort_qual$log_fish<-(log10(effort_qual$fish)+1)
effort_stock$log_fish<-(log10(effort_stock$fish)+1)
effort_pref$log_fish<-(log10(effort_pref$fish)+1)
#Visual inspection of log CPUE for different size classes
ggdensity(effort_all$log_fish)

shapiro.test(effort_all$log_fish)

    Shapiro-Wilk normality test

data:  effort_all$log_fish
W = 0.97823, p-value = 0.4666

ggdensity(effort_qual$log_fish)

ggdensity(effort_stock$log_fish)

ggdensity(effort_pref$log_fish)

#normality not achieve for all
summary(lm(year~fish,data=effort_pref))

Call:
lm(formula = year ~ fish, data = effort_pref)

Residuals:
    Min      1Q  Median      3Q     Max 
-6.3882 -3.8882  0.3461  3.6118  5.1713 

Coefficients:
             Estimate Std. Error  t value Pr(>|t|)    
(Intercept) 1.999e+03  1.008e+00 1982.187   <2e-16 ***
fish        1.469e-01  9.068e-02    1.619    0.112    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 3.753 on 49 degrees of freedom
Multiple R-squared:  0.0508,    Adjusted R-squared:  0.03143 
F-statistic: 2.623 on 1 and 49 DF,  p-value: 0.1118

summary(lm(year~fish,data=effort_qual))

Call:
lm(formula = year ~ fish, data = effort_qual)

Residuals:
   Min     1Q Median     3Q    Max 
-6.254 -1.539 -0.254  2.392  8.166 

Coefficients:
              Estimate Std. Error  t value Pr(>|t|)    
(Intercept) 2001.39270    0.76164 2627.739   <2e-16 ***
fish          -0.14233    0.05867   -2.426    0.019 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 3.639 on 49 degrees of freedom
Multiple R-squared:  0.1072,    Adjusted R-squared:  0.089 
F-statistic: 5.885 on 1 and 49 DF,  p-value: 0.019

summary(lm(year~fish,data=effort_stock))

Call:
lm(formula = year ~ fish, data = effort_stock)

Residuals:
    Min      1Q  Median      3Q     Max 
-7.3208 -2.1646 -0.1996  2.7748  6.7814 

Coefficients:
              Estimate Std. Error  t value Pr(>|t|)    
(Intercept) 2002.58217    1.29171 1550.338   <2e-16 ***
fish          -0.07008    0.03239   -2.163   0.0354 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 3.68 on 49 degrees of freedom
Multiple R-squared:  0.08719,   Adjusted R-squared:  0.06856 
F-statistic: 4.681 on 1 and 49 DF,  p-value: 0.03541

rbind(effort_stock,effort_pref,effort_qual)->effort2
duration<-data.frame(year=unique(effort2$year),duration=c(2:14))
left_join(effort2,duration,by=c("year"))->effort2
(effort2 %>% mutate(model=factor(model,levels = c("Stock","Quality","Preferred"))) %>% 
ggplot()+theme_classic()+
  geom_point(aes(x=duration,y=log_fish),pch=1)+
  geom_smooth(aes(x=duration,y=log_fish,group=model),linetype="dashed",col="black",method = "lm",se=F)+
    labs(x="Years of Regulation ",y="log10 CPUE (fish / net)")+facet_wrap(~model)+
  theme(text = element_text(family = "serif",size = 12),
    strip.background = element_blank())+
  scale_x_continuous(breaks=seq(4,12,4))->abundance_duration_categories)
  ggsave(abundance_duration_categories,file="E:/Grad School/Spring 2020/Quantitative Fisheries Analysis/Population Assessment/tiffs/abundance_duration_categories.tiff",dpi=500,width = 9,height = 6)

  
effort %>% ggplot(aes(x=year,y=CPUE))+theme_classic()+geom_point(aes(pch=model))+geom_line(aes(group=model))

Condition

ggdensity(combined$FishLength)

ggdensity(combined$FishWeight)

wsVal(("Walleye") ,simplify=TRUE,unit="metric")->wsWA
combined %>% filter(is.na(LengthGroup)) %>% 
  mutate( logL=log10(FishLength),
          logW=log10(FishWeight),
    Ws=(10^(wsWA[["int"]]+wsWA[["slope"]]*logL)),
         Wr=FishWeight/Ws*100)->condition
ggdensity(condition$Wr)

shapiro.test(condition$Wr)

    Shapiro-Wilk normality test

data:  condition$Wr
W = 0.99609, p-value = 7.379e-05

condition$year<-as.factor(condition$year)
condition_lm<-lm(Wr~year,data=condition)
summary(lm(Wr~as.numeric(year),data=subset(condition,psd=="preferred")))

Call:
lm(formula = Wr ~ as.numeric(year), data = subset(condition, 
    psd == "preferred"))

Residuals:
    Min      1Q  Median      3Q     Max 
-43.575  -5.638  -0.637   5.897  28.174 

Coefficients:
                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)       97.7692     0.9828  99.477  < 2e-16 ***
as.numeric(year)   0.4897     0.1146   4.275 2.36e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 9.75 on 429 degrees of freedom
Multiple R-squared:  0.04085,   Adjusted R-squared:  0.03862 
F-statistic: 18.27 on 1 and 429 DF,  p-value: 2.361e-05

summary(condition_lm)

Call:
lm(formula = Wr ~ year, data = condition)

Residuals:
    Min      1Q  Median      3Q     Max 
-39.705  -5.302  -0.416   5.461  35.844 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  99.4015     0.6577 151.136  < 2e-16 ***
year1995      6.5635     0.9549   6.873 8.48e-12 ***
year1996     -0.8240     0.8924  -0.923   0.3559    
year1997     -7.0235     0.8945  -7.852 6.79e-15 ***
year1998    -11.4655     0.8424 -13.610  < 2e-16 ***
year1999     -2.0954     1.1438  -1.832   0.0671 .  
year2000      1.1543     0.9584   1.204   0.2286    
year2001     -2.7048     1.1130  -2.430   0.0152 *  
year2002      1.0816     1.1346   0.953   0.3405    
year2003      1.1279     0.9912   1.138   0.2553    
year2004      5.5812     1.0271   5.434 6.21e-08 ***
year2005      4.0983     0.9637   4.253 2.21e-05 ***
year2006      0.1176     0.9154   0.128   0.8978    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 8.423 on 1902 degrees of freedom
  (2 observations deleted due to missingness)
Multiple R-squared:  0.2876,    Adjusted R-squared:  0.2831 
F-statistic: 63.99 on 12 and 1902 DF,  p-value: < 2.2e-16

cond_tuk <- multcomp::glht(condition_lm,mcp(year="Tukey"))#no difference in age
(cld(cond_tuk))

Completion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > absepsCompletion with error > abseps

1994 1995 1996 1997 1998 1999 2000 2001 2002 2003 2004 2005 2006 
"cd"  "f" "cd"  "b"  "a" "cd" "de"  "c" "ce" "ce"  "f" "ef" "cd" 

condition %>% ggplot(aes(x=as.numeric(year),y=Wr))+theme_classic()+geom_point()+
  geom_smooth(method="lm",se=F)+facet_wrap(~psd)
conditiontable<-condition %>% filter(!is.na(Wr)) %>% group_by(year,psd) %>% dplyr::summarise(mean=mean(Wr),
                   n=n(),
                   sd=sd(Wr),
                   error= qnorm(0.975)*sd/sqrt(n),
                   LCI = mean - error,
                   UCI = mean + error)
conditiontable %<>% filter(n>1)
#looking at psd
conditiontable %>% filter(psd%in%c("preferred","quality","stock")) %>% 
ggplot(aes(x=year,y=mean))+theme_classic()+geom_hline(aes(yintercept=100),linetype="dashed",col="grey50")+
  geom_point()+geom_errorbar(aes(ymin=LCI,ymax=UCI),width=.5)+facet_wrap(~psd,ncol = 1)+
  geom_path(aes(group=psd))+labs(x="Year",y="")+
  theme(strip.background = element_blank(),
        axis.text.x = element_text(angle = 45,vjust=.5),
        axis.ticks.x = element_line())->condition_psd
#COmbined PSD categories
conditiontable_all<-condition %>% filter(!is.na(Wr),psd%in%c("preferred","quality","stock")) %>% group_by(year) %>% dplyr::summarise(mean=mean(Wr),
                   n=n(),
                   sd=sd(Wr),
                   error= qnorm(0.975)*sd/sqrt(n),
                   LCI = mean - error,
                   UCI = mean + error)
condition %>% filter(!is.na(Wr),psd%in%c("preferred","quality","stock")) %>% group_by() %>% dplyr::summarise(mean=mean(Wr))

conditiontable_all %>% 
ggplot(aes(x=year,y=mean))+theme_classic()+geom_hline(aes(yintercept=100),linetype="dashed",col="grey50")+
  geom_point()+geom_errorbar(aes(ymin=LCI,ymax=UCI),width=.5)+
  geom_path(aes())+labs(y="Mean Relative Weight (Wr)",x="Year")+
  theme(strip.background = element_blank(),
        axis.text.x = element_text(angle = 45,vjust=.5))->condition_all
ggarrange(condition_all,condition_psd,ncol=2)->condition_plot

geom_path: Each group consists of only one observation. Do you need to adjust the
group aesthetic?

 ggsave(condition_plot,file="E:/Grad School/Spring 2020/Quantitative Fisheries Analysis/Population Assessment/tiffs/condition_plot.tiff",dpi=500,width = 9,height = 6)

Age and Growth Von Bertalanffy Growth

  combined %>%filter(!is.na(Age)) %>%  dplyr::select(year,FishLength,Age:Annulus30)->age_data
#visualize data
ggplot(data=age_data,aes(y=FishLength,x=Age))+theme_classic()+geom_point()+geom_smooth(method = "loess",se=F,col="black")

#Von B Growth curve
(svTyp1 <- vbStarts(FishLength~Age,data=age_data))#start ages

$`Linf`
[1] 772.2778

$K
[1] 0.1845135

$t0
[1] -1.634973

(vbTyp <- function(age,Linf,K,t0) Linf*(1-exp(-K*(age-t0))))#making VBLF function

function(age,Linf,K,t0) Linf*(1-exp(-K*(age-t0)))

# Model Formula
(fitTyp <- nls(FishLength~vbTyp(Age,Linf,K,t0), data = age_data, start = svTyp1))

Nonlinear regression model
  model: FishLength ~ vbTyp(Age, Linf, K, t0)
   data: age_data
    Linf        K       t0 
731.1103   0.1932  -1.8463 
 residual sum-of-squares: 2064737

Number of iterations to convergence: 4 
Achieved convergence tolerance: 5.821e-07

# extract parameters from nls() object with coef()
coef(fitTyp)#linf=731.1103119   k=0.1932223  t0=-1.8462618 

       Linf           K          t0 
731.1103119   0.1932223  -1.8462618 

# fit confidence intervals with from nls() with confint()
confint(fitTyp)

Waiting for profiling to be done...

            2.5%       97.5%
Linf 703.1967490 764.8466416
K      0.1714312   0.2156038
t0    -2.0682640  -1.6471211

bootTyp <-nlsBoot(fitTyp)
headtail(bootTyp$coefboot,n=2)

           Linf         K        t0
[1,]   724.0829 0.1976146 -1.842576
[2,]   737.1911 0.1882552 -1.903502
[998,] 727.8511 0.2001721 -1.750622
[999,] 735.6324 0.1942933 -1.798989

confint(bootTyp, plot=TRUE)

         95% LCI     95% UCI
Linf 703.2287734 764.1821776
K      0.1713741   0.2157366
t0    -2.0760166  -1.6578644

summary(bootTyp)

------
Bootstrap statistics
        Estimate  Std. error
Linf 732.3180293 15.58353070
K      0.1929125  0.01110442
t0    -1.8533778  0.10537776

------
Median of bootstrap estimates and percentile confidence intervals
          Median        2.5%       97.5%
Linf 731.5494841 703.2287734 764.1821776
K      0.1931985   0.1713741   0.2157366
t0    -1.8462313  -2.0760166  -1.6578644

CIall<-confint(bootTyp)
nd <- data.frame(Age=c(0:13))
nd$Length<-predict(fitTyp,nd)
residPlot(fitTyp)
vonb  <- data.frame(vbTyp(0:13,Linf = 731.1103119, K = 0.1932223, t0 = -1.8462618 ))
colnames(vonb)[1] <- "vonb"
vonb$max<-vbTyp(0:13,Linf = 764.8466416, K = 0.2156038, t0 = -1.6471211 )
vonb$Age<-c(0:13)
vonb$min<-vbTyp(0:13,Linf = 703.1967490, K = 0.1714312, t0 =-2.0682640 )
ggplot(data=age_data,aes(y=FishLength,x=Age))+theme_classic()+geom_point(pch=1)+
  scale_y_continuous(breaks = seq(150,750,100))+
  geom_line(data=vonb,aes(x=Age,y=vonb),cex=1)+
  scale_x_continuous(breaks = seq(1,13,1))+
  labs(x="Age",y="Total Length (mm)")+
  geom_hline(aes(yintercept=457),linetype="dashed",col="grey30")+
  geom_hline(aes(yintercept=635),linetype="dashed",col="grey30")->growth_curve
ggsave(growth_curve,file="E:/Grad School/Spring 2020/Quantitative Fisheries Analysis/Population Assessment/tiffs/growth_curve.tiff",dpi=500,width = 9,height = 6)

Age-at-Length

#Length Categories
age_data %<>%  mutate(Annulus0=ifelse(is.na(Annulus1),Edge,NA)) %>% mutate(lencat10=lencat(FishLength,10)) 
#Frequency of Length cat by age, creating Matrix for proportion of age categories within lencat categories. 
alk.freq <- xtabs(~lencat10+Age,data=age_data) 
rowSums(alk.freq)

170 180 190 200 210 220 240 250 260 270 280 290 300 310 320 330 340 350 360 370 380 390 
  1   2   8   3   1   1   2   5   8  13  22  33  25  23  21  17  30  32  32  40  41  40 
400 410 420 430 440 450 460 470 480 490 500 510 520 530 540 550 560 570 580 590 600 610 
 42  59  52  53  59  56  55  66  65  52  47  49  59  54  42  27  29  29  16  15  14   7 
620 630 640 650 660 670 680 690 700 710 720 730 
 10  10   8   5   4   5   4   2   6   1   6   1 

alk <- prop.table(alk.freq,margin=1)
cc.mlr <- nnet::multinom(Age~lencat10,data=age_data)

# weights:  39 (24 variable)
initial  value 3614.013645 
iter  10 value 3038.674952
iter  20 value 2077.005888
iter  30 value 1649.602111
iter  40 value 1438.217383
iter  50 value 1409.010800
iter  60 value 1408.866317
iter  70 value 1407.389981
iter  80 value 1406.734571
iter  90 value 1406.190467
iter 100 value 1406.098352
final  value 1406.098352 
stopped after 100 iterations

lens <- seq(170,740,10)
alk.sm <- predict(cc.mlr,data.frame(lencat10=lens),type="probs")
row.names(alk.sm) <- lens
round(alk.sm,3)

        0     1     2     3     4     5     6     7     8     9    10    11    13
170 0.999 0.001 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
180 0.997 0.003 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
190 0.991 0.009 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
200 0.972 0.028 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
210 0.916 0.084 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
220 0.773 0.227 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
230 0.516 0.484 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
240 0.250 0.749 0.001 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
250 0.094 0.902 0.003 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
260 0.031 0.962 0.007 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
270 0.010 0.977 0.013 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
280 0.003 0.972 0.025 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
290 0.001 0.951 0.047 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
300 0.000 0.912 0.087 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
310 0.000 0.844 0.155 0.001 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
320 0.000 0.737 0.260 0.003 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
330 0.000 0.591 0.402 0.007 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
340 0.000 0.427 0.557 0.016 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
350 0.000 0.276 0.691 0.033 0.001 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
360 0.000 0.161 0.776 0.061 0.001 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
370 0.000 0.087 0.804 0.106 0.003 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
380 0.000 0.044 0.779 0.170 0.008 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
390 0.000 0.021 0.707 0.256 0.016 0.001 0.000 0.000 0.000 0.000 0.000 0.000 0.000
400 0.000 0.009 0.599 0.359 0.031 0.002 0.000 0.000 0.000 0.000 0.000 0.000 0.000
410 0.000 0.004 0.469 0.467 0.056 0.004 0.000 0.000 0.000 0.000 0.000 0.000 0.000
420 0.000 0.001 0.338 0.558 0.094 0.008 0.000 0.000 0.000 0.000 0.000 0.000 0.000
430 0.000 0.000 0.224 0.614 0.144 0.017 0.000 0.000 0.000 0.000 0.000 0.000 0.000
440 0.000 0.000 0.138 0.625 0.205 0.031 0.001 0.000 0.000 0.000 0.000 0.000 0.000
450 0.000 0.000 0.079 0.593 0.271 0.055 0.002 0.000 0.000 0.000 0.000 0.000 0.000
460 0.000 0.000 0.042 0.527 0.336 0.090 0.005 0.000 0.000 0.000 0.000 0.000 0.000
470 0.000 0.000 0.021 0.439 0.390 0.139 0.009 0.001 0.000 0.000 0.000 0.000 0.000
480 0.000 0.000 0.010 0.343 0.426 0.201 0.017 0.002 0.001 0.000 0.000 0.000 0.000
490 0.000 0.000 0.004 0.251 0.435 0.273 0.030 0.003 0.002 0.001 0.000 0.000 0.000
500 0.000 0.000 0.002 0.173 0.418 0.349 0.048 0.006 0.003 0.001 0.000 0.000 0.000
510 0.000 0.000 0.001 0.112 0.378 0.418 0.073 0.011 0.005 0.002 0.000 0.000 0.000
520 0.000 0.000 0.000 0.068 0.323 0.474 0.106 0.017 0.008 0.003 0.000 0.000 0.000
530 0.000 0.000 0.000 0.040 0.262 0.511 0.144 0.026 0.013 0.005 0.000 0.000 0.000
540 0.000 0.000 0.000 0.022 0.202 0.524 0.187 0.037 0.019 0.008 0.001 0.000 0.000
550 0.000 0.000 0.000 0.012 0.150 0.515 0.234 0.052 0.026 0.011 0.001 0.000 0.000
560 0.000 0.000 0.000 0.006 0.107 0.486 0.280 0.069 0.035 0.016 0.002 0.000 0.000
570 0.000 0.000 0.000 0.003 0.073 0.443 0.323 0.088 0.046 0.022 0.003 0.000 0.000
580 0.000 0.000 0.000 0.001 0.048 0.390 0.360 0.110 0.058 0.029 0.005 0.000 0.000
590 0.000 0.000 0.000 0.001 0.031 0.333 0.390 0.132 0.070 0.037 0.007 0.000 0.000
600 0.000 0.000 0.000 0.000 0.019 0.277 0.411 0.155 0.083 0.046 0.009 0.000 0.000
610 0.000 0.000 0.000 0.000 0.012 0.224 0.422 0.177 0.096 0.056 0.013 0.000 0.000
620 0.000 0.000 0.000 0.000 0.007 0.178 0.424 0.198 0.109 0.066 0.017 0.000 0.000
630 0.000 0.000 0.000 0.000 0.004 0.138 0.419 0.217 0.121 0.078 0.023 0.000 0.000
640 0.000 0.000 0.000 0.000 0.002 0.106 0.406 0.234 0.132 0.089 0.029 0.001 0.000
650 0.000 0.000 0.000 0.000 0.001 0.080 0.388 0.249 0.143 0.101 0.037 0.001 0.000
660 0.000 0.000 0.000 0.000 0.001 0.059 0.366 0.261 0.151 0.113 0.047 0.002 0.000
670 0.000 0.000 0.000 0.000 0.000 0.044 0.340 0.270 0.158 0.125 0.058 0.004 0.000
680 0.000 0.000 0.000 0.000 0.000 0.032 0.312 0.276 0.164 0.136 0.072 0.008 0.001
690 0.000 0.000 0.000 0.000 0.000 0.022 0.282 0.277 0.167 0.146 0.086 0.015 0.004
700 0.000 0.000 0.000 0.000 0.000 0.016 0.250 0.273 0.166 0.153 0.102 0.028 0.013
710 0.000 0.000 0.000 0.000 0.000 0.011 0.213 0.258 0.159 0.154 0.116 0.051 0.038
720 0.000 0.000 0.000 0.000 0.000 0.007 0.169 0.228 0.142 0.145 0.123 0.086 0.101
730 0.000 0.000 0.000 0.000 0.000 0.004 0.117 0.175 0.111 0.119 0.113 0.125 0.236
740 0.000 0.000 0.000 0.000 0.000 0.002 0.065 0.109 0.069 0.079 0.084 0.147 0.446

len.n <- xtabs(~lencat10,data=age_data)
tmp <- sweep(alk,MARGIN=1,FUN="*",STATS=len.n)
ad1 <- colSums(tmp)
round(prop.table(ad1),3)

    0     1     2     3     4     5     6     7     8     9    10    11    13 
0.012 0.128 0.196 0.223 0.164 0.160 0.069 0.023 0.013 0.008 0.003 0.001 0.001 

alkAgeDist(alk,lenA.n=rowSums(alk.freq),len.n=len.n)

alkMeanVar(alk,FishLength~lencat10+Age,data=age_data,len.n=len.n)->age_key
#Applying to all-unaged individuals. Final step will be to combined aged fish (done), unaged measured fish(A), and unaged unmeasured fish (B, LengthGroup)
#A. Make Dataframe of all unaged fish lengths (no.ages.final)
combined %>%filter(is.na(LengthGroup),is.na(Age)) %>% group_by( year,lencat10) %>% 
  dplyr::summarise(n=n())->no.ages.measured#508 measured, non-aged fish
combined %>%filter(!is.na(LengthGroup),is.na(Age)) %>%  group_by(LengthGroup,FishCount,year) %>% 
  dplyr::summarise()->no.ages.LG
full_join(no.ages.measured,no.ages.LG,by=c("year"="year","lencat10"="LengthGroup","n"="FishCount"))->no.ages.final
sum(no.ages.final$n)#508+32=540 correct

[1] 540

uncount(no.ages.final,n)->no.ages.final
no.ages.final$Age<-as.numeric(NA)
assign.ages <- alkIndivAge(alk,Age~lencat10,data=no.ages.final)
#Df of aged fish
combined %>% filter(!is.na(Age)) %>% dplyr::select(year,lencat10,Age)-> aged.fish
#combined ages all fish (n=1,949)
rbind(aged.fish,assign.ages)->all.age

Description: Here I created an age-at-length key using a proportional age-length matrix. Essentially, this matrix samples the population of aged fish I have and states the proportion of fish with a certain length category that fit into a given age category. Using this matrix I was able to create an age_length key that fit my VonB very very closely. The next step was to assign ages to fish that were not aged using the age-length matrix and procedure described by Iserman & knight (2005?). THe last step was to combine the dataframes with originally aged fish and the df with predicted ages. THhs will allow comparisons of mortality

Back Calculations

#formatting data for back-calculations
age_data %>% dplyr::select(- Annulus14:-Annulus30) %>% #delete unneeded columns 
  mutate(Annulus0=ifelse(is.na(Annulus1),Edge,NA)) %>% #presevering age 0 in new dataframe
  pivot_longer(cols = Annulus1:Annulus0,names_to = "Annuli",values_to = "radius") %>% #Format long data
  separate(Annuli,c("delete","Annuli"),"s") %>% #cleanup 
  dplyr::select(-delete) %>% #cleanup 
  filter(!is.na(radius)) %>% #cleanup
  mutate(Annuli=as.numeric(Annuli),#creating formula for Dahl Lea calculations
         dahl_lea= round(FishLength*(radius/Edge),2),
         year_at=(year-Age+Annuli),
         nat_year=year-Age)->backcalc
lm1 <- lm(dahl_lea~radius,data=backcalc)
anova(lm1)

Analysis of Variance Table

Response: dahl_lea
            Df   Sum Sq  Mean Sq F value    Pr(>F)    
radius       1 52041432 52041432  4868.2 < 2.2e-16 ***
Residuals 4859 51943571    10690                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

a <- coef(lm1)[1]               # get the correction factor (first coefficient in lm1)
backcalc$Fraser_lee <- a + (backcalc$FishLength-a)*(backcalc$radius/backcalc$Edge)
left_join(backcalc,duration,by=c("year"))->backcalc
backcalc %<>% mutate(era=ifelse(duration<=5,"0-5",ifelse(duration>10,"10+","5-10")))
#visuals
ggplot(data=backcalc,aes(x=year_at,y=Fraser_lee))+geom_point()+facet_wrap(~nat_year)

ggplot(data=subset(FL_sumbackcalc,Annuli!=14),aes(x=Annuli,y=mean))+theme_classic()+
  geom_hline(aes(yintercept=508))+
  geom_vline(aes(xintercept=4.3),linetype="dashed")+
  geom_hline(aes(yintercept=508))+
  geom_line(data=vonb,aes(x=Age,y=vonb),linetype="dashed",col="grey30",cex=1)+
  geom_point(col="black")+scale_x_continuous(breaks = seq(1,13,1))+
  labs(y="Mean Back-Calculated Length (mm)",x="Age")+
  geom_errorbar(aes(ymin=LCI,ymax=UCI),width=.5,col="black")

Error in subset(FL_sumbackcalc, Annuli != 14) : 
  object 'FL_sumbackcalc' not found

Takeaways: 1. No evidence for lee’s Phenonomenon occuring

Traditional, Catch-Curve Mortality Estimates Across all Years Excerise material only, do not include in report

Takeaway: 1. Catch curve analysis indicated little difference in mortality across the 3 distinct eras. This is shown in the models lm13,lm12,lm23 results and mort_plot_era. 2. Catch_curve analysis however, is likely not fitting as a key assumption is that recruitement stays consistant across all other covariates like time or location. We know that walleye can have highly variable recruitment caused by either natural variations or stocking efforts. Therefore, catch curves are likely not the most useful too to investigate the mortality of this population.

Year-Class Strength

all.age %>% mutate(yearcl=year-Age)  %>% filter(Age%in%c(1,2,3,4,5,6),yearcl%in%c(1992:2001)) %>% 
  dplyr::select(-lencat10,-duration) %>%  group_by(year,yearcl,Age) %>% dplyr::summarise(n=n(),logn=log(n))->year_class

Error in eval(lhs, parent, parent) : object 'all.age' not found