7-Maturity_Calibration.Rmd 13.3 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
---
title: "Maturity exploration"
author: "Camille POULET"
date: "08/03/2021"
output: word_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
rm(list = ls())
```

```{r library, include=FALSE}

library(dplyr)
library(tidyr)
library(ggplot2)
library(readr)
library(forcats)
library(knitr)
library(officedown)
library(flextable)
library(stringr)
```

```{r data, include=FALSE}

#load("growthCurves.RData")
load('SOS.rdata')
load("growth_optimal.RData")
source("../GR3D_Rdescription/GR3Dfunction.R")

```

# 1. Estimation of American shad maturity schedules (ASMFC, 2020)

Maturity schedules or ogives were calculated by determining the proportion of mature and immature fish as a function of age for both male and female based on the capture of fish of a species.
When age data were not available, length was used as a surrogate. 
In either case, this methods required satisfying the assumption of equal mixing of mature and immature fish and equal recruitment to gear.For most anadromous fish, this condition was not satisfied, and the assumption was then made for 'in-river' captures when all encountered fish were reproductively mature (ASMFC, 2020)

To estimate the age at maturity and identify previous spawning activity, two methods have been used: scales and otholiths.
However, Maki et al. (2001) developped a new approah to estimate maturity probability using scale-derived age in conjonction with the scale-derived patterns of spawning (Maki et al. 2001)
Consistent with this method, the ASMFC used at sea mortality estimates to inform a maturity ogive, by assuming that all fish observed returning to river system were reproductively mature.

## 1.1. Maturity ogives

Otholiths data were used as age informations and scales as a binary variable discriminating first time spawners from multiple spawners.
Following this, a matrix of age-at-capture and age_at first spawn was build for both male and female from each of the 3 regional metapopulations (semelparous, sourthern iteroparous, northern iteroparous)

Assuming a proportional capture of fish across years, and that all fish are susceptible to capture, the proportion of fish mature at age a (p_a) would then calculated as follow:

$$ 
p_a = \frac{observed_a}{(observed_a + immature_a)} 
$$
where p_a is the proportion of fish mature at age a, observed_a are mature fish observed at age a, with a assumed to be between years 1 and 9, and immature a are the immature fish observed each age.
As it was considered that returning fish were exclusively mature, there is no observations of immature fish at age. 
Then, the proportion of immature fish was estimated based on all the returns greater than age a and the estimated mortality (M) at sea as below: 

$$ 
immature_a = \sum_{j=a+1}^{9}\frac{(observed_j)}{\exp^{-M(j-a)}} 
$$
## 1.2.Natural mortality


The natural mortality at sea (M) was calculated based on maximum ages following the method described in Then et al. (2015) for each region. Estimates of natural mortality by region and sex were provided in Table 10, and Table 27 from the ASMFC. 

```{r echo = FALSE, fig.cap= "Table 1. Natural mortality estimates by region and sex. Maximum age data for the northern iteroparous (NI) region were borrowed from the sourthern iteroparous (SI) region due to more prevalent and long term impacts of man-made barriers in the NI region", out.width='100%'}

knitr::include_graphics("American_Shad_natural_mortality.PNG")

```


## 1.3 Results 

By considering both sex and the mortality correction approach, female have a later maturity-at-age than males in each region, as commonly observed in other shad groups(Taverny, 1991).
Both sexes begin to mature at 3 years,with males maturing on average in 4.3 years and females maturing on average in 4.6 years (ASMFC, 2020).

```{r echo = FALSE, fig.cap= "Table 2. American shad maturity-at-age by region and sex using the mortality correcion approach", out.width='100%'}

knitr::include_graphics("American_Shad_maturity_at_age.PNG")

```
#2. Application for the GR3D model 

## 2.1. Age-at-maturity

We computed the number of mature (pA) and immature individuals at-age for both sex in the 3 regions following the methodology describe above (see 1.).
We then computed the relative number of mature fish-at-age regarding the number of mature fish for all age-class combined.
The mean age-at maturity was then computed as a weighted mean following the relative number of mature fish-at-age (Table 3)

```{r explore the relationship between age and lentgh at maturity from growth curves, echo = FALSE, include = TRUE }

#build a dataframe to store the age at maturity compute from the ogives 
growthCurvesBySex <- expand_grid(dataVerif, sex= c("M","F")) %>% 
  group_by(metapop, sex) %>% 
  mutate(age_at_maturity = case_when (
    metapop == "semelparous" & sex =="F" ~ 5.097,
    metapop == "semelparous" & sex =="M" ~ 4.247,
    metapop == "southern iteroparous" & sex =="F" ~ 5.591,
    metapop == "southern iteroparous" & sex =="M" ~ 5.065,
    metapop == "northern iteroparous" & sex =="F" ~ 5.153,
    metapop == "northern iteroparous" & sex =="M" ~ 4.611)) 

```


```{r echo = FALSE, fig.cap= "Table 3. American shad mean age at maturity by region and sex using the maturity schedule from the ASMFC report (2020)", out.width='100%' }

growthCurvesBySex %>% 
  distinct(metapop,sex,age_at_maturity) %>% 
  flextable() %>% 
  autofit()

```

```{r interpolation on age based on the mean length from Stich et al. (2020), echo = FALSE , include = TRUE, fig.cap= "Table 4. American shad mean age at maturity by region and sex based on the growth cruve from Stich et al. (2020)", out.width='100% }

#mean_size by metapop 
growthCurvesBySex %>% 
  group_by(metapop) %>% 
  mutate(mean_LStich = mean(LStich)) %>% 
  distinct(metapop,mean_LStich)

#Interpoplation on age 
growthCurvesBySex %>% 
  group_by(metapop) %>% 
  mutate(mean_LStich = mean(LStich)) %>% 
  transmute(age_at_maturity = approx(LStich, age, xout = mean_LStich , method = "linear", ties="ordered")$y) %>%  
  distinct(metapop,age_at_maturity) %>% 
  flextable() %>% 
  autofit()

```
Age-at-maturity derived from Stich growth curve are lowest than those estimates from the maturity ogive if we considered a mean size of approximately 38cm.

## 2.2. Length at maturity

In the GR3D model, females are mature for the first time at the length of 45cm (TL), and 40cm (TL) for males. 

The mean length-at-age was computed for each river in the 3 region based on the VBF growth functions derived from Stich et al. (2020) described in "6-Growth_Calibration.Rmd". 

```{r growth curves , echo = FALSE, include = TRUE }

dataVerif %>% 
  group_by(metapop, age) %>% 
  summarise(L = mean(L), LStich = mean(LStich)) %>% 
  ungroup() %>% 
  mutate(metapop = factor(metapop, levels = c('northern iteroparous',  'southern iteroparous', 'semelparous'))) %>% 
  ggplot(aes(x = age, color = metapop)) + 
  geom_line(aes(y = LStich)) +
  geom_line(aes(y = L)) 


```

```{r plot growth curve, echo = FALSE , include = TRUE }

#plot
growthCurvesBySex %>% 
  group_by(metapop, age,sex) %>% 
  summarise(L = mean(L), 
            LStich = mean(LStich), 
            age_at_maturity = mean(age_at_maturity)) %>% 
  ungroup() %>% 
  mutate(metapop = factor(metapop, levels = c('northern iteroparous',  'southern iteroparous', 'semelparous'))) %>% 
  ggplot(aes(x = age, color = metapop)) + 
  geom_line(aes(y = LStich)) +
  geom_line(aes(y = L)) +
  geom_vline(aes(xintercept = age_at_maturity, color = sex)) +
  labs(x = 'Age (year)', y = 'Fork length (cm)', color = "Metapopulation (Stich)") +
  facet_wrap(.~ metapop, ncol =2)

```


```{r interpolation on size based on the mean age-at-maturity, echo = FALSE , include = TRUE, fig.cap= "Table 4. American shad mean age at maturity by region and sex based on the maturity ogives from the ASFCM report (2020)", out.width='100% }

#Interpolation on size 
growthCurvesBySex %>% 
  group_by(metapop) %>% 
  mutate(meanLStich = mean(LStich)) %>% 
  ungroup() %>% 
  group_by(metapop,sex) %>% 
  mutate(size_at_maturity = approx(age, LStich, xout = age_at_maturity , method = "linear", ties="ordered")$y) %>%  
  distinct(metapop,sex,size_at_maturity) %>% 
  flextable() %>% 
  autofit()


```


```{r explore the mean length-at-maturity from growth curves derived from Stich et al (2020), echo = FALSE, include = TRUE }

#mean length at age for each metapopulation
df<-dataVerif %>% 
  filter(season == "spring") %>% 
  group_by(age, metapop) %>% 
  summarize("mean_length_GR3D" = mean(L), 
            "mean_length_Stich" = mean(LStich)) 
  

df<-split(df,df$metapop) 

#TODO: optimize 
#semelparous region
df$semelparous%>% 
  arrange(age) %>% 
  select("age","mean_length_GR3D","mean_length_Stich") %>% 
  flextable() %>% 
  autofit()

#southern iteroaprous region
df$`southern iteroparous`%>% 
  arrange(age) %>% 
  select("age","mean_length_GR3D","mean_length_Stich") %>% 
  flextable() %>% 
  autofit()

#northern iteroparous region
df$`northern iteroparous`%>% 
  arrange(age) %>% 
  select("age","mean_length_GR3D","mean_length_Stich") %>% 
  flextable() %>% 
  autofit()

#mean lentgh at age - coastwide
 dataVerif %>% 
  filter(season == "spring") %>% 
  group_by(age) %>% 
  summarize("mean_length_GR3D" = mean(L), 
            "mean_length_Stich" = mean(LStich)) %>% 
  flextable() %>% 
  autofit()

```


```{r specific river analysis}

growthInBasin %>% 
  filter(basin_name =="St. John" |basin_name =="Miramichi") %>% 
  group_by(basin_name, age) %>% 
  summarise(LStich = mean(LStich)) %>% 
  ungroup() %>% 
  ggplot(aes(x = age, color = basin_name)) + 
  geom_line(aes(y = LStich))+
  geom_hline(yintercept = 40.5, color = "lightblue", linetype = 1)+ #St John
  geom_vline(xintercept = 3.76, linetype = 1, color = "lightblue")+ #St John
  geom_hline(yintercept = 44.35, color = "pink", linetype = 1)+#St John
  geom_vline(xintercept = 4.25, linetype = 1, color = "pink")+#St John
 
  geom_hline(yintercept = 42.75, color = "blue", linetype = 1)+ #Miramichi
  geom_vline(xintercept = 4.23, linetype = 1, color = "blue")+ #Miramichi
  geom_hline(yintercept = 46.87, color = "red", linetype = 1)+#Miramichi
  geom_vline(xintercept = 4.68, linetype = 1, color = "red")#Miramichi
  

```

```{r compute ogives schedule with R function - equivalent to the xls file }
#Compute ages from ogives 

# Ogive_df_semelparous <- data.frame(age = seq (1, 9, 1),
#                                    pa = c(0,0,0.01,0.09,0.33,0.63,0.92,1,1,
#                                           0,0,0.07,0.29,0.60,0.86,1,1,1),
#                                    gender = c(rep ("female",9),rep("male",9)),
#                                    metapop = "semelparous") 
# 
# 
# Ogive_df_SI <- data.frame(age = seq (1, 9, 1),
#                           pa = c(0,0,0,0.04,0.27,0.64,0.81,0.9,1,
#                                  0,0,0,0.12,0.50,0.79,0.84,0.80,1),
#                           gender = c(rep ("female",9),rep("male",9)),
#                           metapop = "southern iteroparous") 
# 
# 
# Ogive_df_NI <- data.frame(age = seq (1, 9, 1),
#                           pa = c(0,0,0,0.04,0.69,0.69,0.9,1,1,
#                                  0,0,0.03,0.21,0.82,1,1,1,1),
#                           gender = c(rep ("female",9),rep("male",9)),
#                           metapop = "northern iteroparous") 
# 
# 
# Ogives_df <- rbind(Ogive_df_semelparous, Ogive_df_SI, Ogive_df_NI)
# 
# # Ogives_df<-expand_grid(metapop = c("semelparous", "southern iteroparous", "northern iteroparous"),gender = c("female","male"), age = seq(1,9,1)) 
# #   
# 
# S= 1000
# 
# Ogives <- Ogives_df %>% 
#   ungroup() %>% 
#   #rename("proportion_of_mature" = "pa") %>% 
#   mutate(M = ifelse(
#     metapop =="semelparous", 0.65, 0.47))
# 
# 
# 
# miaou = function (Zsea, prop_mature){
#   
#   Ogives <- Ogives %>% 
#     #arrange(metapop, gender, age) %>% 
#     mutate(stockAtSea = if_else(age == 1, S, 0)) 
#   
#   ages <- Ogives %>% ungroup() %>% 
#     select(age) %>% distinct() %>% arrange(age) %>% unlist(use.names = FALSE)
#   
#   for (i in 2:length(ages)) {
#     
#     previousAge  = ages[i - 1]
#     currentAge = ages[i]
#     
#     previousStockAtSea = Ogives  %>% 
#     distinct() %>% 
#     filter(age == previousAge) %>%
#     select(stockAtSea) %>% unlist(use.names = FALSE) 
#     # 
#     
#     # pa= Ogives %>% 
#     #   filter(age == currentAge) %>% 
#     #   select(pa) %>% unlist(use.names = FALSE)
#     # 
#     # M= Ogives %>% 
#     #   filter(age== currentAge) %>% 
#     #   select(M) %>% unlist(use.names = FALSE)
#     # 
#     # Ogives = Ogives %>% 
#     #   mutate(stock = previousStockAtSea*exp(-Zsea),
#     #          mature = stock*pa)
#     
# 
#     stock = stockAtSea[previousAge]*exp(-Zsea)
#     mature = stock*prop_mature
#     
#     currentStockAtSea = (previousStockAtSea -mature)*exp(-Zsea)
#     
#     Ogives= Ogives%>% 
#       mutate(stockAtSea = replace(stockAtSea, age == ages[i], CurrentStockAtSea))
#     
#   }
#   
#   return (mature)
# }
# 
# miaou(0.65,0.01)
# 
# 
# #function
# 
# Ogives= Ogives %>% 
#   group_by(age,gender,metapop) %>% 
#   mutate(matureFish = surviving_fish*(pa)) %>% 
#   mutate(eff_surviving_fish = surviving_fish - matureFish) %>% 
#   mutate(eff_matureFish = eff_surviving_fish *(pa*100))
# 
# 
# growthCurves %>% 
#   filter(season == "spring") %>% 
#   group_by(age, metapop) %>% 
#   summarize("mean_GR3D" = mean(L), 
#             "mean_Stich" = mean(LStich)) %>% 
#   arrange(age,mean_GR3D) %>% 
#   filter(mean_GR3D >=40) %>% #males
#   filter(mean_GR3D >=45) #females

```