default.itemps: Default Sigmoidal, Harmonic and Geometric Temperature Ladders In tgp: Bayesian Treed Gaussian Process Models

Description

Parameterized by the minimum desired inverse temperature, this function generates a ladder of inverse temperatures `k[1:m]` starting at `k[1] = 1`, with `m` steps down to the final temperature `k[m] = k.min` progressing sigmoidally, harmonically or geometrically. The output is in a format convenient for the `b*` functions in the tgp package (e.g. `btgp`), including stochastic approximation parameters c0 and n0 for tuning the uniform pseudo-prior output by this function

Usage

 ```1 2 3``` ```default.itemps(m = 40, type = c("geometric", "harmonic","sigmoidal"), k.min = 0.1, c0n0 = c(100, 1000), lambda = c("opt", "naive", "st")) ```

Arguments

 `m` Number of temperatures in the ladder; `m=1` corresponds to importance sampling at the temperature specified by `k.min` (in this case all other arguments are ignored) `type` Choose from amongst two common defaults for simulated tempering and Metropolis-coupled MCMC, i.e., geometric (default) or harmonic, or a sigmoidal ladder (default) that concentrates more inverse temperatures near 1 `k.min` Minimum inverse temperature desired `c0n0` Stochastic approximation parameters used to tune the simulated tempering pseudo-prior (`\$pk`) to get a uniform posterior over the inverse temperatures; must be a 2-vector of positive integers `c(c0, n0)`; see the Geyer \& Thompson reference below `lambda` Method for combining the importance samplers at each temperature. Optimal combination (`"opt"`) is the default, weighting the IS at each temperature k by lambda[k] = sum(w[k,]))^2/sum(w[k,]^2). Setting `lambda = "naive"` allows each temperature to contribute equally (λ[k] = 1, or equivalently ignores delineations due to temperature when using importance weights. Setting `lambda = "st"` allows only the first (cold) temperature to contribute to the estimator, thereby implementing simulated tempering

Details

The geometric and harmonic inverse temperature ladders are usually defined by an index i = 1:m and a parameter delta > 0. The geometric ladder is defined by

k[i] = (1 + delta)^(1-i),

k[i] = (1 + delta*(i-1))^(-1).

Alternatively, specifying the minimum temperature k.min in the ladder can be used to uniquely determine delta. E.g., for the geometric ladder

delta = k.min^(1/(1-m))-1,

and for the harmonic

delta = (k.min^(-1)-1)/(m-1).

In a similar spirit, the sigmoidal ladder is specified by first situating m indices j[i] in Re so that k[1] = k(j[1]) = 1 and k[m] = k(j[m]) = k.min under

k(j[i]) = 1.01 - 1/(1+exp(-j[i])).

The remaining j[2:(m-1)] are spaced evenly between j[i] and j[m] to fill out the ladder k[2:(m-1)] = k(j[2:(m-1)]).

For more details, see the Importance tempering paper cited below and a full demonstration in `vignette("tgp2")`

Value

The return value is a `list` which is compatible with the input argument `itemps` to the `b*` functions (e.g. `btgp`), containing the following entries:

 `c0n0 ` A copy of the `c0n0` input argument `k ` The generated inverse temperature ladder; a vector with `length(k) = m` containing a decreasing sequence from `1` down to `k.min` `pk ` A vector with `length(pk) = m` containing an initial pseudo-prior for the temperature ladder of `1/m` for each inverse temperature `lambda` IT method, as specified by the input argument

Author(s)

Robert B. Gramacy, [email protected], and Matt Taddy, [email protected]

References

Gramacy, R.B., Samworth, R.J., and King, R. (2007) Importance Tempering. ArXiV article 0707.4242 http://arxiv.org/abs/0707.4242. To appear in Statistics and Computing.

For stochastic approximation and simulated tempering (ST):

Geyer, C.~and Thompson, E.~(1995). Annealing Markov chain Monte Carlo with applications to ancestral inference. Journal of the American Statistical Association, 90, 909–920.

Neal, R.M.~(2001) Annealed importance sampling. Statistics and Computing, 11, 125–129

Justifying geometric and harmonic defaults:

Liu, J.S.~(1002) Monte Carlo Strategies in Scientific Computing. New York: Springer. Chapter 10 (pages 213 \& 233)

`btgp`
 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38``` ```## comparing the different ladders geo <- default.itemps(type="geometric") har <- default.itemps(type="harmonic") sig <- default.itemps(type="sigmoidal") par(mfrow=c(2,1)) matplot(cbind(geo\$k, har\$k, sig\$k), pch=21:23, main="inv-temp ladders", xlab="indx", ylab="itemp") legend("topright", pch=21:23, c("geometric","harmonic","sigmoidal"), col=1:3) matplot(log(cbind(sig\$k, geo\$k, har\$k)), pch=21:23, main="log(inv-temp) ladders", xlab="indx", ylab="itemp") ## Not run: ## using Importance Tempering (IT) to improve mixing ## on the motorcycle accident dataset library(MASS) out.it <- btgpllm(X=mcycle[,1], Z=mcycle[,2], BTE=c(2000,22000,2), R=3, itemps=default.itemps(), bprior="b0", trace=TRUE, pred.n=FALSE) ## compare to regular tgp w/o IT out.reg <- btgpllm(X=mcycle[,1], Z=mcycle[,2], BTE=c(2000,22000,2), R=3, bprior="b0", trace=TRUE, pred.n=FALSE) ## compare the heights explored by the three chains: ## REG, combining all temperatures, and IT p <- out.it\$trace\$post L <- length(p\$height) hw <- suppressWarnings(sample(p\$height, L, prob=p\$wlambda, replace=TRUE)) b <- hist2bar(cbind(out.reg\$trace\$post\$height, p\$height, hw)) par(mfrow=c(1,1)) barplot(b, beside=TRUE, xlab="tree height", ylab="counts", col=1:3, main="tree heights encountered") legend("topright", c("reg MCMC", "All Temps", "IT"), fill=1:3) ## End(Not run) ```