findMaxTemper: Find the maximum temperature for parallel MCMC chains

In EMC: Evolutionary Monte Carlo (EMC) algorithm

Description

Multiple MCMC chains based algorithms (e.g., parallel tempering, evolutionary Monte Carlo) need a temperature ladder. This function finds the maximum temperature for constructing the ladder.

Below sampDim refers to the dimension of the sample space, temperLadderLen refers to the length of the temperature ladder, and levelsSaveSampForLen refers to the length of levelsSaveSampFor. Note, this function calls evolMonteCarlo, so some of the arguments below have the same name and meaning as the corresponding ones for evolMonteCarlo. See details below for explanation on the arguments.

Usage

findMaxTemper(nIters,
              statsFuncList,
              startingVals,
              logTarDensFunc,
              MHPropNewFunc,
              logMHPropDensFunc = NULL, 
              temperLadder      = NULL,
              temperLimits      = NULL,
              ladderLen         = 10,
              scheme            = 'exponential',
              schemeParam       = 0.5,
              cutoffDStats      = 1.96,
              cutoffESS         = 50,
              guideMe           = TRUE,
              levelsSaveSampFor = NULL,
              saveFitness       = FALSE,
              doFullAnal        = TRUE,
              verboseLevel      = 0,
              ...)   

Arguments

nIters	integer > 0.
statsFuncList	list of functions of one argument each, which return the value of the statistic evaluated at one MCMC sample or draw.
startingVals	double matrix of dimension temperLadderLen x sampDim or vector of length sampDim, in which case the same starting values are used for every temperature level.
logTarDensFunc	function of two arguments (draw, ...) that returns the target density evaluated in the log scale.
MHPropNewFunc	function of four arguments (temperature, block, currentDraw, ...) that returns new Metropolis-Hastings proposals. See details below on the argument block.
logMHPropDensFunc	function of five arguments (temperature, block, currentDraw, proposalDraw, ...) that returns the proposal density evaluated in the log scale. See details below on the argument block.
temperLadder	double vector with all positive entries, in decreasing order.
temperLimits	double vector with two positive entries.
ladderLen	integer > 0.
scheme	character.
schemeParam	double > 0.
cutoffDStats	double > 0.
cutoffESS	double > 0.
guideMe	logical.
levelsSaveSampFor	integer vector with positive entries.
saveFitness	logical.
doFullAnal	logical.
verboseLevel	integer, a value >= 2 produces a lot of output.
...	optional arguments to be passed to logTarDensFunc, MHPropNewFunc and logMHPropDensFunc.

Details

This function is based on the method to find the temperature range introduced in section 4.1 of Goswami and Liu (2007).

statsFuncList

The user specifies this list of functions, each of which is known to be sensitive to the presence of modes. For example, if both dimension 1 and 3 are sensitive to presence of modes, then one could use:

    coord1        <- function (xx) { xx[1] }
    

    coord3        <- function (xx) { xx[3] }
    

    statsFuncList <- list(coord1, coord3)
    

temperLadder

This is the temperature ladder needed for the first stage preliminary run. One can either specify a temperature ladder via temperLadder or specify temperLimits, ladderLen, scheme and schemeParam. For details on the later set of parameters, see below. Note, temperLadder overrides temperLimits, ladderLen, scheme and schemeParam.

temperLimits

temperLimits = c(lowerLimit, upperLimit) is a two-tuple of positive numbers, where the lowerLimit is usually 1 and upperLimit is a number in [100, 1000]. If stochastic optimization (via sampling) is the goal, then lowerLimit is taken to be in [0, 1].

ladderLen, scheme and schemeParam

These three parameters are required (along with temperLimits) if temperLadder is not provided. We recommend taking ladderLen in [15, 30]. The allowed choices for scheme and schemeParam are:

scheme	schemeParam
========	=============
linear	NA
log	NA
geometric	NA
mult-power	NA
add-power	>= 0
reciprocal	NA
exponential	>= 0
tangent	>= 0

We recommended using scheme = 'exponential' and schemeParam in [0.3, 0.5].

cutoffDStats

This cutoff comes from Normal_1(0, 1), the standard normal distribution (Goswami and Liu, 2007); the default value 1.96 is a conservative cutoff. Note if you have more than one statistic in statsFuncList, which is usually the case, using this cutoff may result in different suggested maximum temperatures (as can be seen by calling the print function on the result of findMaxTemper). A conservative recommendation is that you choose the maximum of the suggested temperatures as the final maximum temperature for use in placeTempers and later in parallelTempering or evolMonteCarlo.

cutoffESS

a cutoff for the effective sample size (ESS) of the underlying Markov chain ergodic estimator and the importance sampling estimators.

guideMe

If guideMe = TRUE, then the function suggests different modifications to alter the setting towards a re-run, in case there are problems with the underlying MCMC run.

doFullAnal

If doFullAnal = TRUE, then the search for the maximum temperature is conducted among all the levels of the temperLadder. In case this switch is turned off, the search for maximum temperature is done in a greedy (and faster) manner, namely, search is stopped as soon as all the statistic(s) in the statsFuncList find some maximum temperature(s). Note, the greedy search may result in much higher maximum temperature (and hence sub-optimal) than needed, so it is not recommended.

levelsSaveSampFor

This is passed to evolMonteCarlo for the underlying MCMC run.

Value

This function returns a list with the following components:

temperLadder	the temperature ladder used for the underlying MCMC run.
DStats	the D-statistic (Goswami and Liu, 2007) values used to find the maximum temperature.
cutoffDStats	the cutoffDStats argument.
nIters	the post burn-in nIters.
levelsSaveSampFor	the levelsSaveSampFor argument.
draws	array of dimension nIters x sampDim x levelsSaveSampForLen, if saveFitness = FALSE. If saveFitness = TRUE, then the returned array is of dimension nIters x (sampDim + 1) x levelsSaveSampForLen; i.e., each of the levelsSaveSampForLen matrices contain the fitness values in their last column.
startingVals	the startingVals argument.
intermediate statistics	a bunch of intermediate statistics used in the computation of DStats, namely, MCEsts, MCVarEsts, MCESS, ISEsts, ISVarEsts, ISESS, each being computed for all the statistics provided by statsFuncList argument.
time	the time taken by the run.

Note

The effect of leaving the default value NULL for some of the arguments above are as follows:

logMHPropDensFunc	the proposal density MHPropNewFunc is deemed symmetric.
temperLadder	valid temperLimits, ladderLen, scheme and schemeParam
	are provided, which are used to construct the temperLadder.
temperLimits	a valid temperLadder is provided.
levelsSaveSampFor	temperLadderLen.

Author(s)

Gopi Goswami [e-mail]

References

Gopi Goswami and Jun S. Liu (2007). On learning strategies for evolutionary Monte Carlo. Statistics and Computing 17:1:23-38.

See Also

placeTempers, parallelTempering, evolMonteCarlo

Examples

## Not run: 
coord1        <- function (xx) { xx[1] }
ss            <- function (xx) { sum(xx) }
pp            <- function (xx) { prod(xx) }
statsFuncList <- list(coord1, ss, pp)
maxTemperObj  <-
    with(VShapedFuncGenerator(-13579),
         findMaxTemper(nIters            = 15000,
                       statsFuncList     = statsFuncList,
                       temperLadder      = c(20, 15, 10, 5, 1),
                       startingVals      = c(0, 0),
                       logTarDensFunc    = logTarDensFunc,
                       MHPropNewFunc     = MHPropNewFunc,
                       levelsSaveSampFor = seq_len(5),
                       doFullAnal        = TRUE,
                       verboseLevel      = 1))
print(maxTemperObj)
print(names(maxTemperObj))
with(maxTemperObj,
 {
     par(mfcol = c(3, 3))
     for (ii in seq_along(levelsSaveSampFor)) {
         main <- paste('temper:', round(temperLadder[levelsSaveSampFor[ii]], 3))
         plot(draws[ , , ii],
              xlim = c(-10, 25),
              ylim = c(-10, 10),
              pch  = '.',
              ask  = FALSE,
              main = as.expression(main),
              xlab = as.expression(substitute(x[xii], list(xii = 1))),
              ylab = as.expression(substitute(x[xii], list(xii = 2))))   
     }
 })

## End(Not run)

EMC documentation built on May 2, 2019, 5:11 a.m.