default.itemps | R Documentation |
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
default.itemps(m = 40, type = c("geometric", "harmonic","sigmoidal"), k.min = 0.1, c0n0 = c(100, 1000), lambda = c("opt", "naive", "st"))
m |
Number of temperatures in the ladder; |
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 ( |
lambda |
Method for combining the importance samplers at each
temperature. Optimal combination ( lambda[k] = sum(w[k,]))^2/sum(w[k,]^2). Setting |
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),
and the harmonic ladder by
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")
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 |
k |
The generated inverse temperature ladder; a vector
with |
pk |
A vector with |
lambda |
IT method, as specified by the input argument |
Robert B. Gramacy, rbg@vt.edu, and Matt Taddy, mataddy@amazon.com
Gramacy, R.B., Samworth, R.J., and King, R. (2010) Importance Tempering. ArXiV article 0707.4242 Statistics and Computing, 20(1), pp. 1-7; https://arxiv.org/abs/0707.4242.
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.
For the geometric temperature ladder:
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)
https://bobby.gramacy.com/r_packages/tgp/
btgp
## 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)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.