MC2Pack Function Reference

Chapter contents

Tables

Table 1: Sampling methods

Table 2: Candidate densities

Table 3: Marginal likelihood computation methods

MC2Pack function members

GetDraws

MC2Pack::GetDraws(const amTheta, const avW)

amTheta

out: Pointer to matrix of size k x nTheta with sampled parameter values.

avW

out: if !0 on input, 1 x nTheta vector with weights, for samplers using importance sampling, scaled such that the weights sum to 1. For non-importance samplers, equal weights are returned.

Return value

1 if a sample was recently run, 0 otherwise.

Return the draws of the last run of the sampler

GetLocationScale

MC2Pack::GetLocationScale()

Return value

Returns two arrays avMu and amS2 with both the first, initial, and the updated versions of locations and scale matrices of intermediate rotations.

Extracts the location and scale matrices of successive rotations.

GetPackageName

MC2Pack::GetPackageName()

Return value

String with name of package, in this case "MC2Pack"

GetPackageVersion

MC2Pack::GetPackageVersion()

Return value

Integer, version number of "MC2Pack" times 100

GetParNames

MC2Pack::GetParNames()

Return value

asNames: array of strings, parameter names

Return the names of the parameters being sampled

GetRNE

MC2Pack::GetRNE()

Return value

mRNE: iK x nTap matrix with relative numerical efficiencies, or 0 when no sample is available

Return the relative numerical efficiency according to [Kim et al. (1998)], with

^ R   Bm  = 1 +  2Bm Bm-1 Bm å i=1  K(i/Bm) ^ r   (i)

with K(j) the Parzen kernel, B_m the bandwidth and [^(r)](i) the ith order autocorrelation. The bandwidth, i.e. the number of drawings used in the computations, is set with SetTaper.

GetTaper

MC2Pack::GetTaper()

Return value

vTaper: 1 x nTap vector with tapering perunages

Return the values used in tapering, in computing the relative numerical efficiencies

SetCandidate

MC2Pack::SetCandidate(const iCand, const xCand);

iCand

in: Either MC_CNORM, MC_CSTUD, MC_CRW or MC_CUSER, see also table Candidate densities

xCand

in: Meaning depends on cForm. If cForm is

MC_CNORM

double dFact

MC_CSTUD

vector with dFact and degrees of freedom of Student-t density

MC_CRW

double dFact

MC_CUSER

function of format

            fnRndCand(const iR, const vTheta, const vMu, const mCS2,
                      const amCand, const avLnDens)
            

with inputs

iR

in: integer, number of vectors to be drawn

vTheta

in: iK x 1 vector of last draw

vMu

in: iK x 1 vector of location parameter

mCS2

in: iK x iK matrix with choleski decomposition of scale matrix

and outputs

amCand

out: iK x iR matrix of random vectors drawn from candidate density

avLnDens

out: 1 x iR vector with logarithm of density in amCand.

No return value
Select the functional form of the candidate density to be used. The value dFact specifies the multiplication factor for the choleski decomposition of the last covariance matrix. For MC_CNORM, MC_CSTUD a value slightly larger than 1 would be logical, for MC_CRW a small value e.g. 0.05 is better.

Note that the matrix mCS2 may be singular, when only a subset of parameters should be sampled as indicated through the SetConditional(vInd, -1) function.

SetConditional

MC2Pack::SetConditional(const vInd, const fnRanCondGG)
MC2Pack::SetConditional(const vInd, const fnRanCondGG, const fnDensCondGG)

vInd

in: iI vector with indices into the parameter vector, indicating which parameters are sampled from this conditional density.

fnRanCondGG

in: function which samples from the conditional density for elements vInd. Function should be of the format

ir= fnRanCondGG(const vInd, const avP)

The function should replace the elements avP[0][vInd] with newly sampled elements.

If, instead of a function, a number (e.g. -1) is passed along to SetConditional, the parameter is skipped in the simulation.

fnLnDensCondGG

in: (optional) function which computes the logarithm of the conditional density for elements vInd, conditional on the other elements. The function should be of the format

dLnDens= fnLnDensCondGG(const vInd, const vP)

The function should compute the logarithm of the conditional density of the elements vP[vInd], conditional on the other elements. Note that the conditional density should contain all integrating constants.

The function is only used for computing the marginal likelihood using the Gibbs method, see MargLik(MC_MLGG, ...).

Return value

1 if everything went correctly, 0 in case of problems

Indicate the sampling function from the analytic conditional density for the (Griddy) Gibbs sampler.

SetDebug

MC2Pack::SetDebug(const iDebug)

iDebug

in: integer, level of debug information.

No return value

Set the level of debugging information to be generated. Default is 0, no extra information.

SetEps

MC2Pack::SetEps(const dEpsAbs, const dEpsRel)

dEpsAbs

in: double, absolute precision used by QuadPack in MC2Pack, default 0

dEpsRel

in: double, relative precision used by QuadPack in MC2Pack, default 1e-2

No return value

Set the level of precision during the integration routines.

SetGraphs

MC2Pack::SetGraphs(const bShowGraphs)

bShowGraphs

in: boolean, indicating if graphs should be shown (default=TRUE)

No return value

Indicate if intermediate graphics should be shown on the screen

SetInfo

MC2Pack::SetInfo(const iInfoRep)

iInfoRep

in: integer, number of iterations between printing info. Default= 1000 (accepted) drawings

No return value

Change the setting of the number of iterations between printing a report on the expected time until termination.

SetIntegration

MC2Pack::SetIntegration(const iAdaptive, ...)
MC2Pack::SetIntegration(const iAdaptive, const iBounds, const bSMax, 
                         const bSeparate)

iAdaptive

in: integer, indicating the method of integration

0: non-adaptive QNG (default)

1: adaptive QAG

2: iterative Simpson

iBounds

in: integer, indicating change of bounds, with

0: No change

1: Change both bounds (default)

2: Use integration from -¥ to +¥

bSMax

in: boolean, search maximum? (default= TRUE)

bSeparate

in: boolean, indicating if integration should be separated (fixing p(h|r = 0) at 0) or not. Default is true, only used for APS/APIS sampling.

No return value

Specify the method of integration.

SetLikelihood

MC2Pack::SetLikelihood(const fnLikelihood)

fnLikelihood

in: function, of same format as in MaxBFGS

No return value

Specify the function which returns the average loglikelihood, where the averaging is done with respect to the number of observations specified in SetLocationScale.

SetLimits

MC2Pack::SetLimits(const mLUBounds)
MC2Pack::SetLimits(const mLUBounds, const mBC)
MC2Pack::SetLimits(const mLUBounds, const mBC, const bChangeLU)

mLUBounds

in: iK x 2 matrix with bounds on the parameters.

mBC

in: iR x iK+1 matrix, consisting of iR x iK matrix mB with a iR x 1 vector mC added, indicating the linear restraints mB*vTheta <= mC

bChange

in: Boolean, if TRUE the bounds mLUBounds are adapted between rotations to the average between the (empirical bounds +- 10 bounds, in such a way that the original bounds are never enlarged. Default= FALSE.

No return value

Set the limits for the sampling using APS/APIS/Griddy Gibbs

SetLocationScale

MC2Pack::SetLocationScale(const avMu, const amS2, const iN, 
                           const bOptimize, ...)
MC2Pack::SetLocationScale(const avMu, const amS2, const iN, 
                           const bOptimize, const asNames)

avMu

in: (column) vector, initial value of mean

amS2

in: matrix, initial value of covariance

iN

in: Integer, number of observations contained in the posterior.

bOptim

in: boolean, indicating if initial values should be used as a starting point for an optimization to find the posterior mode.

asNames

in: (optional) array of strings with names of parameters

avMu

out: vector, initial value of mean used for sampling (unchanged if bOptim was false)

amS2

out: matrix, initial value of covariance used for sampling (unchanged if bOptim was false)

Return value.

Value according to MaxBFGS, or MAX_CONV if bOptim was false

Provide a starting point for the sampler

SetMahalanobisFraction

MC2Pack::SetMahalanobisFraction(const dFrac)

dFrac

in: double, the fraction of improvement in the Mahalanobis distance that warrants repeating a rotation, default= 0.5

No return value

SetMethod

MC2Pack::SetMethod(const iMethod)

iMethod

in: Method, one of MC_MH, MC_IS, MC_APS, MC_APIS, MC_GG, see table MCMC sampling Methods

No return value

Choose the method of simulation

SetOutput

MC2Pack::SetOutput(const sResbase)
MC2Pack::SetOutput(const sResbase, const bRewrite)

sResbase

in: string, base-name of output files

bRewrite

in: boolean, optional argument, indicates if output file is rewritten (default is FALSE)

No return value

Indicate the base-name of the output files. By default, output is written to files with name mcmc<method>. If bRewrite=TRUE (default), textual output file is rewritten.

SetParNames

MC2Pack::SetParNames(const asNames)

asNames

in: array of strings, parameter names

No return value

Set the names of the parameters being sampled

SetPosterior

MC2Pack::SetPosterior(const fnPosterior)

fnPosterior

in: function, of same format as in MaxBFGS

No return value

Specify the function which returns the average log-posterior, where the averaging is done with respect to the number of observations specified in SetLocationScale.

SetPrior

MC2Pack::SetPrior(const fnPrior)

fnPrior

in: function, of same format as in MaxBFGS

No return value

Specify the function which returns the average log-prior, where the averaging is done with respect to the number of observations specified in SetLocationScale.

SetRanPrior

MC2Pack::SetRanPrior(const fnRanPrior)

fnRanPrior

in: function, of format as in

      mU= fnRanPrior(iR);

Return value

iK x iR matrix with iR draws from the prior.

Specify the function which returns iR draws from the prior

SetSample

MC2Pack::SetSample(const vRep, ...)

vRep

in: vector of size iRot, or scalar, number of repetitions of drawing sampled parameter values in each of the rotations

vBurn

in: vector of size iRot, or scalar, number of drawings which are disregarded for a burn-in period

vSkip

in: vector of size iRot, or scalar, number of drawings from which one drawing is saved (ie, if vSkip = 1, all drawings are saved)

vSerRej

in: vector of size iRot, or scalar, number of iterated rejections of the Metropolis-Hastings step, before an acceptance is forced.

vDir

in: vector of size iRot, or scalar, number of directions which are sampled in the case of APS/APIS sampling

No return value

Set the size of the sample

SetTaper

MC2Pack::SetTaper(const vTaper)

vTaper

1 x nTap vector with tapering perunages, or number of parameter vectors to use. Default= <0.05, 0.15>

No return value

Set the values used in tapering, in computing the relative numerical efficiencies

CuSum

MC2Pack::Cusum(const bAdapt)

bAdapt

in: Boolean, indicating if Cusum plots should be adapted to the variance of the sample

No return value

Prepares and shows a Cusum plot according to [Yu and Mykland (1998)].

CumPlot

MC2Pack::CumPlot(const viCum)

viCum

in: integer, or vector of integers, with length of past to account for in cumulative mean plots. Default (if NaN is passed along) is 100.

No return value

Prepares and shows a cumulative mean plot.

Marglik

Marglik(const iMethod, ...)

The actual inputs depend on the method which is chosen for computing the marginal likelihood. Possible inputs are:

iMethod

in: One of MC_MLKERN, MC_MLLP, MC_MLHM, MC_MLGG, ML_MLINT, MC_MLPR, MC_MLIS (see also table 3), indicating the method to be used.

vMu

in: iK x iP matrix with iP location vectors for computation. If not provided, the mean of the sampled drawings is used.

mS2

in: k x k covariance matrix. If not provided (or M_NAN) the numerical approximation to the covariance is used.

nRep

in: Number of repetitions to be used in the Gibbs chain, in the prior sampling, or in the importance sampling. If not provided, the number of repetitions in the last iteration is used. Do not put this number too high, as computing the marginal likelihood using the Gibbs method may be a lengthy operation.

nInt

in: Number of intervals using in numerically integrating the posterior density, when using the method MC_MLNUM. Default= 20. Put this number very low, as the numerical integral in iK dimension may take a tremendously long time.

dFrac

in: Double, fraction of sample to use in Harmonic Mean computations, or number of parameter vectors (if > 1). Default= 1.

The following table gives the arguments for each of the computation methods, and functions which have to be set before calling Marglik:

Method	Arguments	Functions used
MC_MLKERN	vMu	Posterior
MC_MLLP	vMu, mS2	Posterior
MC_MLHM	dFrac	Likelihood
MC_MLGG	vMu, nRep	Posterior, Conditionals (if available)
MC_MLINT	nInt	Posterior
MC_MLPR	nRep	Likelihood, RanPrior
MC_MLIS	vMu, mS2, nRep	Prior, Likelihood

Return value

1 x iP vector with logarithm of marginal likelihoods at locations vMu.
Compute the marginal likelihood of a model, using either a kernel approximation, a LaPlace approximation, the harmonic mean method, the Gibbs approximation, brute force numerical integration or sampling from the prior density. See [Bos (2002)] for an overview and further references.

Simulate

MC2Pack::Simulate()

Return value

Array with components vMu, mS2, the mean and covariance of the final sample.

Initiates the sampling routines. Output is in the text file as indicated by SetOutput with extension .out, and in a data file with extension .fmt. The data file is of size nDraws x k or nDraws x k+1, with the extra column containing the log-weights of the importance sampling algorithms.

Bibliography

[Bauwens et al. (2004)]

Bauwens, L., C. S. Bos, H. K. Van Dijk, and R. D. Van Oest (2004). Adaptive radial-based direction sampling: Some flexible and robust Monte Carlo integration methods. Journal of Econometrics , #25. Forthcoming.

[Bos (2002)]

Bos, C. S. (2002). A comparison of marginal likelihood computation methods. In W. Härdle and B. Ronz (Eds.), COMPSTAT 2002 - Proceedings of Computational Statistics, pp. 111-117.

[Chib and Greenberg (1995)]

Chib, S. and E. Greenberg (1995). Understanding the Metropolis-Hastings algorithm. The American Statistician  49 (4), 327-335.

[Kim et al. (1998)]

Kim, S., N. Shephard, and S. Chib (1998). Stochastic volatility: Likelihood inference and comparison with ARCH models. Review of Economic Studies  64, 361-393.

[Kloek and Van Dijk (1978)]

Kloek, T. and H. K. Van Dijk (1978). Bayesian estimates of equation system parameters: An application of integration by Monte Carlo. Econometrica  46, 1-20.

[Ritter and Tanner (1992)]

Ritter, C. and M. A. Tanner (1992, September). Facilitating the Gibbs sampler: The Gibbs stopper and the Griddy-Gibbs sampler. Journal of the American Statistical Association   87 (419), 861-868.

[Yu and Mykland (1998)]

Yu, B. and P. Mykland (1998). Looking at Markov samplers through cusum path plots: A simple diagnostic idea. Statistics and Computing  8 (3), 275-286.