ctsem (version 3.2.1)

stanoptimis: Optimize / importance sample a stan or ctStan model.

Description

Optimize / importance sample a stan or ctStan model.

Usage

stanoptimis(
  standata,
  sm,
  init = "random",
  initsd = 0.01,
  sampleinit = NA,
  deoptim = FALSE,
  estonly = FALSE,
  tol = 1e-12,
  decontrol = list(),
  stochastic = "auto",
  nopriors = FALSE,
  carefulfit = TRUE,
  subsamplesize = 0.5,
  finitediff = FALSE,
  datadrivenpars = c(),
  plot = FALSE,
  is = FALSE,
  isloopsize = 1000,
  finishsamples = 1000,
  tdf = 10,
  chancethreshold = 100,
  finishmultiply = 5,
  verbose = 0,
  cores = 2
)

Arguments

standata

list object conforming to rstan data standards.

sm

compiled stan model object.

init

vector of unconstrained parameter values, or character string 'random' to initialise with random values very close to zero.

initsd

positive numeric specifying sd of normal distribution governing random sample of init parameters, if init='random' .

sampleinit

either NA, or an niterations * nparams matrix of samples to initialise importance sampling.

deoptim

Do first pass optimization using differential evolution? Slower, but better for cases with multiple minima / difficult optimization.

estonly

if TRUE,just return point estimates under $rawest subobject.

tol

objective tolerance.

decontrol

List of control parameters for differential evolution step, to pass to DEoptim.control.

stochastic

Logical. Use stochastic gradient descent instead of mize (bfgs) optimizer. Still experimental, worth trying for either robustness checks or problematic, high dimensional, nonlinear, problems.

nopriors

logical. If TRUE, a nopriors integer is set to 1 (TRUE) in the standata object -- only has an effect if the stan model uses this value.

carefulfit

Logical. If TRUE, priors are always used for a rough first pass to obtain starting values when nopriors=TRUE.

subsamplesize

value between 0 and 1 representing proportion of subjects to include in first pass fit.

finitediff

Either 'ask', TRUE, or FALSE. Whether to use the slow finite difference calculations for the Hessian (used for confidence intervals) if other approaches do not give a positive definite result.

datadrivenpars

vector of integers denoting which parameters should begin fixed at zero, and only freed if expected to substantially contribute to model fit.

plot

Logical. If TRUE, plot iteration details. Probably slower.

is

Logical. Use importance sampling, or just return map estimates?

isloopsize

Number of samples of approximating distribution per iteration of importance sampling.

finishsamples

Number of samples to draw (either from hessian based covariance or posterior distribution) for final results computation.

tdf

degrees of freedom of multivariate t distribution. Higher (more normal) generally gives more efficent importance sampling, at risk of truncating tails.

chancethreshold

drop iterations of importance sampling where any samples are chancethreshold times more likely to be drawn than expected.

finishmultiply

Importance sampling stops once available samples reach finishsamples * finishmultiply , then the final samples are drawn without replacement from this set.

verbose

Integer from 0 to 2. Higher values print more information during model fit -- for debugging.

cores

Number of cpu cores to use, should be at least 2.

Value

list containing fit elements

Examples

# NOT RUN {
library(rstan)
scode <- "
parameters {
  real y[2];
}
model {
  y[1] ~ normal(0, 1);
  y[2] ~ double_exponential(0, 2);
}
"

sm <- stan_model(model_code=scode)
fit <- sampling(sm, iter = 10000)
summary(fit)$summary

## extract samples as a list of arrays
e <- ctExtract(fit, permuted = TRUE)

#for ml or map estimates
optimis <- stanoptimis(standata = list(),sm = sm,finishsamples = 3000,cores=2)
optimis$optimfit

#for posterior distributions
optimis <- stanoptimis(standata = list(),sm = sm,finishsamples = 3000,cores=2,tdf=5)

apply(optimis$rawposterior,2,mean)
apply(optimis$rawposterior,2,sd)
isdiag(optimis)

plot(density(optimis$rawposterior[,2],bw=.05))
points(density(e$y[,2],bw=.05),type='l',col=2)
# }