diff --git a/exploration/GR3D_Rdescription/GR3Dfunction.R b/exploration/GR3D_Rdescription/GR3Dfunction.R index 2903bd6e788b16218a91e8829d79a79be63880db..1327205183329906c69c35ba218f47ce3749100e 100644 --- a/exploration/GR3D_Rdescription/GR3Dfunction.R +++ b/exploration/GR3D_Rdescription/GR3Dfunction.R @@ -39,6 +39,7 @@ vonBertalanffyInverse = function(L, L0, Linf, K){ # von Bertalanffy increment # pas cohérent avec la temperature effet sur le coeff de croissance mais ca marche + vonBertalanffyIncrement = function(nStep, L0, Linf, K, timeStepDuration, sigma, withTempEffect=FALSE, TrefAtSea = c(9.876946, 13.489854, 15.891487, 11.554104) ){ tempEffect = temperatureEffect( TrefAtSea , 3, 17, 26) L = matrix(nrow = nStep + 1) @@ -156,8 +157,24 @@ vonBertalanffyWithRandomVector = function(L, Linf, K, timeStepDuration, randomVe return(L) } +computeOgive = function(lengthTrajectories, lengthAtMaturity){ + ogive <- lengthTrajectories %>% + group_by(age) %>% + summarise(nTotal = n()) %>% + left_join(lengthTrajectories %>% group_by(age) %>% + filter(L >= lengthAtMaturity) %>% + summarise(taller = n(), .groups = 'drop'), + by = 'age') %>% + replace_na(list(taller = 0)) %>% + mutate(mature = c(0,diff(taller)), + immature = nTotal - cumsum(mature), + p = if_else(mature + immature > 0 , mature / (mature + immature), 1)) %>% + select(age, p) + return(ogive) +} -computeOgive = function(lengthTrajectories, lengthAtMaturity) { +# deprecated +computeOgive3 = function(lengthTrajectories, lengthAtMaturity) { ogive <- lengthTrajectories %>% group_by(age) %>% diff --git a/exploration/NEA_calibration_offline/maturationPL.Rmd b/exploration/NEA_calibration_offline/maturationPL.Rmd index 9253bb80e5190b43439594f46957eece2f98e84e..9987933c33d03b7aa9ba139c37718985e518a353 100644 --- a/exploration/NEA_calibration_offline/maturationPL.Rmd +++ b/exploration/NEA_calibration_offline/maturationPL.Rmd @@ -67,130 +67,131 @@ temperaturePattern <- growthInBasin %>% ``` ```{r local function} - -# generate a cohort of length trajectories -computeMultipleLengthTrajectories = function(temperaturePattern, Nind = 10, growPar){ - ages = temperaturePattern$age - - res <- expand_grid(ind = seq.int(Nind), age = ages) %>% - inner_join(temperaturePattern, by = 'age') %>% - arrange(age, ind) %>% - mutate(temperatureEffect = temperatureEffect(temperature, - growPar$tempMinGrow, - growPar$tempOptGrow, - growPar$tempMaxGrow), - L = if_else(age == 0, growPar$lengthAtHatching, 0)) - - for (i in 2:length(ages)) { - previousAge = ages[i - 1] - currentAge = ages[i] - tempEffect = res %>% filter(age == currentAge) %>% select(temperatureEffect) %>% unlist(use.names = FALSE) - previousL <- res %>% filter(age == previousAge) %>% select(L) %>% unlist(use.names = FALSE) - - currentL <- vonBertalanffyWithNextIncrement(L = previousL, - Linf = growPar$Linf, - K = growPar$kOpt, - timeStepDuration = currentAge - previousAge, - sigma = growPar$sigmaDeltaLVonBert, - tempEffect = tempEffect) - res = res %>% mutate(L =replace(L, age == ages[i], currentL) ) - } - return(res) -} - -computeMultipleLengthTrajectoriesWithRandomSeed = function(temperaturePattern, - Nind = 10, - growPar, - RNGseed =1){ - set.seed(RNGseed) - ages = temperaturePattern %>% - distinct(age) %>% - unlist(use.names = FALSE) - - res <- expand_grid(ind = seq.int(Nind), age = ages) %>% - mutate(random = rnorm(Nind * length(ages))) %>% - inner_join(temperaturePattern, by = 'age') %>% - arrange(age, ind) %>% - mutate(temperatureEffect = temperatureEffect(temperature, - growPar$tempMinGrow, - growPar$tempOptGrow, - growPar$tempMaxGrow), - L = if_else(age == 0, growPar$lengthAtHatching, 0)) - - for (i in 2:length(ages)) { - previousAge = ages[i - 1] - currentAge = ages[i] - tempEffect <- res %>% filter(age == currentAge) %>% - select(temperatureEffect) %>% unlist(use.names = FALSE) - - previousL <- res %>% filter(age == previousAge) %>% - select(L) %>% unlist(use.names = FALSE) - rnd <- res %>% filter(age == currentAge) %>% select(random) %>% unlist(use.names = FALSE) - currentL <- vonBertalanffyWithRandomVector(L = previousL, - Linf = growPar$Linf, - K = growPar$kOpt, - timeStepDuration = currentAge - previousAge, - randomVector = rnd, - sigma = growPar$sigmaDeltaLVonBert, - tempEffect = tempEffect) - res = res %>% mutate(L = replace(L, age == ages[i], currentL) ) - } - return(res) -} - -#computeMultipleLengthTrajectoriesWithRandomSeed(temperaturePattern, Nind = 10, -#growPar = growParUnisex, -# RNGseed = 1) - -vonBertalanffyWithRandomVector = function(L, Linf, K, timeStepDuration, randomVector, sigma, tempEffect ){ - if (sigma == 0) { - mu = if_else(L < Linf, log((Linf - L) * (1 - exp(-K * tempEffect * timeStepDuration))),-Inf) - #mu = log((Linf - L) * (1 - exp(-K * tempEffect * timeStepDuration))) - increment = exp(mu) - } - else { - mu = if_else(L < Linf, log((Linf - L) * (1 - exp(-K * tempEffect * timeStepDuration))) - (sigma * sigma) / 2, -Inf) - # mu = log((Linf - L) * (1 - exp(-K * tempEffect * timeStepDuration))) - (sigma * sigma) / 2 - increment = exp(randomVector * sigma + mu) - } - L = pmin(Linf, L + increment) - return(L) -} - - -computeOgive = function(lengthTrajectories, lengthAtMaturity) { - ogive <- - lengthTrajectories %>% - group_by(age) %>% - summarise(nTotal=n()) %>% - left_join(lengthTrajectories %>% group_by(age) %>% - filter(L > lengthAtMaturity) %>% - summarise(mature = n(), .groups = 'drop'), - by = 'age') %>% - replace_na(list(mature = 0)) %>% - mutate(p = mature/nTotal) %>% - select(age, p) - - return(ogive) -} - -#Deprecated -computeOgive2 = function(lengthTrajectories, lengthAtMaturity) { - ogive <- lengthTrajectories %>% - distinct(age) %>% arrange(age) %>% - left_join(lengthTrajectories %>% filter(L > lengthAtMaturity) %>% - group_by(ind) %>% - summarise(age = min(age), .groups = 'drop') %>% - group_by(age) %>% - summarise(mature = n(), .groups = 'drop'), - by = 'age') %>% - replace_na(list(mature = 0)) %>% - mutate(immature = sum(mature) - cumsum(mature), - p = if_else(mature + immature == 0, 0, mature/(mature + immature))) %>% - select(age, p) - - return(ogive) -} +# TRENSFERED IN GR3D functiond + +# # generate a cohort of length trajectories +# computeMultipleLengthTrajectories = function(temperaturePattern, Nind = 10, growPar){ +# ages = temperaturePattern$age +# +# res <- expand_grid(ind = seq.int(Nind), age = ages) %>% +# inner_join(temperaturePattern, by = 'age') %>% +# arrange(age, ind) %>% +# mutate(temperatureEffect = temperatureEffect(temperature, +# growPar$tempMinGrow, +# growPar$tempOptGrow, +# growPar$tempMaxGrow), +# L = if_else(age == 0, growPar$lengthAtHatching, 0)) +# +# for (i in 2:length(ages)) { +# previousAge = ages[i - 1] +# currentAge = ages[i] +# tempEffect = res %>% filter(age == currentAge) %>% select(temperatureEffect) %>% unlist(use.names = FALSE) +# previousL <- res %>% filter(age == previousAge) %>% select(L) %>% unlist(use.names = FALSE) +# +# currentL <- vonBertalanffyWithNextIncrement(L = previousL, +# Linf = growPar$Linf, +# K = growPar$kOpt, +# timeStepDuration = currentAge - previousAge, +# sigma = growPar$sigmaDeltaLVonBert, +# tempEffect = tempEffect) +# res = res %>% mutate(L =replace(L, age == ages[i], currentL) ) +# } +# return(res) +# } +# +# computeMultipleLengthTrajectoriesWithRandomSeed = function(temperaturePattern, +# Nind = 10, +# growPar, +# RNGseed =1){ +# set.seed(RNGseed) +# ages = temperaturePattern %>% +# distinct(age) %>% +# unlist(use.names = FALSE) +# +# res <- expand_grid(ind = seq.int(Nind), age = ages) %>% +# mutate(random = rnorm(Nind * length(ages))) %>% +# inner_join(temperaturePattern, by = 'age') %>% +# arrange(age, ind) %>% +# mutate(temperatureEffect = temperatureEffect(temperature, +# growPar$tempMinGrow, +# growPar$tempOptGrow, +# growPar$tempMaxGrow), +# L = if_else(age == 0, growPar$lengthAtHatching, 0)) +# +# for (i in 2:length(ages)) { +# previousAge = ages[i - 1] +# currentAge = ages[i] +# tempEffect <- res %>% filter(age == currentAge) %>% +# select(temperatureEffect) %>% unlist(use.names = FALSE) +# +# previousL <- res %>% filter(age == previousAge) %>% +# select(L) %>% unlist(use.names = FALSE) +# rnd <- res %>% filter(age == currentAge) %>% select(random) %>% unlist(use.names = FALSE) +# currentL <- vonBertalanffyWithRandomVector(L = previousL, +# Linf = growPar$Linf, +# K = growPar$kOpt, +# timeStepDuration = currentAge - previousAge, +# randomVector = rnd, +# sigma = growPar$sigmaDeltaLVonBert, +# tempEffect = tempEffect) +# res = res %>% mutate(L = replace(L, age == ages[i], currentL) ) +# } +# return(res) +# } +# +# #computeMultipleLengthTrajectoriesWithRandomSeed(temperaturePattern, Nind = 10, +# #growPar = growParUnisex, +# # RNGseed = 1) +# +# vonBertalanffyWithRandomVector = function(L, Linf, K, timeStepDuration, randomVector, sigma, tempEffect ){ +# if (sigma == 0) { +# mu = if_else(L < Linf, log((Linf - L) * (1 - exp(-K * tempEffect * timeStepDuration))),-Inf) +# #mu = log((Linf - L) * (1 - exp(-K * tempEffect * timeStepDuration))) +# increment = exp(mu) +# } +# else { +# mu = if_else(L < Linf, log((Linf - L) * (1 - exp(-K * tempEffect * timeStepDuration))) - (sigma * sigma) / 2, -Inf) +# # mu = log((Linf - L) * (1 - exp(-K * tempEffect * timeStepDuration))) - (sigma * sigma) / 2 +# increment = exp(randomVector * sigma + mu) +# } +# L = pmin(Linf, L + increment) +# return(L) +# } +# +# +# computeOgive = function(lengthTrajectories, lengthAtMaturity) { +# ogive <- +# lengthTrajectories %>% +# group_by(age) %>% +# summarise(nTotal=n()) %>% +# left_join(lengthTrajectories %>% group_by(age) %>% +# filter(L > lengthAtMaturity) %>% +# summarise(mature = n(), .groups = 'drop'), +# by = 'age') %>% +# replace_na(list(mature = 0)) %>% +# mutate(p = mature/nTotal) %>% +# select(age, p) +# +# return(ogive) +# } +# +# #Deprecated +# computeOgive2 = function(lengthTrajectories, lengthAtMaturity) { +# ogive <- lengthTrajectories %>% +# distinct(age) %>% arrange(age) %>% +# left_join(lengthTrajectories %>% filter(L > lengthAtMaturity) %>% +# group_by(ind) %>% +# summarise(age = min(age), .groups = 'drop') %>% +# group_by(age) %>% +# summarise(mature = n(), .groups = 'drop'), +# by = 'age') %>% +# replace_na(list(mature = 0)) %>% +# mutate(immature = sum(mature) - cumsum(mature), +# p = if_else(mature + immature == 0, 0, mature/(mature + immature))) %>% +# select(age, p) +# +# return(ogive) +# } ``` @@ -632,52 +633,64 @@ SSEforOneMetapop(fullPar = enframe(c(par_0, fixedPar)) %>% pivot_wider(), RNGseed = 1) toc() -fixedPar0 = c(tempMinGrow = 1.6, #min temperature experienced by fish in the 11 RB from Stich et al, 2020 - tempMaxGrow = 27.9,#max temperature experienced by fish in the 11 RB from Stich et al, 2020 - lengthAtHatching = 2.8) -# avirer -parFemale = enframe(c(resA_par, fixedPar0)) %>% pivot_wider() %>% select(tempMinGrow, tempOptGrow, tempMaxGrow, - lengthAtHatching, - Linf = linfVonBert, - kOpt = kOptForFemale, - sigmaDeltaLVonBert, - lengthAtMaturity = lFirstMaturityForFemale) - - -#================================ ICI ======================================= -rm(lengthTrajectories) -lengthTrajectories = computeMultipleLengthTrajectoriesWithRandomSeed(temperaturePattern = temperaturePattern %>% filter(metapop == "semelparous"), - Nind = 1000, - growPar = parFemale, - RNGseed = 1) %>% - arrange(ind,age) %>% - filter(season == "spring") - - -lengthTrajectories %>% - ggplot(aes(x = age, y= L, color = as.factor(ind))) + -geom_path( show.legend = FALSE) - - -lengthTrajectories %>% - group_by(age) %>% - summarise(L = mean(L)) %>% - print(n=Inf) - -lengthTrajectories %>% - computeOgive(lengthAtMaturity = parFemale$lengthAtMaturity) %>% - print(n=Inf) - - - - -vonBertalanffyWithRandomVector(L = c(2.8, 2.8), - Linf = parFemale$Linf, - timeStepDuration = 0.25, - K = parFemale$kOpt, - randomVector = c(0.184, -0.253), - sigma = parFemale$sigmaDeltaLVonBert, - tempEffect = c(0.0639, 0.0639)) +# +# fixedPar +# res_Abis$par +# # fixedPar0 = c(tempMinGrow = 1.6, #min temperature experienced by fish in the 11 RB from Stich et al, 2020 +# # tempMaxGrow = 27.9,#max temperature experienced by fish in the 11 RB from Stich et al, 2020 +# # lengthAtHatching = 2.8) +# # avirer +# parFemale = enframe(c(res_Abis$par, fixedPar)) %>% pivot_wider() %>% select(tempMinGrow, tempOptGrow, tempMaxGrow, +# lengthAtHatching, +# Linf = linfVonBert, +# kOpt = kOptForFemale, +# sigmaDeltaLVonBert, +# lengthAtMaturity = lFirstMaturityForFemale) +# +# +# #================================ ICI ======================================= +# rm(lengthTrajectories) +# lengthTrajectories = computeMultipleLengthTrajectoriesWithRandomSeed(temperaturePattern = temperaturePattern %>% filter(metapop == "semelparous"), +# Nind = 1000, +# growPar = parFemale, +# RNGseed = 1) %>% +# arrange(ind,age) %>% +# filter(season == "spring") +# +# +# lengthTrajectories %>% +# ggplot(aes(x = age, y= L, color = as.factor(ind))) + +# geom_path( show.legend = FALSE) +# +# +# lengthTrajectories %>% +# group_by(age) %>% +# summarise(L = mean(L)) %>% +# print(n=Inf) +# +# +# lengthTrajectories <- computeMultipleLengthTrajectoriesWithRandomSeed(temperaturePattern = temperaturePattern %>% filter(metapop == "semelparous"), +# Nind = 100000, +# growPar = parFemale %>% mutate(sigmaDeltaLVonBert =.6), +# RNGseed = 1) %>% +# arrange(ind,age) %>% +# filter(season == "spring") +# +# +# lengthTrajectories %>% +# computeOgive(lengthAtMaturity = parFemale$lengthAtMaturity) %>% +# print(n=Inf) +# +# +# +# +# vonBertalanffyWithRandomVector(L = c(2.8, 2.8), +# Linf = parFemale$Linf, +# timeStepDuration = 0.25, +# K = parFemale$kOpt, +# randomVector = c(0.184, -0.253), +# sigma = parFemale$sigmaDeltaLVonBert, +# tempEffect = c(0.0639, 0.0639)) ``` ```{r set parameters values for calibration, warning = "FALSE"} @@ -721,14 +734,15 @@ par_0 = c(tempOptGrow = 5, #from optimisation on growth curves kOptForFemale = 0.32, #XML - Alosa alosa kOptForMale = .21, #XML - Alosa alosa lFirstMaturityForMale = 40, #XML - Alosa sapidissima - lFirstMaturityForFemale = 45)#XML - Alosa sapidissima + lFirstMaturityForFemale = 45,#XML - Alosa sapidissima + linfVonBert = 76, #XML - Alosa alosa + sigmaDeltaLVonBert = .2) #XML - Alosa alosa fixedPar = c(tempMinGrow = 1.6, #min temperature experienced by fish in the 11 RB from Stich et al, 2020 tempMaxGrow = 27.9,#max temperature experienced by fish in the 11 RB from Stich et al, 2020 - lengthAtHatching = 2.8,# XML - Alosa sapidissima - linfVonBert = 76, #XML - Alosa alosa - sigmaDeltaLVonBert = .2) #XML - Alosa alosa + lengthAtHatching = 2.8)# XML - Alosa sapidissima + #optimisation # XML - Alosa sapidissia)