Commit 842dedce authored by patrick.lambert's avatar patrick.lambert
Browse files

with automatic dermination of process number

parent f2c8629a
......@@ -18,21 +18,23 @@ rm(list = ls())
# load the function to calculate the metapopulations
source("exploration/GR3D_Rdescription/GR3Dmetapopulation.R")
# load the xml
# ---------------- load the xml
fishXML <- xmlToList(xmlParse("data/input/northeastamerica/fishRIOBasin_Sapidissima_Rjava.xml"))
no_processesEachStep = tibble( process = names(fishXML[["species.DiadromousFishGroup"]][["processes"]][["processesEachStep"]][])) %>%
mutate(no_R = row_number()) %>%
filter(process != "comment") %>%
mutate(no_java = row_number()-1) %>%
no_processesEachStep = tibble( process = names(fishXML[["species.DiadromousFishGroup"]][["processes"]][["processesEachStep"]][])) %>%
mutate(no_R = row_number()) %>%
filter(process != "comment") %>%
mutate(no_java = row_number()-1) %>%
select(no_R, no_java, process)
no_processesAtEnd = tibble( process = names(fishXML[["species.DiadromousFishGroup"]][["processes"]][["processesAtEnd"]][])) %>%
mutate(no_R = row_number()) %>%
filter(process != "comment") %>%
mutate(no_java = row_number()-1) %>%
no_processesAtEnd = tibble( process = names(fishXML[["species.DiadromousFishGroup"]][["processes"]][["processesAtEnd"]][])) %>%
mutate(no_R = row_number()) %>%
filter(process != "comment") %>%
mutate(no_java = row_number()-1) %>%
select(no_R, no_java, process)
# prepare Sensitivity analysis
source("exploration/NEA_sensitivity_analysis/prepareSensitivityAnalysis.R")
# ====================================== LOCAL FUNCTIONS =====================================
......@@ -59,7 +61,7 @@ convertFromSALtoR = function(java_out){
# .jarray(parametersNamesENV), .jarray(thetasENV))
# }
exploreEasyRun = function(simulationName){
exploreEasyRun = function(simulationName, no_processesEachStep, no_processesAtEnd){
# call GR3D method to get model outputs ----------------------------------------------
# 00. simulation name
......@@ -68,21 +70,35 @@ exploreEasyRun = function(simulationName){
# "simulationName")
# 0. abundance (caution spawner abundance of the last reproduction (no averaging) )
noProcess <- no_processesEachStep %>%
filter(process == "analysis.AnalyseFishDistribution") %>%
select(no_java) %>%
unlist(use.names = FALSE)
lastSpawnerAdundance = .jcall("analysis.EasyRun","D", "getSingleValueFromAquanismGroupProcess",
"processes.processesEachStep.10.getSpawnerEffective")
paste0("processes.processesEachStep.", noProcess, ".getSpawnerEffective"))
# 1. northernmost populated basin_id
# 2. centroid latitude
noProcess <- no_processesEachStep %>%
filter(process == "analysis.AnalyseFishDistribution") %>%
select(no_java) %>%
unlist(use.names = FALSE)
out_range = .jcall("analysis.EasyRun","[[Ljava/lang/String;", "getValuesFromAquanismGroupProcess",
"processes.processesEachStep.10.exportToR")
paste0("processes.processesEachStep.", noProcess,".exportToR"))
range = suppressMessages(convertFromSALtoR(out_range))
# 3. log-likelihood ----
# the proportion of reproductions leading to a recruitment > 50 over the 1920-1950 period was used
# as a proxy of predicted probability of presence in a basin.
# Then a log-likelihood of a binomial distribution is computed
noProcess <- no_processesAtEnd %>%
filter(process == "analysis.AnalyseLikelihoodOfPresence") %>%
select(no_java) %>%
unlist(use.names = FALSE)
logL = .jcall("analysis.EasyRun", "D", "getSingleValueFromAquanismGroupProcess",
"processes.processesAtEnd.0.getLogLikelihood");
paste0("processes.processesAtEnd.", noProcess,"processes.processesAtEnd.0.getLogLikelihood"))
# 4. proportion of first-time spawners
# annual % of first-time spawners according to latitude (or basin_id)
......@@ -92,8 +108,13 @@ exploreEasyRun = function(simulationName){
# 5. mean age of female and age
# annual mean age over the period (1920 and 1950) in each basin,
# and average over basin
noProcess <- no_processesEachStep %>%
filter(process == "analysis.AnalyseSpawnerFeatures") %>%
select(no_java) %>%
unlist(use.names = FALSE)
out_spawners_features = .jcall("analysis.EasyRun","[[Ljava/lang/String;","getValuesFromAquanismGroupProcess",
"processes.processesEachStep.9.exportToR")
paste0("processes.processesEachStep.", noProcess,".exportToR"))
spawners_features = suppressMessages (convertFromSALtoR(out_spawners_features))
slopePrimiparous = lm(pct_primiparous~latitude, data = spawners_features)$coeff[2]
......@@ -103,14 +124,19 @@ exploreEasyRun = function(simulationName){
# ------------------------------------------------------------ #
# 6. number of metapopulation based on the exchange matrix ----
# mean annual exchange matrix between 1920 and 1950
noProcess <- no_processesEachStep %>%
filter(process == "analysis.ExportExchangesBetweenCatchments") %>%
select(no_java) %>%
unlist(use.names = FALSE)
out_exchanges = .jcall("analysis.EasyRun","[[Ljava/lang/String;","getValuesFromAquanismGroupProcess",
"processes.processesEachStep.16.exportToR")
paste0("processes.processesEachStep.", noProcess,".exportToR"))
exchanges = suppressMessages(convertFromSALtoR(out_exchanges))
# keep only basin with a totalRun (sum over migrationBasin) >0
selectedBasin <- exchanges %>% group_by(basin = migrationBasin) %>%
summarise(totalRun = sum(effective), .groups = 'drop') %>%
filter(totalRun > 0) %>% select(basin)
if (nrow(selectedBasin > 0)) {
metapopulations <- exchanges %>% inner_join(selectedBasin, by = c('migrationBasin' = 'basin')) %>%
inner_join(selectedBasin, by=c('natalBasin' = 'basin')) %>%
......@@ -122,9 +148,10 @@ exploreEasyRun = function(simulationName){
summarise(n = n()) %>%
unlist(use.names = FALSE)
} else
nbMetapopulation =NA
nbMetapopulation = NA
return(tibble(simulationName = simulationName,
lastSpawnerAdundance=lastSpawnerAdundance,
lastSpawnerAdundance = lastSpawnerAdundance,
northern_basin_id = range$northern_basin_id,
effective_centroid_latitude = range$effective_centroid_latitude,
logL = logL,
......@@ -236,7 +263,11 @@ exploreEasyRun = function(simulationName){
# titi= c(0, 0.75)
# dim(titi)
#==================== run simulation 2 =====
runSimulation2 = function(simulationParameter, simSDuration = 600) {
runSimulation2 = function(simulationParameter, no_processesEachStep, no_processesAtEnd, simSDuration = 600) {
# ----------------------------------------------------------------------
cat(simulationParameter$name, ' = ', simulationParameter$longName,'\n')
# ---------------------------------------------------------- start the VM
jarfile = "target/GR3D-3.2-SNAPSHOT.jar"
......@@ -269,14 +300,14 @@ runSimulation2 = function(simulationParameter, simSDuration = 600) {
#toc()
#tic("simulationExploration")
result = exploreEasyRun(simulationParameter$name)
result = exploreEasyRun(simulationParameter$name, no_processesEachStep, no_processesAtEnd)
#toc()
return(result)
}
runSimulation2(simulationList[[35]])
runSimulation2(simulationList[[35]], no_processesEachStep = no_processesEachStep, no_processesAtEnd = no_processesAtEnd)
#resSimulation = lapply(simulationList, runSimulation2)
......@@ -289,7 +320,8 @@ runSimulation2(simulationList[[35]])
#registerDoParallel(cores = 1)
stopImplicitCluster()
cl <- makeCluster(3)
#cl <- makeCluster(3)
cl <- detectCores() %>% -1 %>% makeCluster
registerDoParallel(cl)
rm(resSimulations)
......@@ -297,10 +329,12 @@ rm(resSimulations)
tic('parralel')
resSimulations = foreach(simul = c(1,7,13,19),
.combine = rbind, .verbose = FALSE, .multicombine = TRUE,
.packages = c('tidyverse', 'rJava', 'jdx')
packages = c('tidyverse', 'rJava', 'jdx'),
.inorder = FALSE) %dopar% {
res = runSimulation2(simulationList[[simul]],
)
no_processesEachStep,
no_processesAtEnd)
)
#res = runSimulation3(simul)
}
toc()
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment