Spatial stochasticity analysis pipeline

Data cleaning directory:

Produces usable data files joined from raw csv files. Designed to be executed from the spatial_infection_kernel directory

Run ./data_cleaning/incidence.R and ./data_cleaning/larvae_density.R to produce all necessary files in ./data_shared/field_output Run ./data_cleaning/verification.R to ensure data matches across joins and tables for use in writing code

Model exec package:

Files to write, run, and store output of the spatial epizootic models Designed to be executed from the ./model_exec_pkg directory

To generate branch-, tree- and site- level linesearch, mcmc and lppd algorithm codes: run ./model_exec_pkg/models_metacode/wrap_writer.R with grain of 0, 1, and 2 respectively.

To run a given model "model.c":

Run the following code into the cmd line to make the executable

module load openmpi gcc gsl
mpicc ./models/model.c -o ./execs/model `gsl-config --cflags --libs`
mpiexec -np $N$ ./execs/model $1 $2 $3 $4

where

• $N$ is equal to (number of sites * reps per site) and reps per site is of the form 2^c for pairwise averaging lhood scores
• Input $1 is grid depth d (distance from center to perimeter) and grid cells in the lattice x = 3d^2 - 3d + 1
• Input $2 is the model name ( see ./model_exec_pkg/models_metacode/model_key.R for more information )
• Input $3 is the likelihood function ( 0 for kernel-binomial, 1 for kbinom * pois(larval density), more functions avaibile and described in lppd.h )
• Input $4 is reps per process, for total number of realizations per likelihood score evaluation to be ($N$ * $4 / number_of_sites)

The code is set up to only take in realizations per lhood score that are a power of 2 to reduce rounding error in calculation

Processing model fits

Note that these files are run from the project directory, not from the model_exec_pkg folder

After running ./model_exec_pkg/execs/linearch_grain, run

module load R
Rscript ./model_output_mgmt/linesearch_pca.R

to generate the PCA objects used in generating the proposal distribution for MCMC. This dynamically selects the number of linesearch iterations to include by taking the best log posterior probability and finding its difference to the median log probability

After running ./model_output_mgmt/execs/mcmc_$grain, run

module load R
Rscript ./model_output_mgmt/iat_informed_target_pca.R

to store the MCMC chains in the particular target proposal cycle and recalibrate the proposal distribution for next iteration. Set cycle_files <- F in order to do troubleshooting without overwriting file management and directories. This is done by selecting the top performing samples using ELPD WAIC in waic_informed_target_pca.R

We

• (1) calculate ACF,
• (2) identify the number of MCMC iterations that consitute one independent sample, and
• (3) calculate an ELPD based on site-level likelihood scores for all iterations in that sample

With slight abuse of model comparison,

• we select iterations in the top q% of ELPD scores or
• keep all scores within n points of the best score,
• whichever retains the most samples

In iat_informed_processing.R, we skip the first two independent samples using median IAT (integrated autocorrelation time). Median is taken across chains and parameter values but is model specific.

For both (MCMC and LS) processing code scripts, no information is overwritten. Checks made to ensure that PCA is only done for models in current stage of processing

Software versions

Utils replicable versions on RCC Midway:

• gsl - default 2.7

• openmpi - default 2.0.1

• gcc - default 6.1.0 or 6.2.0

• intel - default 16.0 or 16.917

• pgi - default 2017

Utils replicable versions on CRI Randi

• gcc - 12.1.0

• openmpi - 4.1.4

• gsl - 2.7.1

• gcc - 11.3.0

• mpich - 4.0.3-srun

• gsl - 2.7.1