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 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77
|
temper <- function(obj, initial, neighbors, nbatch, blen = 1,
nspac = 1, scale = 1, outfun, debug = FALSE, parallel = FALSE, ...)
UseMethod("temper")
temper.tempering <- function(obj, initial, neighbors, nbatch, blen = 1,
nspac = 1, scale = 1, outfun, debug = FALSE, parallel = FALSE, ...)
{
# like metrop, ignore initial argument
initial <- obj$final
# makes no (at least little?) sense to change neighbor structure
neighbors <- obj$neighbors
if (missing(nbatch)) nbatch <- obj$nbatch
if (missing(blen)) blen <- obj$blen
if (missing(nspac)) nspac <- obj$nspac
if (missing(scale)) scale <- obj$scale
if (missing(debug)) debug <- obj$debug
# makes no sense to change from parallel to serial or vice versa
# size and shape of state wouldn't even be the same
parallel <- obj$parallel
assign(".Random.seed", obj$final.seed, .GlobalEnv)
if (missing(outfun)) {
if (is.null(obj$outfun)) {
temper.function(obj$lud, initial, neighbors, nbatch, blen,
nspac, scale, debug = debug, parallel = parallel, ...)
} else {
temper.function(obj$lud, initial, neighbors, nbatch, blen,
nspac, scale, obj$outfun, debug = debug, parallel = parallel,
...)
}
} else {
temper.function(obj$lud, initial, neighbors, nbatch, blen,
nspac, scale, outfun, debug = debug, parallel = parallel, ...)
}
}
temper.function <- function(obj, initial, neighbors, nbatch, blen = 1,
nspac = 1, scale = 1, outfun, debug = FALSE, parallel = FALSE, ...)
{
if (! exists(".Random.seed")) runif(1)
saveseed <- .Random.seed
func1 <- function(state) obj(state, ...)
env1 <- environment(fun = func1)
if (missing(outfun)) {
func2 <- NULL
env2 <- NULL
outfun <- NULL
} else if (is.function(outfun)) {
func2 <- function(state) outfun(state, ...)
env2 <- environment(fun = func2)
}
stopifnot(is.numeric(initial))
storage.mode(initial) <- "double"
if (is.list(scale)) {
for (i in 1:length(scale)) {
stopifnot(is.numeric(scale[[i]]))
storage.mode(scale[[i]]) <- "double"
}
} else {
stopifnot(is.numeric(scale))
storage.mode(scale) <- "double"
}
out.time <- system.time(
out <- .Call(C_temper, func1, initial, neighbors, nbatch, blen, nspac,
scale, func2, debug, parallel, env1, env2)
)
result <- structure(c(list(lud = obj,
neighbors = neighbors, nbatch = nbatch, blen = blen, nspac = nspac,
scale = scale, outfun = outfun, debug = debug, parallel = parallel,
initial.seed = saveseed, final.seed = .Random.seed, time = out.time),
out), class = c("mcmc", "tempering"))
return(result)
}
|