#!/usr/bin/env Rscript ############################################################################### # script pour lancer calibration sur Datarmor avec PBS qsub # toutes les variables d���environnement sont d��finies dans le script PBS ############################################################################### # normalement 2 packages charg��s: {Rmpi} et {snow} qui mettent �� dispo: # `makeCluster`, `clusterEvalQ`, `clusterExport`, `stopCluster` sessionInfo() library(doParallel, include.only = c("registerDoParallel")) # pour ��viter d�����craser `makeCluster` # la fa��on dont {calibrar} utilise {foreach} perturbe {doParallel} # il faut toujours charger tous les pkgs n��cessaires ici et dans `clusterEvalQ` dessous # si vous le faites pas en double, erreur ���function not found��� directe! library(calibrar, lib.loc = Sys.getenv("NICO_REPO")) library(dplyr, warn.conflicts = FALSE) library(data.table, include.only = c("fread")) library(stringi, include.only = c("stri_extract_last_regex", "stri_join")) options(dplyr.summarise.inform = FALSE) cl <- makeCluster() registerDoParallel(cl) clusterEvalQ(cl, { library(calibrar, lib.loc = Sys.getenv("NICO_REPO")) library(dplyr, warn.conflicts = FALSE) library(data.table, include.only = c("fread")) library(stringi, include.only = c("stri_extract_last_regex", "stri_join")) options(dplyr.summarise.inform = FALSE) }) %>% # cacher des messages encombrants suppressMessages() %>% suppressWarnings() %>% capture.output() %>% invisible() ############################################################################### runModel <- function(par) { # "par" provided by the calibrate function # run ISIS ====================================================== # ATTENTION we are in RUN_*/i* folder # if necessary, manually save q of each generation current_time <- format(Sys.time(), format = "%F %Hh%M", tz = "CET", usetz = TRUE) current_gen <- 1L + length(list.files(pattern = "^params .+\\.csv$")) info_suffix <- sprintf(" gen %03d start %s", current_gen, current_time) # get back the capturability parameters (par) and write them in the format the model is able to read (csv, semantic) textLines <- readLines(Sys.getenv("PARAMS_FILE")) for (i in seq_along(par)) { id <- 2 + i tmp <- strsplit(textLines[id], ";")[[1]] textLines[id] <- paste0(tmp[1], ";", par[i]) } write(textLines, file = Sys.getenv("PARAMS_FILE")) file.copy(from = Sys.getenv("PARAMS_FILE"), to = paste0("params", info_suffix, ".csv")) # launch ISIS from R, cf le script PBS current_ind <- basename(getwd()) # current individual of the population nomSimul_i <- Sys.getenv("SIMUL_BASENAME") %>% paste(., current_ind, sep = "_") # debug_file <- paste0("debugISIS", info_suffix, ".log") system2("myIsis", args = nomSimul_i, stdout = "debugISIS_raja2015.log", stderr = "debugISIS_raja2015.log") dossier_simul <- file.path(Sys.getenv("ISIS_BD"), "isis-database/simulations", nomSimul_i, "resultExports") # read the model outputs and tranform them into a list ========== # catch at age export1 <- file.path(dossier_simul, "CapturesPoids.csv") %>% fread(data.table = FALSE) %>% mutate(names = paste0("gp",group,"_Raja2015")) %>% select(names, value) # combine (if several inputs)======================================================= # export <- rbind(export1, export2) export <- export1 output <- as.list(export$value) names(output) <- export$names return(output) } ############################################################################### #' @title residual sum of squares function #' @description name of this fx must match what is declared in calibrationInfo.csv #' @param obs,sim numeric values of LENGTH 1 userFunction <- function(obs, sim) { # obs peut ��tre NA, c-a-d pas de valeurs observ��es mais valeurs simul��es existent dans les exports # d�� �� des combinaisons irr��elles entre flottille avec GSA zone pop res <- if (is.na(obs)) 0 else if (obs == 0) sim else (obs - sim) / obs return(res ^ 2) } # calibrationInfo calib_set <- calibration_setup( file = "data_calib/data_calib_Raja2015/calibrationInfo.csv", control = list(col_skip = 1) # pour compatibilit�� avec ancienne version ) # Observed data calib_obs <- calibration_data(setup = calib_set, path = "data_calib/data_calib_Raja2015") # creation of the objective function obj_fn <- calibration_objFn(model = runModel, setup = calib_set, observed = calib_obs) # aggFn = calibrar:::.weighted.sum # by default clusterExport(cl, c("userFunction", "obj_fn", "calib_set", "calib_obs")) ############################################################################### ncpus <- as.integer(Sys.getenv("mpiproc")) # control argument of the calibrate function control <- list( master = "master/master_raja2015", # a folder run = paste0("RUN_", Sys.getenv("SIMUL_BASENAME")), # a folder restart.file = paste0("ckpt_", Sys.getenv("SIMUL_BASENAME")), convergence = 0.001, maxgen = 500, popsize = if (ncpus < 50) 2 * ncpus else ncpus, REPORT = 1, trace = 5, verbose = TRUE, parallel = TRUE, nCores = ncpus ) # N param��tres calib <- calibrate( fn = obj_fn, method = "default", control = control, par = 9.87e-3, lower = 1e-5, upper = 1e-1 ) stopCluster(cl)