libspatialSEIR
0.1
Bayesian Spatial SEIR Modeling
|
#include <ModelContext.hpp>
Public Member Functions | |
ModelContext () | |
~ModelContext () | |
void | populate (InitData *_A0, covariateArgs *xArgs, covariateArgs *xPrsArgs, double *offset, int *Y, compartmentArgs *S_starArgs, compartmentArgs *E_starArgs, compartmentArgs *I_starArgs, compartmentArgs *R_starArgs, scaledDistanceArgs *scaledDistArgs, double *rho, double *phi, double *beta, double *gamma_ei, double *gamma_ir, double *betaPrs, int *N, sliceParameters *sliceWidths, priorControl *priorInformation, modelConfiguration _config) |
void | buildModel () |
void | configureIterationTasks () |
void | calculateS_CPU () |
void | calculateS_CPU (int startLoc, int startTime) |
void | calculateS_givenE_CPU () |
void | calculateS_givenE_CPU (int startLoc, int startTime) |
void | calculateS_OCL () |
void | calculateE_CPU () |
void | calculateE_CPU (int startLoc, int startTime) |
void | calculateE_givenI_CPU () |
void | calculateE_givenI_CPU (int startLoc, int startTime) |
void | calculateE_OCL () |
void | calculateI_CPU () |
void | calculateI_CPU (int startLoc, int startTime) |
void | calculateI_givenR_CPU () |
void | calculateI_givenR_CPU (int startLoc, int startTime) |
void | calculateI_OCL () |
void | calculateR_CPU () |
void | calculateR_CPU (int startLoc, int startTime) |
void | calculateR_givenS_CPU () |
void | calculateR_givenS_CPU (int startLoc, int startTime) |
void | calculateR_OCL () |
void | calculateGenericCompartment_CPU (CompartmentalModelMatrix *comp, int *comp0, CompartmentalModelMatrix *compStarAdd, CompartmentalModelMatrix *compStarSub) |
void | calculateGenericCompartment_CPU (CompartmentalModelMatrix *comp, int *comp0, CompartmentalModelMatrix *compStarAdd, CompartmentalModelMatrix *compStarSub, int startLoc, int startTime) |
void | calculateGenericCompartment_OCL (int *comp, int *comp0, int *compStarAdd, int *compStarSub) |
int | checkCompartmentBounds () |
void | printFCValues () |
void | setRandomSeed (unsigned int seedValue) |
void | simulationIter (bool verbose, bool debug) |
void | runSimulation (int nIterations, bool verbose, bool debug) |
void | updateSamplingParameters (double desiredRatio, double targetWidth, double proportionChange) |
void | setCompartmentSamplingMode (int mode) |
int | getCompartmentSamplingMode () |
void | setParameterSamplingMode (int mode) |
int | getParameterSamplingMode () |
void | cacheP_SE_Calculation () |
void | calculateP_SE_CPU () |
void | calculateP_SE_CPU (int startLoc, int startTime) |
void | calculateP_SE_OCL () |
void | calculateP_EI_CPU () |
void | calculateP_IR_CPU () |
void | calculateP_RS_CPU () |
int | totalS () |
int | totalE () |
int | totalI () |
int | totalR () |
int | totalS (int tpt) |
int | totalE (int tpt) |
int | totalI (int tpt) |
int | totalR (int tpt) |
int | totalS_star () |
int | totalE_star () |
int | totalI_star () |
int | totalR_star () |
int | totalS_star (int tpt) |
int | totalE_star (int tpt) |
int | totalI_star (int tpt) |
int | totalR_star (int tpt) |
double | avgP_SE () |
double | avgP_RS () |
double | avgP_SE (int tpt) |
double | estimateR0 () |
double * | estimateR0 (int t) |
double | estimateEffectiveR0 () |
double * | estimateEffectiveR0 (int t) |
double * | calculateR0Components (int t) |
double * | calculateEffectiveR0Components (int t) |
double * | calculateG (int t) |
double * | calculateIntegratedG (int t) |
The model context class provides the central interface of libspatialSEIR. After instantiation, the populate method must be called before any meaningful work can be done. Once populated, ModelContext holds pointers to all of the required FullConditional instances as well as the data they use. In addition, ModelContext provides calculation functions and holds pointers to additional utility classes such as OCLProvider, IOProvider, and RandomNumberProvider.
SpatialSEIR::ModelContext::ModelContext | ( | ) |
SpatialSEIR::ModelContext::~ModelContext | ( | ) |
double SpatialSEIR::ModelContext::avgP_RS | ( | ) |
Calculates the average of all components of p_rs.
double SpatialSEIR::ModelContext::avgP_SE | ( | ) |
Calculates the average of all components of p_se.
double SpatialSEIR::ModelContext::avgP_SE | ( | int | tpt | ) |
Calculates the average of all components of p_se for a particular time point..
void SpatialSEIR::ModelContext::buildModel | ( | ) |
void SpatialSEIR::ModelContext::cacheP_SE_Calculation | ( | ) |
cacheP_SE_Calculation places the pre-requisites for calculating p_se in p_se_components. This is useful when partial updates can be performed by calculateP_SE_CPU(int startLoc, int startTime)
void SpatialSEIR::ModelContext::calculateE_CPU | ( | ) |
The zero parameter overload of calculateE_CPU calculates the E compartment from A0, I_star, and E_star.
void SpatialSEIR::ModelContext::calculateE_CPU | ( | int | startLoc, |
int | startTime | ||
) |
The two parameter overload of calculateE_CPU again calculates the E compartment from A0, I_star and E_star, but assumes that only values of I_star and E_star at column (location) startLoc and after row (time) startTime have changed since E was last calculated. When applicable, this saves considerable time.
void SpatialSEIR::ModelContext::calculateE_givenI_CPU | ( | ) |
The zero parameter overload of calculate E_givenI_CPU takes advantage of the fact that N=S+E+I+R to update E when I changes.
void SpatialSEIR::ModelContext::calculateE_givenI_CPU | ( | int | startLoc, |
int | startTime | ||
) |
The two parameter overload of calculate E_givenI_CPU takes advantage of the fact that N=S+E+I+R to update E when I changes, while assuming that I has only changed for location startLoc and after time point startTime.
void SpatialSEIR::ModelContext::calculateE_OCL | ( | ) |
calculateE_OCL works identically to calculateS_CPU using calls to oclProvider to perform computations in parallel.
double * SpatialSEIR::ModelContext::calculateEffectiveR0Components | ( | int | t | ) |
Calculates the location specific components of the effective R0 quantity at time t
double * SpatialSEIR::ModelContext::calculateG | ( | int | t | ) |
Calculates the next generation matrix
void SpatialSEIR::ModelContext::calculateGenericCompartment_CPU | ( | CompartmentalModelMatrix * | comp, |
int * | comp0, | ||
CompartmentalModelMatrix * | compStarAdd, | ||
CompartmentalModelMatrix * | compStarSub | ||
) |
comp | pointer to CompartmentalModelMatrix to be calculated. |
comp0 | pointer to integer array containing starting values for comp. |
compStarAdd | pointer to transition matrix into comp. |
compStarSub | pointer to transition CompartmentalModelMatrix out of comp. |
void SpatialSEIR::ModelContext::calculateGenericCompartment_CPU | ( | CompartmentalModelMatrix * | comp, |
int * | comp0, | ||
CompartmentalModelMatrix * | compStarAdd, | ||
CompartmentalModelMatrix * | compStarSub, | ||
int | startLoc, | ||
int | startTime | ||
) |
The six parameter overload of calculateGenericCompartment_CPU takes advantege of the common calculation pattern for S,E,I, and R to provide a single function mapping starting and transition values to the filan CompartmentalModelMatrix values, assuming that only location startLoc and time points after startTime need be considered.
comp | pointer to CompartmentalModelMatrix to be calculated |
comp0 | pointer to integer array containing starting values for comp |
compStarAdd | pointer to transition matrix into comp |
compStarSub | pointer to transition CompartmentalModelMatrix out of comp |
startLoc | location to update |
startTime | time after which update is needed |
void SpatialSEIR::ModelContext::calculateGenericCompartment_OCL | ( | int * | comp, |
int * | comp0, | ||
int * | compStarAdd, | ||
int * | compStarSub | ||
) |
calculateGenericCompartment_OCL works identically to calculateGenericCompartment_CPU while using calls to oclProvider to perform calculations in parallel.
comp | pointer to CompartmentalModelMatrix to be calculated |
comp0 | pointer to integer array containing starting values for comp |
compStarAdd | pointer to transition matrix into comp |
compStarSub | pointer to transition CompartmentalModelMatrix out of comp |
void SpatialSEIR::ModelContext::calculateI_CPU | ( | ) |
The zero parameter overload of calculateI_CPU calculates the I compartment from A0, I_star, and R_star.
void SpatialSEIR::ModelContext::calculateI_CPU | ( | int | startLoc, |
int | startTime | ||
) |
The two parameter overload of calculateI_CPU again calculates the I compartment from A0, I_star and R_star, but assumes that only values of I_star and R_star at column (location) startLoc and after row (time) startTime have changed since I was last calculated. When applicable, this saves considerable time.
void SpatialSEIR::ModelContext::calculateI_givenR_CPU | ( | ) |
The zero parameter overload of calculate I_givenR_CPU takes advantage of the fact that N=S+E+I+R to update I when R changes.
void SpatialSEIR::ModelContext::calculateI_givenR_CPU | ( | int | startLoc, |
int | startTime | ||
) |
The two parameter overload of calculate I_givenR_CPU takes advantage of the fact that N=S+E+I+R to update I when R changes, while assuming that R has only changed for location startLoc and after time point startTime.
void SpatialSEIR::ModelContext::calculateI_OCL | ( | ) |
calculateI_OCL works identically to calculateI_CPU using calls to oclProvider to perform computations in parallel.
double * SpatialSEIR::ModelContext::calculateIntegratedG | ( | int | t | ) |
Calculates the next generation matrix, integrated forward in time.
void SpatialSEIR::ModelContext::calculateP_EI_CPU | ( | ) |
calculateP_EI_CPU calculates the E to I transition probability vector p_ei
void SpatialSEIR::ModelContext::calculateP_IR_CPU | ( | ) |
calculateP_IR_CPU calculates the I to R transition probability vector p_ir
void SpatialSEIR::ModelContext::calculateP_RS_CPU | ( | ) |
calculateP_RS_CPU calculates the reinfection probability p_rs.
void SpatialSEIR::ModelContext::calculateP_SE_CPU | ( | ) |
The zero parameter overload of calculateP_SE_CPU calculates the exposure probability p_se.
void SpatialSEIR::ModelContext::calculateP_SE_CPU | ( | int | startLoc, |
int | startTime | ||
) |
The two parameter overload of calculateP_SE_CPU calculates the exposure probability starting at startLoc and startTime. cacheP_SE_Calculation must be called first.
void SpatialSEIR::ModelContext::calculateP_SE_OCL | ( | ) |
calculateP_SE_OCL works identically to calculateP_SE_CPU, while making calls to oclProvider to perform calculations in parallel.
double * SpatialSEIR::ModelContext::calculateR0Components | ( | int | t | ) |
Calculates the location specific components of the effective R0 quantity at time t
void SpatialSEIR::ModelContext::calculateR_CPU | ( | ) |
The zero parameter overload of calculateR_CPU calculates the R compartment from A0, R_star, and S_star.
void SpatialSEIR::ModelContext::calculateR_CPU | ( | int | startLoc, |
int | startTime | ||
) |
The two parameter overload of calculateR_CPU again calculates the R compartment from A0, R_star and S_star, but assumes that only values of R_star and S_star at column (location) startLoc and after row (time) startTime have changed since I was last calculated. When applicable, this saves considerable time.
void SpatialSEIR::ModelContext::calculateR_givenS_CPU | ( | ) |
The zero parameter overload of calculate R_givenS_CPU takes advantage of the fact that N=S+E+I+R to update R when S changes.
void SpatialSEIR::ModelContext::calculateR_givenS_CPU | ( | int | startLoc, |
int | startTime | ||
) |
The two parameter overload of calculate R_givenS_CPU takes advantage of the fact that N=S+E+I+R to update R when S changes, while assuming that S has only changed for location startLoc and after time point startTime.
void SpatialSEIR::ModelContext::calculateR_OCL | ( | ) |
calculateR_OCL works identically to calculateR_CPU using calls to oclProvider to perform computations in parallel. The four parameter overload of calculateGenericCompartment_CPU takes advantege of the common calculation pattern for S,E,I, and R to provide a single function mapping starting and transition values to the filan CompartmentalModelMatrix values.
void SpatialSEIR::ModelContext::calculateS_CPU | ( | ) |
The zero parameter overload of calculateS_CPU calculates the S compartment from A0, S_star, and E_star.
void SpatialSEIR::ModelContext::calculateS_CPU | ( | int | startLoc, |
int | startTime | ||
) |
The two parameter overload of calculateS_CPU again calculates the S compartment from A0, S_star and E_star, but assumes that only values of S_star and E_star at column (location) startLoc and after row (time) startTime have changed since S was last calculated. When applicable, this saves considerable time.
void SpatialSEIR::ModelContext::calculateS_givenE_CPU | ( | ) |
The zero parameter overload of calculate S_givenE_CPU takes advantage of the fact that N=S+E+I+R to update S when E changes.
void SpatialSEIR::ModelContext::calculateS_givenE_CPU | ( | int | startLoc, |
int | startTime | ||
) |
The two parameter overload of calculate S_givenE_CPU takes advantage of the fact that N=S+E+I+R to update S when E changes, while assuming that E has only changed for location startLoc and after time point startTime.
void SpatialSEIR::ModelContext::calculateS_OCL | ( | ) |
calculateS_OCL works identically to calculateS_CPU using calls to oclProvider to perform computations in parallel.
int SpatialSEIR::ModelContext::checkCompartmentBounds | ( | ) |
checkCompartmentBounds is a debug mode function which checks for impossible compartment values and prints errors accordingly.
void SpatialSEIR::ModelContext::configureIterationTasks | ( | ) |
configureIterationTasks determines whether there are any conditions which require additional logical tasks to be performed during MCMC sampling.
double SpatialSEIR::ModelContext::estimateEffectiveR0 | ( | ) |
Estimates the effective reproductive rate across all time periods and spatial locations
double * SpatialSEIR::ModelContext::estimateEffectiveR0 | ( | int | t | ) |
Estimates the effective reproductive rate at time t.
double SpatialSEIR::ModelContext::estimateR0 | ( | ) |
Estimates the average basic reproductive number across all time periods and spatial locations.
double * SpatialSEIR::ModelContext::estimateR0 | ( | int | t | ) |
Estimates the basic reproductive number for time point t
int SpatialSEIR::ModelContext::getCompartmentSamplingMode | ( | ) |
getCompartmentSamplingMode returns the current sampling mode as an integer
int SpatialSEIR::ModelContext::getParameterSamplingMode | ( | ) |
getParameterSamplingMode returns the current sampling mode for non-compartment parameters as an integer
void SpatialSEIR::ModelContext::populate | ( | InitData * | _A0, |
covariateArgs * | xArgs, | ||
covariateArgs * | xPrsArgs, | ||
double * | offset, | ||
int * | Y, | ||
compartmentArgs * | S_starArgs, | ||
compartmentArgs * | E_starArgs, | ||
compartmentArgs * | I_starArgs, | ||
compartmentArgs * | R_starArgs, | ||
scaledDistanceArgs * | scaledDistArgs, | ||
double * | rho, | ||
double * | phi, | ||
double * | beta, | ||
double * | gamma_ei, | ||
double * | gamma_ir, | ||
double * | betaPrs, | ||
int * | N, | ||
sliceParameters * | sliceWidths, | ||
priorControl * | priorInformation, | ||
modelConfiguration | _config | ||
) |
The populate method is the entry point to working with spatial SEIR models in libspatialSEIR.
_A0 | A0 must be an instance of InitData, which contains the starting values of S0, E0, I0, R0 |
xArgs | xArgs is an instance of covariateArgs which contains the dimensions of the fixed and time varying covariate matrices driving the exposure process. |
xPrsArgs | xPrsArgs is an instance of covariateArgs which contains the dimensions of the fixed and time varying covariate matrices driving the reinfection process. |
offset | Offset allows irregular spacing between time points. |
Y | Y is the actual observed data vector. It must be of the same dimension as the rest of the Compartments. In the default data model, Y is exactly equal to I_star. |
S_starArgs | S_starArgs is an instance of compartmentArgs containing the dimensions and data for S_star, along with the steadyStateConstraintPrecision parameter. |
E_starArgs | E_starArgs is an instance of compartmentArgs containing the dimensions and data for E_star, along with the steadyStateConstraintPrecision parameter. |
I_starArgs | I_starArgs is an instance of compartmentArgs containing the dimensions and data for I_star, along with the steadyStateConstraintPrecision parameter. |
R_starArgs | R_starArgs is an instance of compartmentArgs containing the dimensions and data for R_star, along with the steadyStateConstraintPrecision parameter. |
scaledDistArgs | scaledDistArgs is the scaledDistanceArgs struct containing the data and dimension of the unscaled distance matrices. |
rho | rho gives the starting value of the spatial dependence parameter |
phi | starting value for overdispersion precision |
beta | beta gives the starting vector of regression parameters driving the exposure process. |
gamma_ei | p_ei is the starting value of the parameter driving the exposed to infectious transition probability. |
gamma_ir | p_ir is the starting value of the parameter driving the infectious to removed/recovered transition probability. |
betaPrs | betaPrs is the starting vector of regression parameters driving the reinfection process. |
N | N is the matrix of population sizes, corresponding in dimension to the CompartmentalModelMatrix instances. |
sliceWidths | sliceWidths is an instance of sliceParameters which gives the slice sampling widths/ Metropolis-Hastings tuning parameters for the various FullConditional distributions. |
priorInformation | The priorControl struct gives the prior precisions for the regression parameters, and the prior alpha and beta parameters for the p_ei and p_ir distributions. |
_config | The modelConfiguration struct gives the reinfection mode (SEIRS, SEIR, SSEIR) and the MCMC sampling mode. MCMC sampling modes are under very active development. |
void SpatialSEIR::ModelContext::printFCValues | ( | ) |
printFCValues is a semi-depricated debug mode function which displays the current full conditional likelihood values of the model parameters.
void SpatialSEIR::ModelContext::runSimulation | ( | int | nIterations, |
bool | verbose = false , |
||
bool | debug = false |
||
) |
runSimulation runs nIterations MCMC updates given the current model state. Optionally, verbose and debug level output is available.
void SpatialSEIR::ModelContext::setCompartmentSamplingMode | ( | int | mode | ) |
setCompartmentSamplingMode sets the sampling mode for the disease compartments. This part of the API is in flux.
void SpatialSEIR::ModelContext::setParameterSamplingMode | ( | int | mode | ) |
setParameterSamplingMode sets the sampling mode for non compartment parameters. This part of the API is in flux.
void SpatialSEIR::ModelContext::setRandomSeed | ( | unsigned int | seedValue | ) |
setRandomSeed sets the seed used by the pseudorandom number generator provided by RandomNumberProvider
void SpatialSEIR::ModelContext::simulationIter | ( | bool | verbose = false , |
bool | debug = false |
||
) |
simulationIter runs a single MCMC update given the current model state. Optionally, verbose and debug level output is available.
int SpatialSEIR::ModelContext::totalE | ( | ) |
Calculates the sum of all components of E.
int SpatialSEIR::ModelContext::totalE | ( | int | tpt | ) |
Calculates the sum of all components of E.
int SpatialSEIR::ModelContext::totalE_star | ( | ) |
Calculates the sum of all components of E_star.
int SpatialSEIR::ModelContext::totalE_star | ( | int | tpt | ) |
Calculates the sum of all components of E_star.
int SpatialSEIR::ModelContext::totalI | ( | ) |
Calculates the sum of all components of I.
int SpatialSEIR::ModelContext::totalI | ( | int | tpt | ) |
Calculates the sum of all components of I.
int SpatialSEIR::ModelContext::totalI_star | ( | ) |
Calculates the sum of all components of I_star.
int SpatialSEIR::ModelContext::totalI_star | ( | int | tpt | ) |
Calculates the sum of all components of I_star.
int SpatialSEIR::ModelContext::totalR | ( | ) |
Calculates the sum of all components of R.
int SpatialSEIR::ModelContext::totalR | ( | int | tpt | ) |
Calculates the sum of all components of R.
int SpatialSEIR::ModelContext::totalR_star | ( | ) |
Calculates the sum of all components of R_star.
int SpatialSEIR::ModelContext::totalR_star | ( | int | tpt | ) |
Calculates the sum of all components of R_star.
int SpatialSEIR::ModelContext::totalS | ( | ) |
Calculates the sum of all components of S.
int SpatialSEIR::ModelContext::totalS | ( | int | tpt | ) |
Calculates the sum of all components of S.
int SpatialSEIR::ModelContext::totalS_star | ( | ) |
Calculates the sum of all components of S_star.
int SpatialSEIR::ModelContext::totalS_star | ( | int | tpt | ) |
Calculates the sum of all components of S_star.
void SpatialSEIR::ModelContext::updateSamplingParameters | ( | double | desiredRatio, |
double | targetWidth, | ||
double | proportionChange | ||
) |
updateSamplingParameters moves the various slice sampling widths / Metropolis-Hastings tuning parameters up or down depending on the current acceptance rate, and re-sets the acceptance counters.
desiredRatio | desired MCMC acceptance ratio |
targetWidth | target ratio tolerance (how close is close enough) |
proportionChange | proportion to change the sampling parameters by |
InitData* SpatialSEIR::ModelContext::A0 |
Pointer to InitData instance containing the data for the initial compartment sizes.
double* SpatialSEIR::ModelContext::beta |
Storage for the regression parameters beta
FC_Beta* SpatialSEIR::ModelContext::beta_fc |
Pointer to FullConditional distribution for the regression parameters beta
double* SpatialSEIR::ModelContext::betaPrs |
Storage for the regression parameters betaP_RS
FC_Beta_P_RS* SpatialSEIR::ModelContext::betaPrs_fc |
Pointer to FullConditional distribution for the regression parameters betaP_RS
double* SpatialSEIR::ModelContext::compartmentCache |
ComparmentalModelMatrix sized cache of doubles.
modelConfiguration* SpatialSEIR::ModelContext::config |
Pointer to modelConfiguration instance.
PerformDecorrelationStep* SpatialSEIR::ModelContext::decorrelationStepTask |
Pointer to the decorrelation step task
CompartmentalModelMatrix* SpatialSEIR::ModelContext::E |
Pointer to the CompartmentalModelMatrix data structure storing the E compartment.
FC_E0* SpatialSEIR::ModelContext::E0_fc |
Pointer to FullConditional distribution for the initial exposed population - not used
CompartmentalModelMatrix* SpatialSEIR::ModelContext::E_star |
Pointer to the CompartmentalModelMatrix data structure storing the E_star compartment.
FC_E_Star* SpatialSEIR::ModelContext::E_star_fc |
Pointer to FullConditional distribution for the transition matrix E_star
double* SpatialSEIR::ModelContext::eta |
Storage for the linear predictor, eta
IOProvider* SpatialSEIR::ModelContext::fileProvider |
pointer to an IOProvider instance which writes selected sampler output to a comma delimited text file.
double* SpatialSEIR::ModelContext::gamma |
Storage for the external infection parameters gamma (DEPRICATED)
double* SpatialSEIR::ModelContext::gamma_ei |
Parameter controlling transition probability vector p_ei
FC_Gamma_EI* SpatialSEIR::ModelContext::gamma_ei_fc |
Pointer to FullConditional distribution for the parameter transition probability p_ei
double* SpatialSEIR::ModelContext::gamma_ir |
Parameter controlling transition probability vector p_ir
FC_Gamma_IR* SpatialSEIR::ModelContext::gamma_ir_fc |
Pointer to FullConditional distribution for the parameter controlling transition probability p_ir
CompartmentalModelMatrix* SpatialSEIR::ModelContext::I |
Pointer to the CompartmentalModelMatrix data structure storing the I compartment.
FC_I0* SpatialSEIR::ModelContext::I0_fc |
Pointer to FullConditional distribution for the initial infectious population
CompartmentalModelMatrix* SpatialSEIR::ModelContext::I_star |
Pointer to the CompartmentalModelMatrix data structure storing the I_star compartment.
FC_I_Star_overdispersed* SpatialSEIR::ModelContext::I_star_overdispersed_fc |
Pointer to FullConditional distribution for the I_star transition matrix under overdispersion
int* SpatialSEIR::ModelContext::indexLength |
Length of the indexList array.
int* SpatialSEIR::ModelContext::indexList |
List of indices, used by some sampling functions.
int* SpatialSEIR::ModelContext::isPopulated |
Indicator that ModelContext instance is ready for use.
std::vector<IterationTask*>* SpatialSEIR::ModelContext::iterationTasks |
Iteration tasks are bits of logical code that need to execute periodically during MCMC sampling, but don't fall into the usual FullConditional framework.
std::vector<FullConditional*>* SpatialSEIR::ModelContext::model |
A model in libspatialSEIR is a collection of full conditional distributions which must be sampled from. This collection is assembled during ModelContext.populate based on the reinfection mode, data model (not yet implemented), and the detection of special cases such as the single location case.
int* SpatialSEIR::ModelContext::N |
Matrix of population sizes
int* SpatialSEIR::ModelContext::numIterations |
Number of MCMC samples so far.
OCLProvider* SpatialSEIR::ModelContext::oclProvider |
pointer to an instance of OCLProvider which supplies an interface for selecting active OpenCL devices as well as parallelized computations required for various FullConditional distributions.
double* SpatialSEIR::ModelContext::offset |
Storage for the vector of offsets.
double* SpatialSEIR::ModelContext::p_ei |
Transition probability vector p_ei
double* SpatialSEIR::ModelContext::p_ir |
Transition probability vector p_ir
double* SpatialSEIR::ModelContext::p_rs |
Transition probability vector p_rs
double* SpatialSEIR::ModelContext::p_se |
Storage for the exposure probability matrix, p_se
double* SpatialSEIR::ModelContext::p_se_components |
Storage for cacheable portions of the p_se calculation
PerformHybridSE_EI_UpdateStep* SpatialSEIR::ModelContext::performHybridSE_EI_UpdateTask |
Pointer to the beta+gamma_ei hybrid sampling task.
double* SpatialSEIR::ModelContext::phi |
Storage for the overdispersion parameter
FC_Phi* SpatialSEIR::ModelContext::phi_fc |
Pointer to FullConditional distribution for the overdispersion parameter for I_star
CompartmentalModelMatrix* SpatialSEIR::ModelContext::R |
Pointer to the CompartmentalModelMatrix data structure storing the R compartment.
FC_R0* SpatialSEIR::ModelContext::R0_fc |
Pointer to FullConditional distribution for the initial recovered population - not used
CompartmentalModelMatrix* SpatialSEIR::ModelContext::R_star |
Pointer to the CompartmentalModelMatrix data structure storing the R_star compartment.
FC_R_Star* SpatialSEIR::ModelContext::R_star_fc |
Pointer to FullConditional distribution for the transition matrix R_star
RandomNumberProvider* SpatialSEIR::ModelContext::random |
pointer to a RandomNumberProvider instance which provides random numbers from various distributions as well implementations of several probability density functions.
double* SpatialSEIR::ModelContext::rho |
Storage for the vector of spatial dependence parameters, rho
FC_Rho* SpatialSEIR::ModelContext::rho_fc |
Pointer to FullConditional distribution for the spatial depenence parameter rho
CompartmentalModelMatrix* SpatialSEIR::ModelContext::S |
Pointer to the CompartmentalModelMatrix data structure storing the S compartment.
FC_S0* SpatialSEIR::ModelContext::S0_fc |
Pointer to FullConditional distribution for the initial suscptible population.
CompartmentalModelMatrix* SpatialSEIR::ModelContext::S_star |
Pointer to the CompartmentalModelMatrix data structure storing the S_star compartment.
FC_S_Star* SpatialSEIR::ModelContext::S_star_fc |
Pointer to FullConditional distribution for the transition matrix S_star
std::vector<DistanceMatrix*>* SpatialSEIR::ModelContext::scaledDistMatrices |
Pointer to vector of scaled DistanceMatrix objects.
SetCompartmentSamplingIndicesTask* SpatialSEIR::ModelContext::setSamplingIndicesTask |
Pointer to the sampling indices task.
int* SpatialSEIR::ModelContext::singleLocation |
Indicator that ModelContext is operating in the single location special case.
CompartmentalModelMatrix* SpatialSEIR::ModelContext::tmpContainer |
Extra compartment storage for caching integer computations
CovariateMatrix* SpatialSEIR::ModelContext::X |
Pointer to CovariateMatrix instance containing covariate informating driving the exposure process.
CovariateMatrix* SpatialSEIR::ModelContext::X_pRS |
Pointer to CovariateMatrix instance containing covariate informating driving the reinfection process.
int* SpatialSEIR::ModelContext::Y |
Pointer to the actual data vector: Y. Y is related by a data model to I_star.