Skip to contents

Likelihood-Free Markhov Chain Monte Carlo (LFMCMC)

Usage

LFMCMC(model = NULL)

run_lfmcmc(lfmcmc, params_init_, n_samples_, epsilon_, seed = NULL)

set_observed_data(lfmcmc, observed_data_)

set_proposal_fun(lfmcmc, fun)

use_proposal_norm_reflective(lfmcmc)

set_simulation_fun(lfmcmc, fun)

set_summary_fun(lfmcmc, fun)

set_kernel_fun(lfmcmc, fun)

use_kernel_fun_gaussian(lfmcmc)

set_par_names(lfmcmc, names)

set_stats_names(lfmcmc, names)

get_params_mean(lfmcmc)

get_stats_mean(lfmcmc)

# S3 method for class 'epiworld_lfmcmc'
print(x, ...)

Arguments

model

A model of class epiworld_model or NULL (see details).

lfmcmc

LFMCMC model

params_init_

Initial model parameters, treated as double

n_samples_

Number of samples, treated as integer

epsilon_

Epsilon parameter, treated as double

seed

Random engine seed

observed_data_

Observed data, treated as double

fun

The LFMCMC kernel function

names

The model stats names

x

LFMCMC model to print

...

Ignored

Value

The LFMCMC function returns a model of class epiworld_lfmcmc.

The simulated model of class epiworld_lfmcmc.

The lfmcmc model with the observed data added

The lfmcmc model with the proposal function added

The LFMCMC model with proposal function set to norm reflective

The lfmcmc model with the simulation function added

The lfmcmc model with the summary function added

The lfmcmc model with the kernel function added

The LFMCMC model with kernel function set to gaussian

The lfmcmc model with the parameter names added

The lfmcmc model with the stats names added

The param means for the given lfmcmc model

The stats means for the given lfmcmc model

The lfmcmc model

Details

Performs a Likelihood-Free Markhov Chain Monte Carlo simulation. When model is not NULL, the model uses the same random-number generator engine as the model. Otherwise, when model is NULL, a new random-number generator engine is created.

Examples

## Setup an SIR model to use in the simulation
model_seed <- 122
model_sir <- ModelSIR(name = "COVID-19", prevalence = .1,
  transmission_rate = .9, recovery_rate = .3)
agents_smallworld(
  model_sir,
  n = 1000,
  k = 5,
  d = FALSE,
  p = 0.01
)
verbose_off(model_sir)
run(model_sir, ndays = 50, seed = model_seed)

## Setup LFMCMC
# Extract the observed data from the model
obs_data <- get_today_total(model_sir)

# Define the simulation function
simfun <- function(params) {
  set_param(model_sir, "Recovery rate", params[1])
  set_param(model_sir, "Transmission rate", params[2])
  run(model_sir, ndays = 50)
  res <- get_today_total(model_sir)
  return(res)
}

# Define the summary function
sumfun <- function(dat) {
  return(dat)
}

# Create the LFMCMC model
lfmcmc_model <- LFMCMC(model_sir) |>
  set_simulation_fun(simfun) |>
  set_summary_fun(sumfun) |>
  use_proposal_norm_reflective() |>
  use_kernel_fun_gaussian() |>
  set_observed_data(obs_data)

## Run LFMCMC simulation
# Set initial parameters
par0 <- c(0.1, 0.5)
n_samp <- 2000
epsil <- 1.0

# Run the LFMCMC simulation
run_lfmcmc(
  lfmcmc = lfmcmc_model,
  params_init_ = par0,
  n_samples_ = n_samp,
  epsilon_ = epsil,
  seed = model_seed
)

# Print the results
set_stats_names(lfmcmc_model, get_states(model_sir))
set_par_names(lfmcmc_model, c("Immune recovery", "Infectiousness"))

print(lfmcmc_model)
#> ___________________________________________
#> 
#> LIKELIHOOD-FREE MARKOV CHAIN MONTE CARLO
#> 
#> N Samples : 2000
#> Elapsed t : 1.00s
#> 
#> Parameters:
#>   -Immune recovery :  0.45 [ 0.14,  0.95] (initial :  0.10)
#>   -Infectiousness  :  0.85 [ 0.54,  1.00] (initial :  0.50)
#> 
#> Statistics:
#>   -Susceptible :     0.22 [    0.00,     2.00] (Observed:     0.00)
#>   -Infected    :     0.07 [    0.00,     1.00] (Observed:     0.00)
#>   -Recovered   :   995.71 [  998.00,  1000.00] (Observed:  1000.00)
#> ___________________________________________
#> 

get_stats_mean(lfmcmc_model)
#> [1]   0.2215   0.0655 995.7130
get_params_mean(lfmcmc_model)
#> [1] 0.4481647 0.8513273