plot.PMCMC | R Documentation |
PMCMC
objectsPlot method for PMCMC
objects.
## S3 method for class 'PMCMC' plot( x, type = c("post", "trace"), joint = FALSE, transfunc = NA, ask = TRUE, ... )
x |
A |
type |
Takes the value |
joint |
A logical describing whether joint or marginal distributions are wanted. |
transfunc |
Is a |
ask |
Should the user ask before moving onto next trace plot. |
... |
Not used here. |
A plot of the (approximate) posterior distributions obtained from fitting a particle Markov chain Monte Carlo algorithm, or provides corresponding trace plots.
PMCMC
, print.PMCMC
, predict.PMCMC
, summary.PMCMC
window.PMCMC
## set up data to pass to PMCMC flu_dat <- data.frame( t = 1:14, Robs = c(3, 8, 26, 76, 225, 298, 258, 233, 189, 128, 68, 29, 14, 4) ) ## set up observation process obs <- data.frame( dataNames = "Robs", dist = "pois", p1 = "R + 1e-5", p2 = NA, stringsAsFactors = FALSE ) ## set up model (no need to specify tspan ## argument as it is set in PMCMC()) transitions <- c( "S -> beta * S * I / (S + I + R + R1) -> I", "I -> gamma * I -> R", "R -> gamma1 * R -> R1" ) compartments <- c("S", "I", "R", "R1") pars <- c("beta", "gamma", "gamma1") model <- mparseRcpp( transitions = transitions, compartments = compartments, pars = pars, obsProcess = obs ) ## set priors priors <- data.frame( parnames = c("beta", "gamma", "gamma1"), dist = rep("unif", 3), stringsAsFactors = FALSE) priors$p1 <- c(0, 0, 0) priors$p2 <- c(5, 5, 5) ## define initial states iniStates <- c(S = 762, I = 1, R = 0, R1 = 0) set.seed(50) ## run PMCMC algorithm post <- PMCMC( x = flu_dat, priors = priors, func = model, u = iniStates, npart = 25, niter = 5000, nprintsum = 1000 ) ## plot MCMC traces plot(post, "trace") ## continue for some more iterations post <- PMCMC(post, niter = 5000, nprintsum = 1000) ## plot traces and posteriors plot(post, "trace") plot(post) ## remove burn-in post <- window(post, start = 5000) ## summarise posteriors summary(post)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.