2-RunHBVmodel.Rmd

# Run HBV model 
# =============

Authors: Richard T. Gray, Neil Bretana

This script is used to load a project and run in the HBV model. 

```{r Initialize script}
# Clear workspace
rm(list = ls()) 

# Setup directories after setting working directory to source file 
# directory (so it is easy to move code around)
basePath <- getwd()
dataPath <- file.path(basePath, "data") 
Rcode <- file.path(basePath, "code")
projectFolder <- file.path(basePath, "projects")

# Load useful libraries
source(file.path(Rcode, "LoadLibrary.R"))
source(file.path(Rcode, "DataLibraries.R"))
# source(file.path(Rcode, "PlotOptions.R"))

# Source model files
source(file.path(Rcode, "HBVmodel.R"))

# Current date and time for appending to results
currTime <- format(Sys.time(), "%Y-%m-%d_%H-%M-%S")

```

## User specifications

```{r User specifications}
# Specify project to run
project <- "Australia_HBV_estimates-BM"

# Options
runSamples <- TRUE
saveAsBase <- TRUE # if TRUE doesn't append time to results and overwrites
                    # a base file. Useful for storing main results or
                    # testing

```

## Main script

The first chunk loads all the project specifications and inputs.

```{r Load project}
# This loads pg, inputs, best_estimates, param_sets, best_initial_pop
# param_initial_pops, transitions to the workspace.

projectFile <- file.path(projectFolder, project, paste0(project, ".rda"))
projectVars <- load(projectFile)

```

The next chunk runs the all the simulations in a project.

```{r Run model}

# Start clock
tic <- proc.time()

# Run model on best estimates
bestResults <- HBVmodel(pg, best_estimates, best_initial_pop, pg$pts,
                        transitions, waifw_matrix)

# Run sampled parameter sets
if (runSamples) {
  paramResults <- list()
  
  for (set in 1:pg$parameter_samples){
    paramSet <- paste0("paramSet", toString(set))
    
    paramResults[[paramSet]] <- HBVmodel(pg, param_sets[[paramSet]], 
                                         param_initial_pops[[paramSet]],
                                         pg$pts, transitions, waifw_matrix)
  }
}

# Stop the clock
toc <- proc.time() - tic

```

```{r Save results}

resultsFolder <- file.path(projectFolder, project, "results")

if (saveAsBase) {
  # Create testing folder and results file
  dir.create(file.path(resultsFolder, "results_base"), 
             showWarnings = FALSE)
  resultsFile <- file.path(resultsFolder, "results_base", 
                           "results_base.rda")
} else {
  dir.create(file.path(resultsFolder, paste0("results_", currTime)), 
             showWarnings = FALSE)
  resultsFile <- file.path(resultsFolder, paste0("results_", currTime), 
                           paste0("results_", currTime, ".rda"))
}

# Append the results to the project file - annoyingly R makes this tricky!
if (runSamples) {
  save(list = c("bestResults", "paramResults", projectVars), 
       file = resultsFile)
} else {
  save(list = c("bestResults", projectVars), 
       file = resultsFile)
}

  
```