Skip to content

Commit

Permalink
Multi host2 (#188)
Browse files Browse the repository at this point in the history
* add error messages for checking multi-host data inputs

* add checks for multihost data inputs

* switch from infected and host file to file_list to capture that multiple hosts are being used as inputs.

* replace mortality_on, time_lag, and rate with pest_host_table which contains these plus susceptibility.

* update docs

* update mortality on for pops_model.

* add competecy-table parameter.

* update docs

* accept multi host in pops_model_cpp

* fix typo

* fix typo2

* add pest host table verification and setup for input into cpp

* add pest host table examples for tests

* add competency tables for tests.

* add checks for competency table

* update configuration to create competency table

* update pops_model.R to use multihost

* draw competency table for each model run in functions

* add checks for values in pest_host_table and competency_table

* edit multihost checks

* remove unused variables

* fix typos and lint errors

* update config and helpers to correctly create host pools

* fix helpers function

* make expose_file exposed_file_list

* update mortality_on

* update checks, config, and helpers to fix input to cpp issues.

* use create host pool function in cal, val, pops, and pops_multi

* update testthat for pops

* Store the casted matrix inputs in std::vector to preserve them

* Remove extra reserve

* Read pest-host table from parameters

* add non-mortality version of pesthost table

* update tests to include new tables

* Use the existing matrices used in HostPool for cloning instead of getting them from inputs

* add more test data for mortality tests

* fix exposed issue

* update tests

* fix mortality export in pops.cpp

* update validation to loop through and combine and compare infections for all hosts to infected locations.

* update validation tests

* update validate and tests

* update calibrate and tests

* update popsmultirun to handle multihost api from pops core

* update tests for for pops multirun

* lint

* fix notes in multirun

* add parameters to config tests

* fix mask error in calibrate

* modify order in validate stats

* add check to ensure that all host populations are less than total populations in all cells

* add error messages

* update if to while

* add data for tests

* update tests

* add returns for config errors

* add data for 2 host

* Store references to matrices only after the objects have fixed addresses (i.e., vector size is stable)

* lint

* lint

* update data

* fix config test error

* add uncertainty test data

* add tests for multihost

* add set.seed

* update multihost tests

* export output host pools from all runs

* fix cal error

---------

Co-authored-by: Anna Petrasova <kratochanna@gmail.com>
Co-authored-by: Vaclav Petras <wenzeslaus@gmail.com>
  • Loading branch information
3 people authored Jan 30, 2024
1 parent ac1383a commit b70235a
Show file tree
Hide file tree
Showing 48 changed files with 3,261 additions and 2,186 deletions.
4 changes: 2 additions & 2 deletions R/RcppExports.R
Original file line number Diff line number Diff line change
@@ -1,8 +1,8 @@
# Generated by using Rcpp::compileAttributes() -> do not edit by hand
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393

pops_model_cpp <- function(random_seed, multiple_random_seeds, random_seeds, lethal_temperature, lethal_temperature_month, infected, total_exposed, exposed, susceptible, total_populations, total_hosts, mortality_tracker, mortality, quarantine_areas, quarantine_directions, treatment_maps, treatment_dates, pesticide_duration, resistant, movements, movements_dates, temperature, survival_rates, weather_coefficient, weather_coefficient_sd, bbox, res, rows_cols, soil_reservoirs, reproductive_rate, spatial_indices, season_month_start_end, frequency_config, bool_config, mortality_rate = 0.0, mortality_time_lag = 2L, start_date = "2018-01-01", end_date = "2018-12-31", treatment_method = "ratio", natural_kernel_type = "cauchy", anthropogenic_kernel_type = "cauchy", percent_natural_dispersal = 0.0, natural_distance_scale = 21, anthropogenic_distance_scale = 0.0, natural_dir = "NONE", natural_kappa = 0, anthropogenic_dir = "NONE", anthropogenic_kappa = 0, frequencies_n_config = NULL, model_type_ = "SI", latency_period = 0L, establishment_probability = 0, dispersal_percentage = 0.99, survival_rate_month = 0L, survival_rate_day = 0L, overpopulation_config = NULL, network_config = NULL, network_data_config = NULL, weather_size = 0L, weather_type = "deterministic", dispersers_to_soils_percentage = 0) {
.Call(`_PoPS_pops_model_cpp`, random_seed, multiple_random_seeds, random_seeds, lethal_temperature, lethal_temperature_month, infected, total_exposed, exposed, susceptible, total_populations, total_hosts, mortality_tracker, mortality, quarantine_areas, quarantine_directions, treatment_maps, treatment_dates, pesticide_duration, resistant, movements, movements_dates, temperature, survival_rates, weather_coefficient, weather_coefficient_sd, bbox, res, rows_cols, soil_reservoirs, reproductive_rate, spatial_indices, season_month_start_end, frequency_config, bool_config, mortality_rate, mortality_time_lag, start_date, end_date, treatment_method, natural_kernel_type, anthropogenic_kernel_type, percent_natural_dispersal, natural_distance_scale, anthropogenic_distance_scale, natural_dir, natural_kappa, anthropogenic_dir, anthropogenic_kappa, frequencies_n_config, model_type_, latency_period, establishment_probability, dispersal_percentage, survival_rate_month, survival_rate_day, overpopulation_config, network_config, network_data_config, weather_size, weather_type, dispersers_to_soils_percentage)
pops_model_cpp <- function(random_seed, multiple_random_seeds, random_seeds, lethal_temperature, lethal_temperature_month, host_pools, total_populations, competency_table, pest_host_table, quarantine_areas, quarantine_directions, treatment_maps, treatment_dates, pesticide_duration, movements, movements_dates, temperature, survival_rates, weather_coefficient, weather_coefficient_sd, bbox, res, rows_cols, soil_reservoirs, reproductive_rate, spatial_indices, season_month_start_end, frequency_config, bool_config, start_date = "2018-01-01", end_date = "2018-12-31", treatment_method = "ratio", natural_kernel_type = "cauchy", anthropogenic_kernel_type = "cauchy", percent_natural_dispersal = 0.0, natural_distance_scale = 21, anthropogenic_distance_scale = 0.0, natural_dir = "NONE", natural_kappa = 0, anthropogenic_dir = "NONE", anthropogenic_kappa = 0, frequencies_n_config = NULL, model_type_ = "SI", latency_period = 0L, establishment_probability = 0, dispersal_percentage = 0.99, survival_rate_month = 0L, survival_rate_day = 0L, overpopulation_config = NULL, network_config = NULL, network_data_config = NULL, weather_size = 0L, weather_type = "deterministic", dispersers_to_soils_percentage = 0) {
.Call(`_PoPS_pops_model_cpp`, random_seed, multiple_random_seeds, random_seeds, lethal_temperature, lethal_temperature_month, host_pools, total_populations, competency_table, pest_host_table, quarantine_areas, quarantine_directions, treatment_maps, treatment_dates, pesticide_duration, movements, movements_dates, temperature, survival_rates, weather_coefficient, weather_coefficient_sd, bbox, res, rows_cols, soil_reservoirs, reproductive_rate, spatial_indices, season_month_start_end, frequency_config, bool_config, start_date, end_date, treatment_method, natural_kernel_type, anthropogenic_kernel_type, percent_natural_dispersal, natural_distance_scale, anthropogenic_distance_scale, natural_dir, natural_kappa, anthropogenic_dir, anthropogenic_kappa, frequencies_n_config, model_type_, latency_period, establishment_probability, dispersal_percentage, survival_rate_month, survival_rate_day, overpopulation_config, network_config, network_data_config, weather_size, weather_type, dispersers_to_soils_percentage)
}

# Register entry points for exported C++ functions
Expand Down
133 changes: 53 additions & 80 deletions R/calibrate.R
Original file line number Diff line number Diff line change
Expand Up @@ -123,8 +123,10 @@ calibrate <- function(infected_years_file,
params_to_estimate = c(TRUE, TRUE, TRUE, TRUE, FALSE, FALSE),
number_of_generations = 7,
generation_size = 1000,
infected_file,
host_file,
pest_host_table,
competency_table,
infected_file_list,
host_file_list,
total_populations_file,
temp = FALSE,
temperature_coefficient_file = "",
Expand All @@ -145,9 +147,6 @@ calibrate <- function(infected_years_file,
temperature_file = "",
lethal_temperature = -12.87,
lethal_temperature_month = 1,
mortality_on = FALSE,
mortality_rate = 0,
mortality_time_lag = 0,
mortality_frequency = "year",
mortality_frequency_n = 1,
management = FALSE,
Expand Down Expand Up @@ -183,7 +182,7 @@ calibrate <- function(infected_years_file,
leaving_scale_coefficient = 1,
calibration_method = "ABC",
number_of_iterations = 100000,
exposed_file = "",
exposed_file_list = "",
verbose = TRUE,
write_outputs = "None",
output_folder_path = "",
Expand Down Expand Up @@ -213,8 +212,8 @@ calibrate <- function(infected_years_file,
config$params_to_estimate <- params_to_estimate
config$number_of_generations <- number_of_generations
config$generation_size <- generation_size
config$infected_file <- infected_file
config$host_file <- host_file
config$infected_file_list <- infected_file_list
config$host_file_list <- host_file_list
config$total_populations_file <- total_populations_file
config$temp <- temp
config$temperature_coefficient_file <- temperature_coefficient_file
Expand All @@ -235,9 +234,6 @@ calibrate <- function(infected_years_file,
config$survival_rate_month <- survival_rate_month
config$survival_rate_day <- survival_rate_day
config$survival_rates_file <- survival_rates_file
config$mortality_on <- mortality_on
config$mortality_rate <- mortality_rate
config$mortality_time_lag <- mortality_time_lag
config$management <- management
config$treatment_dates <- treatment_dates
config$treatments_file <- treatments_file
Expand Down Expand Up @@ -272,7 +268,7 @@ calibrate <- function(infected_years_file,
config$leaving_scale_coefficient <- leaving_scale_coefficient
config$calibration_method <- calibration_method
config$number_of_iterations <- number_of_iterations
config$exposed_file <- exposed_file
config$exposed_file_list <- exposed_file_list
# add function name for use in configuration function to skip
# function specific specific configurations namely for validation and
# calibration.
Expand All @@ -296,6 +292,8 @@ calibrate <- function(infected_years_file,
config$use_soils <- use_soils
config$soil_starting_pest_file <- soil_starting_pest_file
config$start_with_soil_populations <- start_with_soil_populations
config$pest_host_table <- pest_host_table
config$competency_table <- competency_table

# call configuration function to perform data checks and transform data into
# format used in pops c++
Expand Down Expand Up @@ -323,49 +321,16 @@ calibrate <- function(infected_years_file,
network_min_distance,
network_max_distance) {

config$random_seed <- sample(1:999999999999, 1, replace = FALSE)
config$random_seed <- as.integer(sample.int(1e9, 1, replace = FALSE))
set.seed(config$random_seed[[1]])
random_seeds <- create_random_seeds(1)
if (config$use_initial_condition_uncertainty) {
config$infected <- matrix_norm_distribution(config$infected_mean, config$infected_sd)
while (any(config$infected < 0)) {
config$infected <- matrix_norm_distribution(config$infected_mean, config$infected_sd)
}
exposed2 <- matrix_norm_distribution(config$exposed_mean, config$exposed_sd)
while (any(exposed2 < 0)) {
exposed2 <- matrix_norm_distribution(config$infected_mean, config$infected_sd)
}
exposed <- config$exposed
exposed[[config$latency_period + 1]] <- exposed2
config$exposed <- exposed
} else {
config$infected <- config$infected_mean
exposed2 <- config$exposed_mean
exposed <- config$exposed
exposed[[config$latency_period + 1]] <- exposed2
config$exposed <- exposed
}

if (config$use_host_uncertainty) {
config$host <- matrix_norm_distribution(config$host_mean, config$host_sd)
while (any(config$host > config$total_populations)) {
config$host <- matrix_norm_distribution(config$host_mean, config$host_sd)
}
} else {
config$host <- config$host_mean
}

susceptible <- config$host - config$infected - exposed2
susceptible[susceptible < 0] <- 0

config$susceptible <- susceptible
config$total_hosts <- config$host
config$total_exposed <- exposed2

if (config$mortality_on) {
mortality_tracker2 <- config$mortality_tracker
mortality_tracker2[[length(mortality_tracker2)]] <- config$infected
config$mortality_tracker <- mortality_tracker2
config <- host_pool_setup(config)
while (any(config$total_hosts > config$total_populations) ||
any(config$total_exposed > config$total_populations) ||
any(config$total_infecteds > config$total_populations)) {
config <- host_pool_setup(config)
}
config$competency_table_list <- competency_table_list_creator(config$competency_table)

data <- pops_model(
random_seed = config$random_seed,
Expand All @@ -377,21 +342,16 @@ calibrate <- function(infected_years_file,
use_survival_rates = config$use_survival_rates,
survival_rate_month = config$survival_rate_month,
survival_rate_day = config$survival_rate_day,
infected = config$infected,
total_exposed = config$total_exposed,
exposed = config$exposed,
susceptible = config$susceptible,
host_pools = config$host_pools,
total_populations = config$total_populations,
total_hosts = config$total_hosts,
competency_table = config$competency_table_list,
pest_host_table = config$pest_host_table_list,
mortality_on = config$mortality_on,
mortality_tracker = config$mortality_tracker,
mortality = config$mortality,
quarantine_areas = config$quarantine_areas,
quarantine_directions = config$quarantine_directions,
treatment_maps = config$treatment_maps,
treatment_dates = config$treatment_dates,
pesticide_duration = config$pesticide_duration,
resistant = config$resistant,
use_movements = config$use_movements,
movements = config$movements,
movements_dates = config$movements_dates,
Expand All @@ -407,8 +367,6 @@ calibrate <- function(infected_years_file,
spatial_indices = config$spatial_indices,
season_month_start_end = config$season_month_start_end,
soil_reservoirs = config$soil_reservoirs,
mortality_rate = config$mortality_rate,
mortality_time_lag = config$mortality_time_lag,
start_date = config$start_date,
end_date = config$end_date,
treatment_method = config$treatment_method,
Expand Down Expand Up @@ -590,16 +548,21 @@ calibrate <- function(infected_years_file,
# the simulation
all_disagreement <-
foreach::foreach(
q = seq_len(length(data$infected)),
q = seq_len(length(data$host_pools[[1]]$infected)),
.combine = rbind,
.packages = c("terra", "PoPS"),
.final = colSums
) %do% {
comparison <- terra::rast(config$infected_file)[[1]]
reference <- terra::rast(config$infected_file)[[1]]
terra::values(comparison) <- data$infected[[q]]
comparison <- terra::rast(config$infected_file_list[[1]])[[1]]
reference <- comparison
mask <- comparison
terra::values(comparison) <- 0
infections <- comparison
for (p in seq_len(length(data$host_pools))) {
terra::values(infections) <- data$host_pools[[p]]$infected[[q]]
comparison <- comparison + infections
}
terra::values(reference) <- config$infection_years2[[q]]
mask <- terra::rast(config$infected_file)[[1]]
terra::values(mask) <- config$mask_matrix
quantity_allocation_disagreement(reference,
comparison,
Expand All @@ -609,7 +572,7 @@ calibrate <- function(infected_years_file,
}

all_disagreement <- as.data.frame(t(all_disagreement))
all_disagreement <- all_disagreement / length(data$infected)
all_disagreement <- all_disagreement / length(data$host_pools[[1]]$infected)
config$quantity <- all_disagreement$quantity_disagreement
config$allocation <- all_disagreement$allocation_disagreement
config$configuration_dis <- all_disagreement$configuration_disagreement
Expand Down Expand Up @@ -928,16 +891,21 @@ calibrate <- function(infected_years_file,

all_disagreement <-
foreach::foreach(
q = seq_len(length(data$infected)),
q = seq_len(length(data$host_pools[[1]]$infected)),
.combine = rbind,
.packages = c("terra", "PoPS"),
.final = colSums
) %do% {
comparison <- terra::rast(config$infected_file)[[1]]
reference <- terra::rast(config$infected_file)[[1]]
terra::values(comparison) <- data$infected[[q]]
comparison <- terra::rast(config$infected_file_list[[1]])[[1]]
reference <- comparison
mask <- comparison
terra::values(comparison) <- 0
infections <- comparison
for (p in seq_len(length(data$host_pools))) {
terra::values(infections) <- data$host_pools[[p]]$infected[[q]]
comparison <- comparison + infections
}
terra::values(reference) <- config$infection_years2[[q]]
mask <- terra::rast(config$infected_file)[[1]]
terra::values(mask) <- config$mask_matrix
quantity_allocation_disagreement(reference,
comparison,
Expand All @@ -947,7 +915,7 @@ calibrate <- function(infected_years_file,
}

all_disagreement <- as.data.frame(t(all_disagreement))
all_disagreement <- all_disagreement / length(data$infected)
all_disagreement <- all_disagreement / length(data$host_pools[[1]]$infected)
config$accuracy <- all_disagreement$accuracy
config$precision <- all_disagreement$precision
config$recall <- all_disagreement$recall
Expand Down Expand Up @@ -1102,26 +1070,31 @@ calibrate <- function(infected_years_file,
# set up comparison
all_disagreement <-
foreach::foreach(
q = seq_len(length(data$infected)),
q = seq_len(length(data$host_pools[[1]]$infected)),
.combine = rbind,
.packages = c("terra", "PoPS"),
.final = colSums
) %do% {
comparison <- terra::rast(config$infected_file)[[1]]
comparison <- terra::rast(config$infected_file_list[[1]])[[1]]
reference <- comparison
mask <- comparison
terra::values(comparison) <- data$infected[[q]]
terra::values(comparison) <- 0
infections <- comparison
for (p in seq_len(length(data$host_pools))) {
terra::values(infections) <- data$host_pools[[p]]$infected[[q]]
comparison <- comparison + infections
}
terra::values(reference) <- config$infection_years2[[q]]
terra::values(mask) <- config$mask_matrix
quantity_allocation_disagreement(reference,
comparison,
use_configuration = FALSE,
use_configuration = config$use_configuration,
mask = mask,
use_distance = config$use_distance)
}

all_disagreement <- as.data.frame(t(all_disagreement))
all_disagreement <- all_disagreement / length(data$infected)
all_disagreement <- all_disagreement / length(data$host_pools[[1]]$infected)
proposed <-
data.frame(all_disagreement[, c("quantity_disagreement", "allocation_disagreement",
"configuration_disagreement", "accuracy", "precision",
Expand Down
Loading

0 comments on commit b70235a

Please sign in to comment.