From f5fa65856b52d8ead1a222d15e69a245b7277b7c Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Mon, 23 Apr 2018 21:28:20 -0700 Subject: [PATCH 01/11] Install working --- install_StanHeaders.R | 14 ++++++++------ rstan/R_Makevars | 1 + rstan/rstan/R/git_rstan_head.R | 1 + 3 files changed, 10 insertions(+), 6 deletions(-) diff --git a/install_StanHeaders.R b/install_StanHeaders.R index 17120a10b..89ac6948c 100644 --- a/install_StanHeaders.R +++ b/install_StanHeaders.R @@ -3,11 +3,13 @@ # including its submodules. This requires the git2r and devtools packages path_rstan <- tempfile(pattern = "git2r-") -git2r::clone("http://github.com/stan-dev/rstan", path_rstan, branch = "develop") -git2r::clone("http://github.com/stan-dev/stan", - file.path(path_rstan, "StanHeaders", "inst", "include", "upstream"), - branch = "develop") # may want to change this branch -git2r::clone("http://github.com/stan-dev/math", - file.path(path_rstan, "StanHeaders", "inst", "include", "mathlib"), +path_stan_dev <- "/home/rgiordan/Documents/git_repos/stan-dev/" +git2r::clone(file.path(path_stan_dev, "rstan"), + path_rstan, branch = "add_hessian4") +git2r::clone(file.path(path_stan_dev, "stan"), + file.path(path_rstan, "StanHeaders", "inst", "include", "upstream"), + branch = "add_hessians") # may want to change this branch +git2r::clone(file.path(path_stan_dev, "math"), + file.path(path_rstan, "StanHeaders", "inst", "include", "mathlib"), branch = "develop") # may want to change this branch devtools::install(file.path(path_rstan, "StanHeaders"), args = "--preclean") diff --git a/rstan/R_Makevars b/rstan/R_Makevars index 889fafb5f..a16a3b1e0 100644 --- a/rstan/R_Makevars +++ b/rstan/R_Makevars @@ -6,3 +6,4 @@ MAKEFLAGS=-j6 # CXXFLAGS = -O3 $(LTO) +CXXFLAGS+=-std=c++14 \ No newline at end of file diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index eb486edd1..5a5fbc7c1 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -29,3 +29,4 @@ git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" +git_head <- function() "50a10cca0de66ee51cf3832d1aa6ef29b3590659" From 284d6c7f4bcce37633c4c0a5cf0ba8c1a837a075 Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Mon, 23 Apr 2018 21:35:49 -0700 Subject: [PATCH 02/11] Expose hessian as stanfit method --- rstan/rstan/NAMESPACE | 30 +- rstan/rstan/R/git_rstan_head.R | 3 + rstan/rstan/R/stanfit-class.R | 721 +++++++++--------- .../rstan/rcpp_module_def_for_rstan.hpp | 2 + rstan/rstan/inst/include/rstan/stan_fit.hpp | 42 +- 5 files changed, 427 insertions(+), 371 deletions(-) diff --git a/rstan/rstan/NAMESPACE b/rstan/rstan/NAMESPACE index f629be884..c3c2bb1d7 100644 --- a/rstan/rstan/NAMESPACE +++ b/rstan/rstan/NAMESPACE @@ -20,7 +20,7 @@ export( stan_model, stanc, stanc_builder, - stan_version, + stan_version, stan, stan_rdump, read_rdump, @@ -34,12 +34,12 @@ export( rstan_options, As.mcmc.list, set_cppo, - stan_plot, - stan_trace, - stan_hist, - stan_dens, + stan_plot, + stan_trace, + stan_hist, + stan_dens, stan_scat, - stan_ac, + stan_ac, stan_diag, stan_par, stan_rhat, @@ -53,8 +53,8 @@ export( cpp_object_initializer, # get_rstan.options check_hmc_diagnostics, - check_divergences, - check_treedepth, + check_divergences, + check_treedepth, check_energy, get_rng, get_stream, @@ -62,21 +62,21 @@ export( ) exportClasses( - stanmodel, stanfit -) + stanmodel, stanfit +) exportMethods( -# print, plot, -# extract, +# print, plot, +# extract, optimizing, vb, - get_cppcode, get_cxxflags, # for stanmodel + get_cppcode, get_cxxflags, # for stanmodel show, sampling, summary, extract, traceplot, plot, get_stancode, get_inits, get_seed, get_cppo_mode, - log_prob, grad_log_prob, + log_prob, grad_log_prob, hessian_log_prob, unconstrain_pars, constrain_pars, get_num_upars, get_seeds, get_adaptation_info, get_sampler_params, - get_logposterior, + get_logposterior, get_posterior_mean, get_elapsed_time, get_stanmodel diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index 5a5fbc7c1..4c4106380 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -30,3 +30,6 @@ git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" git_head <- function() "50a10cca0de66ee51cf3832d1aa6ef29b3590659" +git_head <- function() "57444a10c088dbc52105fb394aadb996cd0c9df0" +git_head <- function() "2165604262c9085f0e47408d60bb04c76e5485c3" +git_head <- function() "2165604262c9085f0e47408d60bb04c76e5485c3" diff --git a/rstan/rstan/R/stanfit-class.R b/rstan/rstan/R/stanfit-class.R index 8b5625132..60d7ba39e 100644 --- a/rstan/rstan/R/stanfit-class.R +++ b/rstan/rstan/R/stanfit-class.R @@ -15,75 +15,75 @@ # along with this program; if not, write to the Free Software # Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. -setMethod("show", "stanfit", +setMethod("show", "stanfit", function(object) { print.stanfit(x = object, pars = object@sim$pars_oi) - }) + }) -print.stanfit <- function(x, pars = x@sim$pars_oi, - probs = c(0.025, 0.25, 0.5, 0.75, 0.975), - digits_summary = 2, include = TRUE, ...) { - if (x@mode == 1L) { +print.stanfit <- function(x, pars = x@sim$pars_oi, + probs = c(0.025, 0.25, 0.5, 0.75, 0.975), + digits_summary = 2, include = TRUE, ...) { + if (x@mode == 1L) { cat("Stan model '", x@model_name, "' is of mode 'test_grad';\n", "sampling is not conducted.\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } else if (x@mode == 2L) { - cat("Stan model '", x@model_name, "' does not contain samples.\n", sep = '') - return(invisible(NULL)) - } + cat("Stan model '", x@model_name, "' does not contain samples.\n", sep = '') + return(invisible(NULL)) + } if(!include) pars <- setdiff(x@sim$pars_oi, pars) - s <- summary(x, pars, probs, ...) + s <- summary(x, pars, probs, ...) if (is.null(s)) return(invisible(NULL)) n_kept <- x@sim$n_save - x@sim$warmup2 cat("Inference for Stan model: ", x@model_name, '.\n', sep = '') - cat(x@sim$chains, " chains, each with iter=", x@sim$iter, - "; warmup=", x@sim$warmup, "; thin=", x@sim$thin, "; \n", - "post-warmup draws per chain=", n_kept[1], ", ", + cat(x@sim$chains, " chains, each with iter=", x@sim$iter, + "; warmup=", x@sim$warmup, "; thin=", x@sim$thin, "; \n", + "post-warmup draws per chain=", n_kept[1], ", ", "total post-warmup draws=", sum(n_kept), ".\n\n", sep = '') - if (!is.null(x@stan_args[[1]]$method) && + if (!is.null(x@stan_args[[1]]$method) && x@stan_args[[1]]$method == "variational") { - print(round(s$summary, digits_summary), ...) - cat("\nApproximate samples were drawn using VB(", x@stan_args[[1]]$algorithm, ") at ", x@date, + print(round(s$summary, digits_summary), ...) + cat("\nApproximate samples were drawn using VB(", x@stan_args[[1]]$algorithm, ") at ", x@date, ".\n", sep = '') message("We recommend genuine 'sampling' from the posterior distribution for final inferences!") return(invisible(NULL)) } - + # round n_eff to integers s$summary[, 'n_eff'] <- round(s$summary[, 'n_eff'], 0) - print(round(s$summary, digits_summary), ...) + print(round(s$summary, digits_summary), ...) sampler <- attr(x@sim$samples[[1]], "args")$sampler_t cat("\nSamples were drawn using ", sampler, " at ", x@date, ".\n", - "For each parameter, n_eff is a crude measure of effective sample size,\n", + "For each parameter, n_eff is a crude measure of effective sample size,\n", "and Rhat is the potential scale reduction factor on split chains (at \n", "convergence, Rhat=1).\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } -setMethod("plot", signature(x = "stanfit", y = "missing"), +setMethod("plot", signature(x = "stanfit", y = "missing"), function(x, ..., plotfun) { if (x@mode == 1L) { cat("Stan model '", x@model_name, "' is of mode 'test_grad';\n", "sampling is not conducted.\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } else if (x@mode == 2L) { - cat("Stan model '", x@model_name, - "' does not contain samples.\n", sep = '') - return(invisible(NULL)) - } + cat("Stan model '", x@model_name, + "' does not contain samples.\n", sep = '') + return(invisible(NULL)) + } if (isTRUE(all.equal(x@sim$n_save, x@sim$warmup2, check.attributes = FALSE, check.names = FALSE))) { - cat("Stan model '", x@model_name, + cat("Stan model '", x@model_name, "' does not contain samples after warmup.\n", sep = '') return(invisible(NULL)) } - - if (missing(plotfun)) fun <- "stan_plot" + + if (missing(plotfun)) fun <- "stan_plot" else { plotters <- paste0("stan_", c("plot", "trace", "scat", "hist", "dens", "ac", "diag", "rhat", "ess", "mcse", "par")) @@ -92,46 +92,46 @@ setMethod("plot", signature(x = "stanfit", y = "missing"), else fun <- match.fun(fun) } do.call(fun, list(object = x, ...)) - }) + }) setGeneric(name = "get_stancode", - def = function(object, ...) { standardGeneric("get_stancode")}) + def = function(object, ...) { standardGeneric("get_stancode")}) -setGeneric(name = "get_cppo_mode", - def = function(object, ...) { standardGeneric("get_cppo_mode") }) +setGeneric(name = "get_cppo_mode", + def = function(object, ...) { standardGeneric("get_cppo_mode") }) -setMethod('get_cppo_mode', signature = "stanfit", - function(object) { +setMethod('get_cppo_mode', signature = "stanfit", + function(object) { cxxf <- get_cxxflags(object@stanmodel) if (identical(cxxf, character(0))) return(NA) l <- get_cxxo_level(cxxf) - if ("" == l) l <- "0" - p <- match(l, c("3", "2", "1", "0")) + if ("" == l) l <- "0" + p <- match(l, c("3", "2", "1", "0")) c("fast", "presentation2", "presentation1", "debug")[p] - }) + }) -setMethod('get_stancode', signature = "stanfit", +setMethod('get_stancode', signature = "stanfit", function(object, print = FALSE) { code <- object@stanmodel@model_code - if (print) cat(code, "\n") + if (print) cat(code, "\n") return(code) - }) + }) -setGeneric(name = 'get_stanmodel', +setGeneric(name = 'get_stanmodel', def = function(object, ...) { standardGeneric("get_stanmodel")}) -setMethod("get_stanmodel", signature = "stanfit", - function(object) { +setMethod("get_stanmodel", signature = "stanfit", + function(object) { return(object@stanmodel) - }) + }) -setGeneric(name = 'get_inits', +setGeneric(name = 'get_inits', def = function(object, ...) { standardGeneric("get_inits")}) -setMethod("get_inits", signature = "stanfit", - function(object, iter = NULL) { +setMethod("get_inits", signature = "stanfit", + function(object, iter = NULL) { if (is.null(iter)) return(object@inits) - stopifnot(is.numeric(iter), iter > 0, + stopifnot(is.numeric(iter), iter > 0, iter <= length(object@sim$samples[[1]][[1]])) inits <- object@inits if (length(inits) == 0) { @@ -146,58 +146,58 @@ setMethod("get_inits", signature = "stanfit", vec <- as.matrix(as.data.frame(object@sim$samples[[c]]))[iter,] inits[[c]] <- relist(vec, skeleton = inits[[c]]) } - return(inits) + return(inits) }) -setGeneric(name = 'get_seed', +setGeneric(name = 'get_seed', def = function(object, ...) { standardGeneric("get_seed")}) -setMethod("get_seed", signature = "stanfit", - function(object) { - if (length(object@stan_args) < 1L) return(NULL) +setMethod("get_seed", signature = "stanfit", + function(object) { + if (length(object@stan_args) < 1L) return(NULL) object@stan_args[[1]]$seed }) -setGeneric(name = 'get_seeds', +setGeneric(name = 'get_seeds', def = function(object, ...) { standardGeneric("get_seeds")}) -setMethod("get_seeds", signature = "stanfit", - function(object) { - if (length(object@stan_args) < 1L) return(NULL) +setMethod("get_seeds", signature = "stanfit", + function(object) { + if (length(object@stan_args) < 1L) return(NULL) sapply(object@stan_args, function(x) x$seed) }) ### HELPER FUNCTIONS -### +### get_kept_samples <- function(n, sim) { # # Args: - # sim: the sim slot in object stanfit - # n: the nth parameter (starting from 1) - # Note: - # samples from different chains are merged. + # sim: the sim slot in object stanfit + # n: the nth parameter (starting from 1) + # Note: + # samples from different chains are merged. - # get chain kept samples (gcks) + # get chain kept samples (gcks) gcks <- function(s, nw, permutation) { - s <- s[[n]][-(1:nw)] - s[permutation] - } + s <- s[[n]][-(1:nw)] + s[permutation] + } ss <- mapply(gcks, sim$samples, sim$warmup2, sim$permutation, - SIMPLIFY = FALSE, USE.NAMES = FALSE) - do.call(c, ss) -} + SIMPLIFY = FALSE, USE.NAMES = FALSE) + do.call(c, ss) +} get_kept_samples2 <- function(n, sim) { - # a different implementation of get_kept_samples - # It seems this one is faster than get_kept_samples - # TODO: to understand why it is faster? + # a different implementation of get_kept_samples + # It seems this one is faster than get_kept_samples + # TODO: to understand why it is faster? lst <- vector("list", sim$chains) - for (ic in 1:sim$chains) { - if (sim$warmup2[ic] > 0) + for (ic in 1:sim$chains) { + if (sim$warmup2[ic] > 0) lst[[ic]] <- sim$samples[[ic]][[n]][-(1:sim$warmup2[ic])][sim$permutation[[ic]]] - else + else lst[[ic]] <- sim$samples[[ic]][[n]][sim$permutation[[ic]]] - } + } do.call(c, lst) } @@ -205,100 +205,100 @@ get_samples <- function(n, sim, inc_warmup = TRUE) { # get chain samples # Args: # n: parameter index (integer) - # sim: the sim list in stanfit - # + # sim: the sim list in stanfit + # # Returns: - # a list of chains for the nth parameter; each chain is an - # element of the list. + # a list of chains for the nth parameter; each chain is an + # element of the list. if (all(sim$warmup2 == 0)) inc_warmup <- TRUE # for the case warmup sample is discarded gcs <- function(s, inc_warmup, nw) { # Args: - # s: samples of all chains - # nw: number of warmup + # s: samples of all chains + # nw: number of warmup if (inc_warmup) return(s[[n]]) - else return(s[[n]][-(1:nw)]) - } - ss <- mapply(gcs, sim$samples, inc_warmup, sim$warmup2, - SIMPLIFY = FALSE, USE.NAMES = FALSE) - ss -} + else return(s[[n]][-(1:nw)]) + } + ss <- mapply(gcs, sim$samples, inc_warmup, sim$warmup2, + SIMPLIFY = FALSE, USE.NAMES = FALSE) + ss +} get_samples2 <- function(n, sim, inc_warmup = TRUE) { - # serves the same purpose with get_samples, but with - # different implementation - # It seems that this one is fast. + # serves the same purpose with get_samples, but with + # different implementation + # It seems that this one is fast. if (all(sim$warmup2 == 0)) inc_warmup <- TRUE # for the case warmup sample is discarded npar <- length(n) lst <- vector("list", sim$chains) for (ic in 1:sim$chains) { - lst[[ic]] <- + lst[[ic]] <- if (inc_warmup) sim$samples[[ic]][[n]] else sim$samples[[ic]][[n]][-(1:sim$warmup2[ic])] } lst -} +} par_traceplot <- function(sim, n, par_name, inc_warmup = TRUE, window = NULL, ...) { # same thin, n_save, warmup2 for all the chains thin <- sim$thin - warmup2 <- sim$warmup2[1] - warmup <- sim$warmup - main <- paste("Trace of ", par_name) + warmup2 <- sim$warmup2[1] + warmup <- sim$warmup + main <- paste("Trace of ", par_name) chain_cols <- rstan_options("rstan_chain_cols") - warmup_col <- rstan_options("rstan_warmup_bg_col") + warmup_col <- rstan_options("rstan_warmup_bg_col") - start_i <- window[1] + start_i <- window[1] window_size <- (window[2] - start_i) %/% thin id <- seq.int(start_i, by = thin, length.out = window_size) start_idx <- (if (warmup2 == 0) (start_i - warmup) else start_i) %/% thin if (start_idx < 1) start_idx <- 1 idx <- seq.int(start_idx, by = 1, length.out = window_size) - yrange <- NULL + yrange <- NULL for (i in 1:sim$chains) - yrange <- range(yrange, sim$samples[[i]][[n]][idx]) + yrange <- range(yrange, sim$samples[[i]][[n]][idx]) if (inc_warmup) { plot(range(id), yrange, type = 'n', bty = 'l', xlab = 'Iterations', ylab = "", main = main, ...) - rect(par("usr")[1], par("usr")[3], sim$warmup, par("usr")[4], + rect(par("usr")[1], par("usr")[3], sim$warmup, par("usr")[4], col = warmup_col, border = NA) - } else { + } else { plot(range(id), yrange, type = 'n', bty = 'l', xlab = 'Iterations (without warmup)', ylab = "", main = main, ...) - } + } - for (i in 1:sim$chains) - lines(id, sim$samples[[i]][[n]][idx], - xlab = '', ylab = '', col = chain_cols[(i-1) %% 6 + 1], ...) -} + for (i in 1:sim$chains) + lines(id, sim$samples[[i]][[n]][idx], + xlab = '', ylab = '', col = chain_cols[(i-1) %% 6 + 1], ...) +} ###### -setGeneric(name = 'get_adaptation_info', +setGeneric(name = 'get_adaptation_info', def = function(object, ...) { standardGeneric("get_adaptation_info")}) -setMethod("get_adaptation_info", +setMethod("get_adaptation_info", definition = function(object) { if (object@mode == 1L) { cat("Stan model '", object@model_name, "' is of mode 'test_grad';\n", "sampling is not conducted.\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } else if (object@mode == 2L) { - cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') - return(invisible(NULL)) - } + cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') + return(invisible(NULL)) + } lai <- lapply(object@sim$samples, function(x) attr(x, "adaptation_info")) - is_empty <- function(x) { - if (is.null(x)) return(TRUE) + is_empty <- function(x) { + if (is.null(x)) return(TRUE) if (is.character(x) && all(nchar(x) == 0)) return(TRUE) FALSE } - if (all(sapply(lai, FUN = is_empty))) return(invisible(NULL)) - return(lai) - }) + if (all(sapply(lai, FUN = is_empty))) return(invisible(NULL)) + return(lai) + }) -setGeneric(name = "get_logposterior", +setGeneric(name = "get_logposterior", def = function(object, ...) { standardGeneric("get_logposterior")}) @@ -307,53 +307,53 @@ setMethod("get_logposterior", "stanfit", if (object@mode == 1L) { cat("Stan model '", object@model_name, "' is of mode 'test_grad';\n", "sampling is not conducted.\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } else if (object@mode == 2L) { - cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') - return(invisible(NULL)) - } + cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') + return(invisible(NULL)) + } - llp <- lapply(object@sim$samples, function(x) x[['lp__']]) + llp <- lapply(object@sim$samples, function(x) x[['lp__']]) if (inc_warmup) return(llp) if (all(object@sim$warmup2 == 0)) return(llp) - return(mapply(function(x, w) x[-(1:w)], + return(mapply(function(x, w) x[-(1:w)], llp, object@sim$warmup2, - SIMPLIFY = FALSE, USE.NAMES = FALSE)) - }) + SIMPLIFY = FALSE, USE.NAMES = FALSE)) + }) -setGeneric(name = 'get_sampler_params', - def = function(object, ...) { standardGeneric("get_sampler_params")}) +setGeneric(name = 'get_sampler_params', + def = function(object, ...) { standardGeneric("get_sampler_params")}) -setMethod("get_sampler_params", +setMethod("get_sampler_params", definition = function(object, inc_warmup = TRUE) { if (object@mode == 1L) { cat("Stan model '", object@model_name, "' is of mode 'test_grad';\n", "sampling is not conducted.\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } else if (object@mode == 2L) { - cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') - return(invisible(NULL)) + cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') + return(invisible(NULL)) } - + if (isTRUE(object@stan_args[[1L]]$method == "variational")) { stop("'get_sampler_params' not available for ", - "meanfield or fullrank algorithms.") + "meanfield or fullrank algorithms.") } - ldf <- lapply(object@sim$samples, - function(x) do.call(cbind, attr(x, "sampler_params"))) - if (all(sapply(ldf, is.null))) return(invisible(NULL)) + ldf <- lapply(object@sim$samples, + function(x) do.call(cbind, attr(x, "sampler_params"))) + if (all(sapply(ldf, is.null))) return(invisible(NULL)) if (inc_warmup) { if (all(object@sim$warmup2 == 0)) warning("warmup samples not saved") return(invisible(ldf)) } else if (all(object@sim$warmup2 == 0)) return(invisible(ldf)) - return(mapply(function(x, w) - if (w > 0) x[-(1:w), , drop = FALSE] else x, - ldf, object@sim$warmup2, - SIMPLIFY = FALSE, USE.NAMES = FALSE)) - }) + return(mapply(function(x, w) + if (w > 0) x[-(1:w), , drop = FALSE] else x, + ldf, object@sim$warmup2, + SIMPLIFY = FALSE, USE.NAMES = FALSE)) + }) setGeneric(name = 'get_elapsed_time', def = function(object, ...) { standardGeneric("get_elapsed_time")}) @@ -378,19 +378,19 @@ setMethod("get_elapsed_time", t }) -setGeneric(name = 'get_posterior_mean', - def = function(object, ...) { standardGeneric("get_posterior_mean")}) +setGeneric(name = 'get_posterior_mean', + def = function(object, ...) { standardGeneric("get_posterior_mean")}) setMethod("get_posterior_mean", signature = "stanfit", definition = function(object, pars) { if (object@mode == 1L) { cat("Stan model '", object@model_name, "' is of mode 'test_grad';\n", "sampling is not conducted.\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } else if (object@mode == 2L) { - cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') - return(invisible(NULL)) - } + cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') + return(invisible(NULL)) + } fnames <- flatnames(object@model_pars, object@par_dims, col_major = TRUE) if (!exists("posterior_mean_4all", envir = object@.MISC, inherits = FALSE)) { mean_pars <- lapply(object@sim$samples, function(x) attr(x, "mean_pars")) @@ -400,7 +400,7 @@ setMethod("get_posterior_mean", signature = "stanfit", if (ncol(m) > 1) { m <- cbind(m, apply(m, 1, mean)) name_allchains <- "mean-all chains" - } + } cids <- sapply(object@stan_args, function(x) x$chain_id) colnames(m) <- c(paste0('mean-chain:', cids), name_allchains) rownames(m) <- fnames @@ -408,99 +408,99 @@ setMethod("get_posterior_mean", signature = "stanfit", } pars <- if (missing(pars)) object@model_pars else check_pars(c(object@model_pars, fnames), pars) pars <- remove_empty_pars(pars, object@sim$dims_oi) - tidx <- pars_total_indexes(object@model_pars, object@par_dims, fnames, pars) + tidx <- pars_total_indexes(object@model_pars, object@par_dims, fnames, pars) tidx <- lapply(tidx, function(x) attr(x, "row_major_idx")) return(object@.MISC$posterior_mean_4all[unlist(tidx), , drop = FALSE]) }) setGeneric(name = "extract", - def = function(object, ...) { standardGeneric("extract") }) + def = function(object, ...) { standardGeneric("extract") }) setMethod("extract", signature = "stanfit", - definition = function(object, pars, permuted = TRUE, + definition = function(object, pars, permuted = TRUE, inc_warmup = FALSE, include = TRUE) { - # Extract the samples in different forms for different parameters. + # Extract the samples in different forms for different parameters. # # Args: - # object: the object of "stanfit" class - # pars: the names of parameters (including other quantiles) + # object: the object of "stanfit" class + # pars: the names of parameters (including other quantiles) # permuted: if TRUE, the returned samples are permuted without - # warming up. And all the chains are merged. - # inc_warmup: if TRUE, warmup samples are kept; otherwise, - # discarded. If permuted is TRUE, inc_warmup is ignored. + # warming up. And all the chains are merged. + # inc_warmup: if TRUE, warmup samples are kept; otherwise, + # discarded. If permuted is TRUE, inc_warmup is ignored. # include: if FALSE interpret pars as those to exclude # # Returns: # If permuted is TRUE, return an array (matrix) of samples with each - # column being the samples for a parameter. + # column being the samples for a parameter. # If permuted is FALSE, return array with dimensions - # (# of iter (with or w.o. warmup), # of chains, # of flat parameters). + # (# of iter (with or w.o. warmup), # of chains, # of flat parameters). if (object@mode == 1L) { cat("Stan model '", object@model_name, "' is of mode 'test_grad';\n", "sampling is not conducted.\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } else if (object@mode == 2L) { - cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') - return(invisible(NULL)) - } + cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') + return(invisible(NULL)) + } if(!include) pars <- setdiff(object@sim$pars_oi, pars) - pars <- if (missing(pars)) object@sim$pars_oi else check_pars_second(object@sim, pars) + pars <- if (missing(pars)) object@sim$pars_oi else check_pars_second(object@sim, pars) pars <- remove_empty_pars(pars, object@sim$dims_oi) - tidx <- pars_total_indexes(object@sim$pars_oi, - object@sim$dims_oi, - object@sim$fnames_oi, - pars) + tidx <- pars_total_indexes(object@sim$pars_oi, + object@sim$dims_oi, + object@sim$fnames_oi, + pars) n_kept <- object@sim$n_save - object@sim$warmup2 fun1 <- function(par_i) { # sss <- sapply(tidx[[par_i]], get_kept_samples2, object@sim) # if (is.list(sss)) sss <- do.call(c, sss) # the above two lines are slower than the following line of code - sss <- do.call(cbind, lapply(tidx[[par_i]], get_kept_samples2, object@sim)) - dim(sss) <- c(sum(n_kept), object@sim$dims_oi[[par_i]]) + sss <- do.call(cbind, lapply(tidx[[par_i]], get_kept_samples2, object@sim)) + dim(sss) <- c(sum(n_kept), object@sim$dims_oi[[par_i]]) dimnames(sss) <- list(iterations = NULL) - sss - } - + sss + } + if (permuted) { - slist <- lapply(pars, fun1) - names(slist) <- pars - return(slist) - } - - tidx <- unlist(tidx, use.names = FALSE) - tidxnames <- object@sim$fnames_oi[tidx] - sss <- lapply(tidx, get_samples2, object@sim, inc_warmup) + slist <- lapply(pars, fun1) + names(slist) <- pars + return(slist) + } + + tidx <- unlist(tidx, use.names = FALSE) + tidxnames <- object@sim$fnames_oi[tidx] + sss <- lapply(tidx, get_samples2, object@sim, inc_warmup) sss2 <- lapply(sss, function(x) do.call(c, x)) # concatenate samples from different chains - sssf <- unlist(sss2, use.names = FALSE) - - n2 <- object@sim$n_save[1] ## assuming all the chains have equal iter - if (!inc_warmup) n2 <- n2 - object@sim$warmup2[1] - dim(sssf) <- c(n2, object@sim$chains, length(tidx)) + sssf <- unlist(sss2, use.names = FALSE) + + n2 <- object@sim$n_save[1] ## assuming all the chains have equal iter + if (!inc_warmup) n2 <- n2 - object@sim$warmup2[1] + dim(sssf) <- c(n2, object@sim$chains, length(tidx)) cids <- sapply(object@stan_args, function(x) x$chain_id) dimnames(sssf) <- list(iterations = NULL, chains = paste0("chain:", cids), parameters = tidxnames) - sssf - }) + sssf + }) -setMethod("summary", signature = "stanfit", - function(object, pars, - probs = c(0.025, 0.25, 0.50, 0.75, 0.975), use_cache = TRUE, ...) { - # Summarize the samples (that is, compute the mean, SD, quantiles, for +setMethod("summary", signature = "stanfit", + function(object, pars, + probs = c(0.025, 0.25, 0.50, 0.75, 0.975), use_cache = TRUE, ...) { + # Summarize the samples (that is, compute the mean, SD, quantiles, for # the samples in all chains and chains individually after removing # warmup samples, and n_eff and split R hat for all the kept samples.) - # - # Returns: + # + # Returns: # A list with elements: - # summary: the summary for all the kept samples - # c_summary: the summary for individual chains. - # - # Note: + # summary: the summary for all the kept samples + # c_summary: the summary for individual chains. + # + # Note: # This function is not straight in terms of implementation as it # saves some standard summaries including n_eff and Rhat in the - # environment of the object. The summaries and created upon - # the first time standard summary is called and resued later if possible. + # environment of the object. The summaries and created upon + # the first time standard summary is called and resued later if possible. # In addition, the indexes complicate the implementation as internally # we use column major indexes for vector/array parameters. But it might # be better to use row major indexes for the output such as print. @@ -508,11 +508,11 @@ setMethod("summary", signature = "stanfit", if (object@mode == 1L) { cat("Stan model '", object@model_name, "' is of mode 'test_grad';\n", "sampling is not conducted.\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } else if (object@mode == 2L) { - cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') - return(invisible(NULL)) - } + cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') + return(invisible(NULL)) + } if (isTRUE(all.equal(object@sim$n_save, object@sim$warmup2, check.attributes = FALSE, check.names = FALSE))) { @@ -520,117 +520,117 @@ setMethod("summary", signature = "stanfit", return(invisible(NULL)) } - if (!exists("summary", envir = object@.MISC, inherits = FALSE) && use_cache) + if (!exists("summary", envir = object@.MISC, inherits = FALSE) && use_cache) assign("summary", summary_sim(object@sim), envir = object@.MISC) - - pars <- if (missing(pars)) object@sim$pars_oi else check_pars_second(object@sim, pars) + + pars <- if (missing(pars)) object@sim$pars_oi else check_pars_second(object@sim, pars) pars <- remove_empty_pars(pars, object@sim$dims_oi) - if (missing(probs)) - probs <- c(0.025, 0.25, 0.50, 0.75, 0.975) + if (missing(probs)) + probs <- c(0.025, 0.25, 0.50, 0.75, 0.975) if (!use_cache) { # not using the cached (and not create cache, which takes time for too many pars) - ss <- summary_sim(object@sim, pars, probs) - s1 <- cbind(ss$msd[, 1, drop = FALSE], ss$sem, ss$msd[, 2, drop = FALSE], + ss <- summary_sim(object@sim, pars, probs) + s1 <- cbind(ss$msd[, 1, drop = FALSE], ss$sem, ss$msd[, 2, drop = FALSE], ss$quan, ss$ess, ss$rhat) colnames(s1) <- c("mean", "se_mean", "sd", colnames(ss$quan), 'n_eff', 'Rhat') - s2 <- combine_msd_quan(ss$c_msd, ss$c_quan) + s2 <- combine_msd_quan(ss$c_msd, ss$c_quan) idx2 <- match(attr(ss, "row_major_idx"), attr(ss, "col_major_idx")) sf <- list(summary = s1[idx2, , drop = FALSE], c_summary = s2[idx2, , , drop = FALSE]) return(sf) - } + } m <- match(probs, default_summary_probs()) - if (any(is.na(m))) { # unordinary quantiles are requested - ss <- summary_sim_quan(object@sim, pars, probs) - col_idx <- attr(ss, "col_major_idx") - ss$ess <- object@.MISC$summary$ess[col_idx, drop = FALSE] - ss$rhat <- object@.MISC$summary$rhat[col_idx, drop = FALSE] - ss$mean <- object@.MISC$summary$msd[col_idx, 1, drop = FALSE] - ss$sd <- object@.MISC$summary$msd[col_idx, 2, drop = FALSE] - ss$sem <- object@.MISC$summary$sem[col_idx] - s1 <- cbind(ss$mean, ss$sem, ss$sd, + if (any(is.na(m))) { # unordinary quantiles are requested + ss <- summary_sim_quan(object@sim, pars, probs) + col_idx <- attr(ss, "col_major_idx") + ss$ess <- object@.MISC$summary$ess[col_idx, drop = FALSE] + ss$rhat <- object@.MISC$summary$rhat[col_idx, drop = FALSE] + ss$mean <- object@.MISC$summary$msd[col_idx, 1, drop = FALSE] + ss$sd <- object@.MISC$summary$msd[col_idx, 2, drop = FALSE] + ss$sem <- object@.MISC$summary$sem[col_idx] + s1 <- cbind(ss$mean, ss$sem, ss$sd, ss$quan, ss$ess, ss$rhat) colnames(s1) <- c("mean", "se_mean", "sd", colnames(ss$quan), 'n_eff', 'Rhat') - s2 <- combine_msd_quan(object@.MISC$summary$c_msd[col_idx, , , drop = FALSE], ss$c_quan) + s2 <- combine_msd_quan(object@.MISC$summary$c_msd[col_idx, , , drop = FALSE], ss$c_quan) idx2 <- match(attr(ss, "row_major_idx"), col_idx) ss <- list(summary = s1[idx2, , drop = FALSE], c_summary = s2[idx2, , , drop = FALSE]) return(ss) } - tidx <- pars_total_indexes(object@sim$pars_oi, - object@sim$dims_oi, - object@sim$fnames_oi, - pars) + tidx <- pars_total_indexes(object@sim$pars_oi, + object@sim$dims_oi, + object@sim$fnames_oi, + pars) tidx <- lapply(tidx, function(x) attr(x, "row_major_idx")) tidx <- unlist(tidx, use.names = FALSE) tidx_len <- length(tidx) - ss <- object@.MISC$summary + ss <- object@.MISC$summary qnames <- colnames(ss$quan)[m] - - if (!is.null(object@stan_args[[1]]$method) && + + if (!is.null(object@stan_args[[1]]$method) && object@stan_args[[1]]$method == "variational") { s1 <- cbind(ss$msd[tidx, 1, drop = FALSE], - ss$msd[tidx, 2, drop = FALSE], + ss$msd[tidx, 2, drop = FALSE], ss$quan[tidx, m, drop = FALSE]) dim(s1) <- c(length(tidx), length(m) + 2L) colnames(s1) <- c("mean", "sd", qnames) } else { - s1 <- cbind(ss$msd[tidx, 1, drop = FALSE], - ss$sem[tidx, drop = FALSE], - ss$msd[tidx, 2, drop = FALSE], - ss$quan[tidx, m, drop = FALSE], + s1 <- cbind(ss$msd[tidx, 1, drop = FALSE], + ss$sem[tidx, drop = FALSE], + ss$msd[tidx, 2, drop = FALSE], + ss$quan[tidx, m, drop = FALSE], ss$ess[tidx, drop = FALSE], ss$rhat[tidx, drop = FALSE]) - dim(s1) <- c(length(tidx), length(m) + 5L) + dim(s1) <- c(length(tidx), length(m) + 5L) colnames(s1) <- c("mean", "se_mean", "sd", qnames, 'n_eff', 'Rhat') } - pars_names <- rownames(ss$msd)[tidx] - rownames(s1) <- pars_names - s2 <- combine_msd_quan(ss$c_msd[tidx, , , drop = FALSE], ss$c_quan[tidx, m, , drop = FALSE]) + pars_names <- rownames(ss$msd)[tidx] + rownames(s1) <- pars_names + s2 <- combine_msd_quan(ss$c_msd[tidx, , , drop = FALSE], ss$c_quan[tidx, m, , drop = FALSE]) # dim(s2) <- c(tidx_len, length(m) + 2, object@sim$chains) - # dimnames(s2) <- list(parameter = pars_names, - # stats = c("mean", "sd", qnames), NULL) - ss <- list(summary = s1, c_summary = s2) + # dimnames(s2) <- list(parameter = pars_names, + # stats = c("mean", "sd", qnames), NULL) + ss <- list(summary = s1, c_summary = s2) return(ss) - }) + }) if (!isGeneric("traceplot")) { setGeneric(name = "traceplot", - def = function(object, ...) { standardGeneric("traceplot") }) -} + def = function(object, ...) { standardGeneric("traceplot") }) +} if (!isGeneric("log_prob")) { - setGeneric(name = "log_prob", - def = function(object, ...) { standardGeneric("log_prob") }) -} + setGeneric(name = "log_prob", + def = function(object, ...) { standardGeneric("log_prob") }) +} if (!isGeneric("unconstrain_pars")) { - setGeneric(name = "unconstrain_pars", - def = function(object, ...) { standardGeneric("unconstrain_pars") }) -} + setGeneric(name = "unconstrain_pars", + def = function(object, ...) { standardGeneric("unconstrain_pars") }) +} -setMethod("unconstrain_pars", signature = "stanfit", +setMethod("unconstrain_pars", signature = "stanfit", function(object, pars) { # pars is a list as in specifying inits for a chain - if (!is_sfinstance_valid(object)) + if (!is_sfinstance_valid(object)) stop("the model object is not created or not valid") object@.MISC$stan_fit_instance$unconstrain_pars(pars) }) if (!isGeneric("constrain_pars")) { - setGeneric(name = "constrain_pars", - def = function(object, ...) { standardGeneric("constrain_pars") }) -} + setGeneric(name = "constrain_pars", + def = function(object, ...) { standardGeneric("constrain_pars") }) +} -setMethod("constrain_pars", signature = "stanfit", +setMethod("constrain_pars", signature = "stanfit", function(object, upars) { - # upars is a vector on the unconstrained space (R^N*), - # where N* is the number of unconstrained parameters. - if (!is_sfinstance_valid(object)) + # upars is a vector on the unconstrained space (R^N*), + # where N* is the number of unconstrained parameters. + if (!is_sfinstance_valid(object)) stop("the model object is not created or not valid") p <- object@.MISC$stan_fit_instance$constrain_pars(upars) idx_wo_lp <- which(object@model_pars != 'lp__') @@ -638,140 +638,152 @@ setMethod("constrain_pars", signature = "stanfit", }) -setMethod("log_prob", signature = "stanfit", +setMethod("log_prob", signature = "stanfit", function(object, upars, adjust_transform = TRUE, gradient = FALSE) { - if (!is_sfinstance_valid(object)) + if (!is_sfinstance_valid(object)) stop("the model object is not created or not valid") return(object@.MISC$stan_fit_instance$log_prob(upars, adjust_transform, gradient)) - }) + }) if (!isGeneric("get_num_upars")) { - setGeneric(name = "get_num_upars", - def = function(object, ...) { standardGeneric("get_num_upars") }) -} + setGeneric(name = "get_num_upars", + def = function(object, ...) { standardGeneric("get_num_upars") }) +} -setMethod("get_num_upars", signature = "stanfit", +setMethod("get_num_upars", signature = "stanfit", function(object) { - if (!is_sfinstance_valid(object)) + if (!is_sfinstance_valid(object)) stop("the model object is not created or not valid") object@.MISC$stan_fit_instance$num_pars_unconstrained() }) if (!isGeneric("grad_log_prob")) { - setGeneric(name = "grad_log_prob", - def = function(object, ...) { standardGeneric("grad_log_prob") }) -} + setGeneric(name = "grad_log_prob", + def = function(object, ...) { standardGeneric("grad_log_prob") }) +} -setMethod("grad_log_prob", signature = "stanfit", +setMethod("grad_log_prob", signature = "stanfit", function(object, upars, adjust_transform = TRUE) { - if (!is_sfinstance_valid(object)) + if (!is_sfinstance_valid(object)) stop("the model object is not created or not valid") - object@.MISC$stan_fit_instance$grad_log_prob(upars, adjust_transform) - }) + object@.MISC$stan_fit_instance$grad_log_prob(upars, adjust_transform) + }) + +if (!isGeneric("hessian_log_prob")) { + setGeneric(name = "hessian_log_prob", + def = function(object, ...) { standardGeneric("hessian_log_prob") }) +} + +setMethod("hessian_log_prob", signature = "stanfit", + function(object, upars, adjust_transform = TRUE) { + if (!is_sfinstance_valid(object)) + stop("the model object is not created or not valid") + object@.MISC$stan_fit_instance$hessian_log_prob(upars, adjust_transform) + }) -setMethod("traceplot", signature = "stanfit", +setMethod("traceplot", signature = "stanfit", function(object, pars, include = TRUE, unconstrain = FALSE, inc_warmup = FALSE, window = NULL, nrow = NULL, ncol = NULL, - ...) { + ...) { if (object@mode == 1L) { cat("Stan model '", object@model_name, "' is of mode 'test_grad';\n", "sampling is not conducted.\n", sep = '') - return(invisible(NULL)) + return(invisible(NULL)) } else if (object@mode == 2L) { - cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') - return(invisible(NULL)) - } - args <- list(object = object, include = include, + cat("Stan model '", object@model_name, "' does not contain samples.\n", sep = '') + return(invisible(NULL)) + } + args <- list(object = object, include = include, unconstrain = unconstrain, inc_warmup = inc_warmup, nrow = nrow, ncol = ncol, window = window, ...) - if (!missing(pars)) { + if (!missing(pars)) { if ("log-posterior" %in% pars) pars[which(pars == "log-posterior")] <- "lp__" args$pars <- pars } do.call("stan_trace", args) - }) + }) is_sf_valid <- function(sf) { # Similar to is_sm_valid, this is only to test whether - # the compiled DSO is loaded - return(is_sm_valid(sf@stanmodel)) -} + # the compiled DSO is loaded + return(is_sm_valid(sf@stanmodel)) +} is_sfinstance_valid <- function(object) { # Args - # object: an instance of S4 class stanfit - exists("stan_fit_instance", envir = object@.MISC, inherits = FALSE) && + # object: an instance of S4 class stanfit + exists("stan_fit_instance", envir = object@.MISC, inherits = FALSE) && !(is_null_ptr(object@.MISC$stan_fit_instance@.xData$.pointer) || - is_null_ptr(object@.MISC$stan_fit_instance@.xData$.module) || - is_null_ptr(object@.MISC$stan_fit_instance@.xData$.cppclass)) -} + is_null_ptr(object@.MISC$stan_fit_instance@.xData$.module) || + is_null_ptr(object@.MISC$stan_fit_instance@.xData$.cppclass)) +} sflist2stanfit <- function(sflist) { # merge a list of stanfit objects into one # Args: - # sflist, a list of stanfit objects, each element of which - # should have equal length of `iter`, `warmup`, and `thin`. - # Returns: - # A new stanfit objects have all the chains in each element of sf_list. - # The date would be where the new object is created. - # Note: - # * method get_seed would not work well for this merged + # sflist, a list of stanfit objects, each element of which + # should have equal length of `iter`, `warmup`, and `thin`. + # Returns: + # A new stanfit objects have all the chains in each element of sf_list. + # The date would be where the new object is created. + # Note: + # * method get_seed would not work well for this merged # stanfit object in that it only returns the seed used - # for the first object. But all the information is still there. - # * When print function is called, the sampler info is obtained - # only from the first chain. - # - sf_len <- length(sflist) + # for the first object. But all the information is still there. + # * When print function is called, the sampler info is obtained + # only from the first chain. + # + sf_len <- length(sflist) if (sf_len == 0) stop("'sflist' should have at least 1 element") - if (!is.list(sflist) || - any(sapply(sflist, function(x) !is(x, "stanfit")))) { + if (!is.list(sflist) || + any(sapply(sflist, function(x) !is(x, "stanfit")))) { stop("'sflist' must be a list of 'stanfit' objects") } non_zero_modes_idx <- which(sapply(sflist, function(x) x@mode) > 0) - if (length(non_zero_modes_idx) > 0) { + if (length(non_zero_modes_idx) > 0) { stop("The following elements of 'sflist' do not contain samples: ", - paste(non_zero_modes_idx, collapse = ', '), ".") - } + paste(non_zero_modes_idx, collapse = ', '), ".") + } if (sf_len == 1) return(sflist[[1]]) - for (i in 2:sf_len) { - if (!identical(sflist[[i]]@sim$pars_oi, sflist[[1]]@sim$pars_oi) || - !identical(sflist[[i]]@sim$dims_oi, sflist[[1]]@sim$dims_oi)) - stop("parameters in element ", i, " (stanfit object) are different from in element 1") - if (sflist[[i]]@sim$n_save[1] != sflist[[1]]@sim$n_save[1] || + for (i in 2:sf_len) { + if (!identical(sflist[[i]]@sim$pars_oi, sflist[[1]]@sim$pars_oi) || + !identical(sflist[[i]]@sim$dims_oi, sflist[[1]]@sim$dims_oi)) + stop("parameters in element ", i, " (stanfit object) are different from in element 1") + if (sflist[[i]]@sim$n_save[1] != sflist[[1]]@sim$n_save[1] || sflist[[i]]@sim$warmup2[1] != sflist[[1]]@sim$warmup2[1]) - stop("all 'stanfit' objects should have equal length of iterations and warmup") - } - n_chains = sum(sapply(sflist, function(x) x@sim$chains)) - sim = list(samples = do.call(c, lapply(sflist, function(x) x@sim$samples)), - chains = n_chains, - iter = sflist[[1]]@sim$iter, - thin = sflist[[1]]@sim$thin, + stop("all 'stanfit' objects should have equal length of iterations and warmup") + } + n_chains = sum(sapply(sflist, function(x) x@sim$chains)) + sim = list(samples = do.call(c, lapply(sflist, function(x) x@sim$samples)), + chains = n_chains, + iter = sflist[[1]]@sim$iter, + thin = sflist[[1]]@sim$thin, warmup = sflist[[1]]@sim$warmup, - n_save = rep(sflist[[1]]@sim$n_save[1], n_chains), - warmup2 = rep(sflist[[1]]@sim$warmup2[1], n_chains), - permutation = do.call(c, lapply(sflist, function(x) x@sim$permutation)), - pars_oi = sflist[[1]]@sim$pars_oi, - dims_oi = sflist[[1]]@sim$dims_oi, + n_save = rep(sflist[[1]]@sim$n_save[1], n_chains), + warmup2 = rep(sflist[[1]]@sim$warmup2[1], n_chains), + permutation = do.call(c, lapply(sflist, function(x) x@sim$permutation)), + pars_oi = sflist[[1]]@sim$pars_oi, + dims_oi = sflist[[1]]@sim$dims_oi, fnames_oi = sflist[[1]]@sim$fnames_oi, - n_flatnames = sflist[[1]]@sim$n_flatnames) - nfit <- new("stanfit", - model_name = sflist[[1]]@model_name, + n_flatnames = sflist[[1]]@sim$n_flatnames) + nfit <- new("stanfit", + model_name = sflist[[1]]@model_name, model_pars = sflist[[1]]@model_pars, - par_dims = sflist[[1]]@par_dims, - mode = 0L, - sim = sim, - inits = do.call(c, lapply(sflist, function(x) x@inits)), - stan_args = do.call(c, lapply(sflist, function(x) x@stan_args)), - stanmodel = sflist[[1]]@stanmodel, - date = date(), - .MISC = new.env(parent = emptyenv())) + par_dims = sflist[[1]]@par_dims, + mode = 0L, + sim = sim, + inits = do.call(c, lapply(sflist, function(x) x@inits)), + stan_args = do.call(c, lapply(sflist, function(x) x@stan_args)), + stanmodel = sflist[[1]]@stanmodel, + date = date(), + .MISC = new.env(parent = emptyenv())) return(nfit) -} +} # sflist2stan(list(l1=ss1, l2=ss2)) names.stanfit <- function(x) dimnames(x)$parameters @@ -792,14 +804,14 @@ names.stanfit <- function(x) dimnames(x)$parameters } as.array.stanfit <- function(x, ...) { - if (x@mode != 0) return(numeric(0)) + if (x@mode != 0) return(numeric(0)) out <- extract(x, permuted = FALSE, inc_warmup = FALSE, ...) # dimnames(out) <- dimnames(x) return(out) -} +} as.matrix.stanfit <- function(x, ...) { - if (x@mode != 0) return(numeric(0)) - e <- extract(x, permuted = FALSE, inc_warmup = FALSE, ...) + if (x@mode != 0) return(numeric(0)) + e <- extract(x, permuted = FALSE, inc_warmup = FALSE, ...) if (is.null(e)) return(e) enames <- dimnames(e) edim <- dim(e) @@ -807,29 +819,29 @@ as.matrix.stanfit <- function(x, ...) { dimnames(e) <- enames[-2] e } - + as.data.frame.stanfit <- function(x, ...) { return( as.data.frame(as.matrix(x, ...)) ) } dim.stanfit <- function(x) { - if (x@mode != 0) return(numeric(0)) + if (x@mode != 0) return(numeric(0)) c(x@sim$n_save[1] - x@sim$warmup2[1], x@sim$chains, x@sim$n_flatnames) } setGeneric("as.mcmc.list", function(object, ...) standardGeneric("as.mcmc.list")) as.mcmc.list.stanfit <- function(object, pars, ...) { - pars <- if (missing(pars)) object@sim$pars_oi else check_pars_second(object@sim, pars) + pars <- if (missing(pars)) object@sim$pars_oi else check_pars_second(object@sim, pars) pars <- remove_empty_pars(pars, object@sim$dims_oi) tidx <- pars_total_indexes(object@sim$pars_oi, object@sim$dims_oi, object@sim$fnames_oi, pars) tidx <- lapply(tidx, function(x) attr(x, "row_major_idx")) tidx <- unlist(tidx, use.names = FALSE) lst <- vector("list", object@sim$chains) - for (ic in 1:object@sim$chains) { + for (ic in 1:object@sim$chains) { x <- do.call(cbind, object@sim$samples[[ic]])[,tidx,drop=FALSE] - warmup2 <- object@sim$warmup2[ic] + warmup2 <- object@sim$warmup2[ic] if (warmup2 > 0) x <- x[-(1:warmup2), ] x <- as.matrix(x) end <- object@sim$iter @@ -837,7 +849,7 @@ as.mcmc.list.stanfit <- function(object, pars, ...) { start <- end - (nrow(x) - 1) * thin class(x) <- 'mcmc' attr(x, "mcpar") <- c(start, end, thin) - lst[[ic]] <- x + lst[[ic]] <- x } class(lst) <- "mcmc.list" return(lst) @@ -853,9 +865,8 @@ As.mcmc.list <- function(object, pars, include = TRUE, ...) { } dimnames.stanfit <- function(x) { - if (x@mode != 0) return(character(0)) + if (x@mode != 0) return(character(0)) cids <- sapply(x@stan_args, function(x) x$chain_id) - list(iterations = NULL, chains = paste0("chain:", cids), parameters = x@sim$fnames_oi) + list(iterations = NULL, chains = paste0("chain:", cids), parameters = x@sim$fnames_oi) } is.array.stanfit <- function(x) return(x@mode == 0) - diff --git a/rstan/rstan/inst/include/rstan/rcpp_module_def_for_rstan.hpp b/rstan/rstan/inst/include/rstan/rcpp_module_def_for_rstan.hpp index 0df4f526d..601a2c53b 100644 --- a/rstan/rstan/inst/include/rstan/rcpp_module_def_for_rstan.hpp +++ b/rstan/rstan/inst/include/rstan/rcpp_module_def_for_rstan.hpp @@ -25,6 +25,8 @@ RCPP_MODULE(stan_fit4%model_name%_mod){ &rstan::stan_fit<%model_name%_namespace::%model_name%, boost::random::ecuyer1988>::param_oi_tidx) .method("grad_log_prob", &rstan::stan_fit<%model_name%_namespace::%model_name%, boost::random::ecuyer1988>::grad_log_prob) + .method("hessian_log_prob", + &rstan::stan_fit<%model_name%_namespace::%model_name%, boost::random::ecuyer1988>::hessian_log_prob) .method("log_prob", &rstan::stan_fit<%model_name%_namespace::%model_name%, boost::random::ecuyer1988>::log_prob) .method("unconstrain_pars", diff --git a/rstan/rstan/inst/include/rstan/stan_fit.hpp b/rstan/rstan/inst/include/rstan/stan_fit.hpp index 70bcff418..c9a9f1bd1 100644 --- a/rstan/rstan/inst/include/rstan/stan_fit.hpp +++ b/rstan/rstan/inst/include/rstan/stan_fit.hpp @@ -997,7 +997,7 @@ class stan_fit { get_all_flatnames(names_oi_, dims_oi_, fnames_oi_, true); // get_all_indices_col2row(dims_, midx_for_col2row); } - + stan_fit(SEXP data, SEXP seed, SEXP cxxf) : data_(data), model_(data_, Rcpp::as(seed), &rstan::io::rcout), @@ -1173,6 +1173,46 @@ class stan_fit { END_RCPP } + + /** + * Expose the hessian_log_prob of the model to stan_fit so R user + * can call this function. + * + * @param upar The real parameters on the unconstrained + * space. + * @param jacobian_adjust_transform TRUE/FALSE, whether + * we add the term due to the transform from constrained + * space to unconstrained space implicitly done in Stan. + */ + SEXP hessian_log_prob(SEXP upar, SEXP jacobian_adjust_transform) { + BEGIN_RCPP + std::vector par_r = Rcpp::as >(upar); + if (par_r.size() != model_.num_params_r()) { + std::stringstream msg; + msg << "Number of unconstrained parameters does not match " + "that of the model (" + << par_r.size() << " vs " + << model_.num_params_r() + << ")."; + throw std::domain_error(msg.str()); + } + std::vector par_i(model_.num_params_i(), 0); + std::vector gradient; + double lp; + // TODO: change to log_prob_hessian. + if (Rcpp::as(jacobian_adjust_transform)) + lp = stan::model::log_prob_grad(model_, par_r, par_i, gradient, &rstan::io::rcout); + else + lp = stan::model::log_prob_grad(model_, par_r, par_i, gradient, &rstan::io::rcout); + Rcpp::NumericVector grad = Rcpp::wrap(gradient); + grad.attr("log_prob") = lp; + SEXP __sexp_result; + PROTECT(__sexp_result = Rcpp::wrap(grad)); + UNPROTECT(1); + return __sexp_result; + END_RCPP + } + /** * Return the number of unconstrained parameters */ From be52c5b610c9e8a7c31819c86213cc94e89c2b94 Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Mon, 23 Apr 2018 22:03:55 -0700 Subject: [PATCH 03/11] Actually use the hessian method --- rstan/rstan/inst/include/rstan/stan_fit.hpp | 17 ++++++++++++----- 1 file changed, 12 insertions(+), 5 deletions(-) diff --git a/rstan/rstan/inst/include/rstan/stan_fit.hpp b/rstan/rstan/inst/include/rstan/stan_fit.hpp index c9a9f1bd1..9e626abef 100644 --- a/rstan/rstan/inst/include/rstan/stan_fit.hpp +++ b/rstan/rstan/inst/include/rstan/stan_fit.hpp @@ -1197,15 +1197,22 @@ class stan_fit { throw std::domain_error(msg.str()); } std::vector par_i(model_.num_params_i(), 0); - std::vector gradient; + Eigen::Matrix grad_f, + Eigen::Matrix hess_f, + double lp; - // TODO: change to log_prob_hessian. if (Rcpp::as(jacobian_adjust_transform)) - lp = stan::model::log_prob_grad(model_, par_r, par_i, gradient, &rstan::io::rcout); + lp = stan::model::log_prob_hessian( + model_, par_r, lp, grad_f, hess_f, &rstan::io::rcout); else - lp = stan::model::log_prob_grad(model_, par_r, par_i, gradient, &rstan::io::rcout); - Rcpp::NumericVector grad = Rcpp::wrap(gradient); + lp = stan::model::log_prob_hessian( + model_, par_r, lp, grad_f, hess_f, &rstan::io::rcout); + Rcpp::NumericVector grad = Rcpp::wrap(grad_f); + Rcpp::NumericMatrix hess = Rcpp::wrap(hess_f); + + // TODO: do this differently grad.attr("log_prob") = lp; + grad.attr("hessian") = hess; SEXP __sexp_result; PROTECT(__sexp_result = Rcpp::wrap(grad)); UNPROTECT(1); From 232f85fb57547338339acfd5e950f4de64d7c9d9 Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Mon, 23 Apr 2018 22:51:34 -0700 Subject: [PATCH 04/11] Make a functor template --- rstan/rstan/R/git_rstan_head.R | 4 ++++ rstan/rstan/inst/include/rstan/stan_fit.hpp | 7 +++++-- 2 files changed, 9 insertions(+), 2 deletions(-) diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index 4c4106380..a7cb8f4e7 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -33,3 +33,7 @@ git_head <- function() "50a10cca0de66ee51cf3832d1aa6ef29b3590659" git_head <- function() "57444a10c088dbc52105fb394aadb996cd0c9df0" git_head <- function() "2165604262c9085f0e47408d60bb04c76e5485c3" git_head <- function() "2165604262c9085f0e47408d60bb04c76e5485c3" +git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" +git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" +git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" +git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" diff --git a/rstan/rstan/inst/include/rstan/stan_fit.hpp b/rstan/rstan/inst/include/rstan/stan_fit.hpp index 9e626abef..914d65e50 100644 --- a/rstan/rstan/inst/include/rstan/stan_fit.hpp +++ b/rstan/rstan/inst/include/rstan/stan_fit.hpp @@ -26,6 +26,9 @@ // void R_CheckUserInterrupt(void); +// TODO: where is this actually supposed to go? +#include + // REF: cmdstan: src/cmdstan/command.hpp #include #include @@ -1197,8 +1200,8 @@ class stan_fit { throw std::domain_error(msg.str()); } std::vector par_i(model_.num_params_i(), 0); - Eigen::Matrix grad_f, - Eigen::Matrix hess_f, + Eigen::Matrix grad_f; + Eigen::Matrix hess_f; double lp; if (Rcpp::as(jacobian_adjust_transform)) From 2d7c20b4be33936c67b88215245a96f4255b3748 Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Mon, 23 Apr 2018 23:16:37 -0700 Subject: [PATCH 05/11] Use templated functor --- rstan/rstan/R/git_rstan_head.R | 3 +++ rstan/rstan/inst/include/rstan/stan_fit.hpp | 10 ++++++---- 2 files changed, 9 insertions(+), 4 deletions(-) diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index a7cb8f4e7..741f28680 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -37,3 +37,6 @@ git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" +git_head <- function() "232f85fb57547338339acfd5e950f4de64d7c9d9" +git_head <- function() "232f85fb57547338339acfd5e950f4de64d7c9d9" +git_head <- function() "232f85fb57547338339acfd5e950f4de64d7c9d9" diff --git a/rstan/rstan/inst/include/rstan/stan_fit.hpp b/rstan/rstan/inst/include/rstan/stan_fit.hpp index 914d65e50..3bc519b5a 100644 --- a/rstan/rstan/inst/include/rstan/stan_fit.hpp +++ b/rstan/rstan/inst/include/rstan/stan_fit.hpp @@ -1189,7 +1189,9 @@ class stan_fit { */ SEXP hessian_log_prob(SEXP upar, SEXP jacobian_adjust_transform) { BEGIN_RCPP - std::vector par_r = Rcpp::as >(upar); + Eigen::Matrix par_r = + Rcpp::as>(upar); + //std::vector par_r = Rcpp::as >(upar); if (par_r.size() != model_.num_params_r()) { std::stringstream msg; msg << "Number of unconstrained parameters does not match " @@ -1199,16 +1201,16 @@ class stan_fit { << ")."; throw std::domain_error(msg.str()); } - std::vector par_i(model_.num_params_i(), 0); + //std::vector par_i(model_.num_params_i(), 0); Eigen::Matrix grad_f; Eigen::Matrix hess_f; double lp; if (Rcpp::as(jacobian_adjust_transform)) - lp = stan::model::log_prob_hessian( + stan::model::log_prob_hessian( model_, par_r, lp, grad_f, hess_f, &rstan::io::rcout); else - lp = stan::model::log_prob_hessian( + stan::model::log_prob_hessian( model_, par_r, lp, grad_f, hess_f, &rstan::io::rcout); Rcpp::NumericVector grad = Rcpp::wrap(grad_f); Rcpp::NumericMatrix hess = Rcpp::wrap(hess_f); From e56060c723e942519335c25d4b41741ec9296af2 Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Mon, 23 Apr 2018 23:41:39 -0700 Subject: [PATCH 06/11] Fix conversion from SEXP --- rstan/rstan/R/git_rstan_head.R | 5 +++++ rstan/rstan/inst/include/rstan/stan_fit.hpp | 6 ++++-- 2 files changed, 9 insertions(+), 2 deletions(-) diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index 741f28680..51acdcf93 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -40,3 +40,8 @@ git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" git_head <- function() "232f85fb57547338339acfd5e950f4de64d7c9d9" git_head <- function() "232f85fb57547338339acfd5e950f4de64d7c9d9" git_head <- function() "232f85fb57547338339acfd5e950f4de64d7c9d9" +git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" +git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" +git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" +git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" +git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" diff --git a/rstan/rstan/inst/include/rstan/stan_fit.hpp b/rstan/rstan/inst/include/rstan/stan_fit.hpp index 3bc519b5a..9e764babe 100644 --- a/rstan/rstan/inst/include/rstan/stan_fit.hpp +++ b/rstan/rstan/inst/include/rstan/stan_fit.hpp @@ -1189,9 +1189,11 @@ class stan_fit { */ SEXP hessian_log_prob(SEXP upar, SEXP jacobian_adjust_transform) { BEGIN_RCPP + std::vector std_par_r = Rcpp::as >(upar); Eigen::Matrix par_r = - Rcpp::as>(upar); - //std::vector par_r = Rcpp::as >(upar); + Eigen::Map>( + std_par_r.data(), std_par_r.size()); + if (par_r.size() != model_.num_params_r()) { std::stringstream msg; msg << "Number of unconstrained parameters does not match " From 9a2bbd719e4eabb86453ffaa8d61712a2c727b7c Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Tue, 24 Apr 2018 09:59:46 -0700 Subject: [PATCH 07/11] Include the Rcppeigen headers --- install_StanHeaders.R | 8 +++- rstan/makefile | 46 ++++++++++----------- rstan/rstan/R/git_rstan_head.R | 1 + rstan/rstan/inst/include/rstan/stan_fit.hpp | 3 +- 4 files changed, 31 insertions(+), 27 deletions(-) diff --git a/install_StanHeaders.R b/install_StanHeaders.R index 89ac6948c..3be94cf9a 100644 --- a/install_StanHeaders.R +++ b/install_StanHeaders.R @@ -4,12 +4,16 @@ path_rstan <- tempfile(pattern = "git2r-") path_stan_dev <- "/home/rgiordan/Documents/git_repos/stan-dev/" + git2r::clone(file.path(path_stan_dev, "rstan"), path_rstan, branch = "add_hessian4") + git2r::clone(file.path(path_stan_dev, "stan"), file.path(path_rstan, "StanHeaders", "inst", "include", "upstream"), - branch = "add_hessians") # may want to change this branch + branch = "add_hessians2") + git2r::clone(file.path(path_stan_dev, "math"), file.path(path_rstan, "StanHeaders", "inst", "include", "mathlib"), - branch = "develop") # may want to change this branch + branch = "develop") + devtools::install(file.path(path_rstan, "StanHeaders"), args = "--preclean") diff --git a/rstan/makefile b/rstan/makefile index 2df2a628c..6ab30f728 100644 --- a/rstan/makefile +++ b/rstan/makefile @@ -16,33 +16,33 @@ ifeq ($(uname_S), Linux) endif ifdef R_HOME - R = $(R_HOME)/bin/R - RSCRIPT = $(R_HOME)/bin/Rscript + R = $(R_HOME)/bin/R + RSCRIPT = $(R_HOME)/bin/Rscript else # R = R --arch x86_64 - R = R - RSCRIPT = Rscript + R = R + RSCRIPT = Rscript endif R_MAKEVARS_USER = $(CURDIR)/R_Makevars export R_MAKEVARS_USER -BUILD_OPTS ?= +BUILD_OPTS ?= RSTANVER := $(shell $(RSCRIPT) -e "cat(read.dcf(file = './rstan/DESCRIPTION')[1, deparse(quote(Version))])") -STANPKG := rstan_$(RSTANVER).tar.gz +STANPKG := rstan_$(RSTANVER).tar.gz GIT = $(strip $(shell $(WH) git)) -# The default target of this makefile +# The default target of this makefile help: .PHONY: help -help: +help: @echo '--------------------------------------------------------------------------------' @echo 'main targets:' @echo ' build: to build rstan source package using R CMD build' - @echo ' install: to install rstan using R CMD INSTALL' + @echo ' install: to install rstan using R CMD INSTALL' @echo ' uninstall: to install rstan using R CMD REMOVE rstan' @echo ' clean: to remove the gzip package file' @echo ' check: use R CMD check on rstan package' @@ -55,7 +55,7 @@ help: .PHONY: build check install clean clean-all install_pre_rpkg test-R example_csv install_pre_rpkg: - @R -q -e "options(repos=structure(c(CRAN = 'http://cran.rstudio.com'))); for (pkg in c('inline', 'Rcpp', 'RcppEigen', 'RUnit', 'BH', 'StanHeaders', 'RInside')) if (!require(pkg, character.only = TRUE)) install.packages(pkg, dep = TRUE); sessionInfo()" + @R -q -e "options(repos=structure(c(CRAN = 'http://cran.rstudio.com'))); for (pkg in c('inline', 'Rcpp', 'RcppEigen', 'RUnit', 'BH', 'StanHeaders', 'RInside')) if (!require(pkg, character.only = TRUE)) install.packages(pkg, dep = TRUE); sessionInfo()" build $(STANPKG): ./rstan/DESCRIPTION ifeq ($(GIT),) @@ -65,26 +65,26 @@ else endif $(R) CMD build rstan --md5 $(BUILD_OPTS) --no-build-vignettes --no-manual -example_csv: +example_csv: ifneq ($(OS), win) @if test -f "example/rstan_doc_ex_1.csv"; then cp example/rstan_doc_ex_*.csv ./rstan/inst/misc && echo "Copy done."; \ else echo "Not found example/rstan_doc_ex_1.csv, run example/examplemodel.R first."; fi endif -check: $(STANPKG) - $(R) CMD check --as-cran $(STANPKG) +check: $(STANPKG) + $(R) CMD check --as-cran $(STANPKG) -install: $(STANPKG) +install: $(STANPKG) $(R) CMD INSTALL $(STANPKG) uninstall: $(R) CMD REMOVE rstan test-R: - cd tests; R -q -f runRunitTests.R + cd tests; R -q -f runRunitTests.R -clean: - rm -f $(STANPKG) +clean: + rm -f $(STANPKG) rm -rf $(shell find tests/cpp -type f -name '*.d') rm -rf $(shell find tests/cpp -type f -name '*.d.*') @@ -92,9 +92,9 @@ clean-all: clean rm -f $(shell find . -name 'rstan_*.tar.gz') rm -rf rstan.Rcheck rm -rf tests/cpp/lib $(GTEST_LIB) - rm -rf $(patsubst %.cpp,%$(EXE),$(CPPTESTS)) + rm -rf $(patsubst %.cpp,%$(EXE),$(CPPTESTS)) -# buildbin: # build a binary version +# buildbin: # build a binary version # R CMD INSTALL -l ./tmp --build rstan @@ -108,7 +108,7 @@ STANLIB := $(STAN_MATH_SUBMODULE)lib STANSRC := $(STAN_SUBMODULE)src GTEST_LIB := tests/cpp/gtest.o -LIBPTHREAD = +LIBPTHREAD = ifeq ($(OS), linux) LIBPTHREAD += -lpthread endif @@ -116,11 +116,11 @@ endif CPPTESTS := $(shell find tests/cpp -type f -name '*_test.cpp') STAN_INSTANTIATION_FILES := $(patsubst rstan/src/%.cpp,tests/cpp/lib/src/%.o,$(wildcard rstan/src/*.cpp)) -## comment this out if you need a different version of R, +## comment this out if you need a different version of R, ## and set set R_HOME accordingly as an environment variable R_HOME := $(shell R RHOME) -## include headers and libraries for R +## include headers and libraries for R RCPPFLAGS := $(shell $(R_HOME)/bin/R CMD config --cppflags) RLDFLAGS := $(shell $(R_HOME)/bin/R CMD config --ldflags) @@ -132,7 +132,6 @@ RCPPLIBS := $(shell echo 'Rcpp:::LdFlags()' | $(R_HOME)/bin/R --slave) RINSIDEINCL := $(shell $(R_HOME)/bin/Rscript -e 'RInside:::CxxFlags()') RINSIDELIBS := $(shell $(R_HOME)/bin/Rscript -e 'RInside:::LdFlags()') - $(STAN_INSTANTIATION_FILES) $(GTEST_LIB) $(patsubst %.cpp,%$(EXE),$(CPPTESTS)) $(patsubst %.cpp,%.d,$(CPPTESTS)): CPPFLAGS += -Wall $(RCPPFLAGS) $(STAN_INSTANTIATION_FILES) $(GTEST_LIB) $(patsubst %.cpp,%$(EXE),$(CPPTESTS)) $(patsubst %.cpp,%.d,$(CPPTESTS)): CXXFLAGS += $(RCPPFLAGS) $(RCPPINCL) $(RINSIDEINCL) $(shell $(R_HOME)/bin/R CMD config CXXFLAGS) -I $(STAN_MATH_SUBMODULE) -I rstan/inst/include -I $(STANSRC) $(addprefix -isystem ,$(wildcard $(STANLIB)/*)) -isystem $(GTESTPATH)/include -isystem $(GTESTPATH) -ftemplate-depth=256 -DBOOST_RESULT_OF_USE_TR1 -DBOOST_NO_DECLTYPE -DBOOST_DISABLE_ASSERTS @@ -174,4 +173,3 @@ endif ifneq (,$(filter $(patsubst %.cpp,%$(EXE),$(CPPTESTS)),$(MAKECMDGOALS))) -include $(patsubst %$(EXE),%.d,$(filter $(patsubst %.cpp,%$(EXE),$(CPPTESTS)),$(MAKECMDGOALS))) endif - diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index 51acdcf93..b911e8e9e 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -45,3 +45,4 @@ git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" +git_head <- function() "e56060c723e942519335c25d4b41741ec9296af2" diff --git a/rstan/rstan/inst/include/rstan/stan_fit.hpp b/rstan/rstan/inst/include/rstan/stan_fit.hpp index 9e764babe..6b3d91057 100644 --- a/rstan/rstan/inst/include/rstan/stan_fit.hpp +++ b/rstan/rstan/inst/include/rstan/stan_fit.hpp @@ -26,8 +26,9 @@ // void R_CheckUserInterrupt(void); -// TODO: where is this actually supposed to go? +// TODO: where are these includes actually supposed to go? #include +#include // REF: cmdstan: src/cmdstan/command.hpp #include From 6aeef54d87f66376dffcac152b392d671bb4b44c Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Tue, 24 Apr 2018 10:12:59 -0700 Subject: [PATCH 08/11] Expose HVP methods --- rstan/rstan/NAMESPACE | 2 +- rstan/rstan/R/git_rstan_head.R | 5 ++ rstan/rstan/R/stanfit-class.R | 12 ++++ .../rstan/rcpp_module_def_for_rstan.hpp | 2 + rstan/rstan/inst/include/rstan/stan_fit.hpp | 57 +++++++++++++++++++ 5 files changed, 77 insertions(+), 1 deletion(-) diff --git a/rstan/rstan/NAMESPACE b/rstan/rstan/NAMESPACE index c3c2bb1d7..367a91978 100644 --- a/rstan/rstan/NAMESPACE +++ b/rstan/rstan/NAMESPACE @@ -71,7 +71,7 @@ exportMethods( get_cppcode, get_cxxflags, # for stanmodel show, sampling, summary, extract, traceplot, plot, get_stancode, get_inits, get_seed, get_cppo_mode, - log_prob, grad_log_prob, hessian_log_prob, + log_prob, grad_log_prob, hessian_log_prob, hessian_times_vector_log_prob, unconstrain_pars, constrain_pars, get_num_upars, get_seeds, get_adaptation_info, diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index b911e8e9e..012e5d570 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -46,3 +46,8 @@ git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" git_head <- function() "e56060c723e942519335c25d4b41741ec9296af2" +git_head <- function() "9a2bbd719e4eabb86453ffaa8d61712a2c727b7c" +git_head <- function() "9a2bbd719e4eabb86453ffaa8d61712a2c727b7c" +git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" +git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" +git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" diff --git a/rstan/rstan/R/stanfit-class.R b/rstan/rstan/R/stanfit-class.R index 60d7ba39e..8c1e6ea50 100644 --- a/rstan/rstan/R/stanfit-class.R +++ b/rstan/rstan/R/stanfit-class.R @@ -681,6 +681,18 @@ setMethod("hessian_log_prob", signature = "stanfit", object@.MISC$stan_fit_instance$hessian_log_prob(upars, adjust_transform) }) +if (!isGeneric("hessian_times_vector_log_prob")) { + setGeneric(name = "hessian_times_vector_log_prob", + def = function(object, ...) { standardGeneric("hessian_times_vector_log_prob") }) +} + +setMethod("hessian_times_vector_log_prob", signature = "stanfit", + function(object, upars, v, adjust_transform = TRUE) { + if (!is_sfinstance_valid(object)) + stop("the model object is not created or not valid") + object@.MISC$stan_fit_instance$hessian_times_vector_log_prob(upars, v, adjust_transform) + }) + setMethod("traceplot", signature = "stanfit", function(object, pars, include = TRUE, unconstrain = FALSE, inc_warmup = FALSE, window = NULL, nrow = NULL, ncol = NULL, diff --git a/rstan/rstan/inst/include/rstan/rcpp_module_def_for_rstan.hpp b/rstan/rstan/inst/include/rstan/rcpp_module_def_for_rstan.hpp index 601a2c53b..0e2b553f7 100644 --- a/rstan/rstan/inst/include/rstan/rcpp_module_def_for_rstan.hpp +++ b/rstan/rstan/inst/include/rstan/rcpp_module_def_for_rstan.hpp @@ -25,6 +25,8 @@ RCPP_MODULE(stan_fit4%model_name%_mod){ &rstan::stan_fit<%model_name%_namespace::%model_name%, boost::random::ecuyer1988>::param_oi_tidx) .method("grad_log_prob", &rstan::stan_fit<%model_name%_namespace::%model_name%, boost::random::ecuyer1988>::grad_log_prob) + .method("hessian_times_vector_log_prob", + &rstan::stan_fit<%model_name%_namespace::%model_name%, boost::random::ecuyer1988>::hessian_times_vector_log_prob) .method("hessian_log_prob", &rstan::stan_fit<%model_name%_namespace::%model_name%, boost::random::ecuyer1988>::hessian_log_prob) .method("log_prob", diff --git a/rstan/rstan/inst/include/rstan/stan_fit.hpp b/rstan/rstan/inst/include/rstan/stan_fit.hpp index 6b3d91057..151e63213 100644 --- a/rstan/rstan/inst/include/rstan/stan_fit.hpp +++ b/rstan/rstan/inst/include/rstan/stan_fit.hpp @@ -28,6 +28,7 @@ // TODO: where are these includes actually supposed to go? #include +#include #include // REF: cmdstan: src/cmdstan/command.hpp @@ -1228,6 +1229,62 @@ class stan_fit { END_RCPP } + /** + * Expose the hessian_log_prob of the model to stan_fit so R user + * can call this function. + * + * @param upar The real parameters on the unconstrained + * space. + * @param jacobian_adjust_transform TRUE/FALSE, whether + * we add the term due to the transform from constrained + * space to unconstrained space implicitly done in Stan. + */ + SEXP hessian_times_vector_log_prob( + SEXP upar, SEXP v, SEXP jacobian_adjust_transform) { + + BEGIN_RCPP + + // Is there a way to convert directly from SEXP to an Eigen matrix? + std::vector std_par_r = Rcpp::as >(upar); + Eigen::Matrix par_r = + Eigen::Map>( + std_par_r.data(), std_par_r.size()); + + std::vector v_r = Rcpp::as >(v); + Eigen::Matrix v_vec = + Eigen::Map>( + v_r.data(), v_r.size()); + + if (par_r.size() != model_.num_params_r()) { + std::stringstream msg; + msg << "Number of unconstrained parameters does not match " + "that of the model (" + << par_r.size() << " vs " + << model_.num_params_r() + << ")."; + throw std::domain_error(msg.str()); + } + //std::vector par_i(model_.num_params_i(), 0); + Eigen::Matrix hess_f_dot_v; + + double lp; + if (Rcpp::as(jacobian_adjust_transform)) + stan::model::log_prob_hessian_times_vector( + model_, par_r, v_vec, lp, hess_f_dot_v, &rstan::io::rcout); + else + stan::model::log_prob_hessian_times_vector( + model_, par_r, v_vec, lp, hess_f_dot_v, &rstan::io::rcout); + Rcpp::NumericVector hess_f_dot_v_r = Rcpp::wrap(hess_f_dot_v); + + // TODO: do this differently + hess_f_dot_v_r.attr("log_prob") = lp; + SEXP __sexp_result; + PROTECT(__sexp_result = Rcpp::wrap(hess_f_dot_v_r)); + UNPROTECT(1); + return __sexp_result; + END_RCPP + } + /** * Return the number of unconstrained parameters */ From e1cbf0d794729327c8bc4f1fed5957e49f231420 Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Tue, 24 Apr 2018 10:49:44 -0700 Subject: [PATCH 09/11] Return the hessian direclty not as an attribute --- rstan/rstan/R/git_rstan_head.R | 1 + rstan/rstan/inst/include/rstan/stan_fit.hpp | 8 ++++---- 2 files changed, 5 insertions(+), 4 deletions(-) diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index 012e5d570..32f219dc3 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -51,3 +51,4 @@ git_head <- function() "9a2bbd719e4eabb86453ffaa8d61712a2c727b7c" git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" +git_head <- function() "6aeef54d87f66376dffcac152b392d671bb4b44c" diff --git a/rstan/rstan/inst/include/rstan/stan_fit.hpp b/rstan/rstan/inst/include/rstan/stan_fit.hpp index 151e63213..3eea5ef28 100644 --- a/rstan/rstan/inst/include/rstan/stan_fit.hpp +++ b/rstan/rstan/inst/include/rstan/stan_fit.hpp @@ -1216,14 +1216,14 @@ class stan_fit { else stan::model::log_prob_hessian( model_, par_r, lp, grad_f, hess_f, &rstan::io::rcout); + Rcpp::NumericVector grad = Rcpp::wrap(grad_f); Rcpp::NumericMatrix hess = Rcpp::wrap(hess_f); - // TODO: do this differently - grad.attr("log_prob") = lp; - grad.attr("hessian") = hess; + hess.attr("log_prob") = lp; + hess.attr("gradient") = grad; SEXP __sexp_result; - PROTECT(__sexp_result = Rcpp::wrap(grad)); + PROTECT(__sexp_result = Rcpp::wrap(hess)); UNPROTECT(1); return __sexp_result; END_RCPP From dd3df4e2fc9671cd886e5d60a9b2376b32be63b3 Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Tue, 20 Nov 2018 15:38:39 -0800 Subject: [PATCH 10/11] More notes --- install_StanHeaders.R | 30 ++++++++++++++++++++++-------- rstan/rstan/R/git_rstan_head.R | 6 ++++++ 2 files changed, 28 insertions(+), 8 deletions(-) diff --git a/install_StanHeaders.R b/install_StanHeaders.R index 3be94cf9a..7fc8fd3c3 100644 --- a/install_StanHeaders.R +++ b/install_StanHeaders.R @@ -3,17 +3,31 @@ # including its submodules. This requires the git2r and devtools packages path_rstan <- tempfile(pattern = "git2r-") +print(path_rstan) path_stan_dev <- "/home/rgiordan/Documents/git_repos/stan-dev/" -git2r::clone(file.path(path_stan_dev, "rstan"), - path_rstan, branch = "add_hessian4") +# git2r::clone(file.path(path_stan_dev, "rstan"), +# path_rstan, branch = "add_hessian4") +# +# git2r::clone(file.path(path_stan_dev, "stan"), +# file.path(path_rstan, "StanHeaders", "inst", "include", "upstream"), +# branch = "add_hessians3") -git2r::clone(file.path(path_stan_dev, "stan"), - file.path(path_rstan, "StanHeaders", "inst", "include", "upstream"), - branch = "add_hessians2") + git2r::clone(file.path(path_stan_dev, "rstan"), + path_rstan, branch = "develop") -git2r::clone(file.path(path_stan_dev, "math"), - file.path(path_rstan, "StanHeaders", "inst", "include", "mathlib"), - branch = "develop") + git2r::clone(file.path(path_stan_dev, "stan"), + file.path(path_rstan, "StanHeaders", "inst", "include", "upstream"), + branch = "develop") + +# Currently broken: +# git2r::clone(file.path(path_stan_dev, "math"), +# file.path(path_rstan, "StanHeaders", "inst", "include", "mathlib"), +# branch = "develop") + +#Currently works: + git2r::clone(file.path(path_stan_dev, "math"), + file.path(path_rstan, "StanHeaders", "inst", "include", "mathlib"), + branch = "master") devtools::install(file.path(path_rstan, "StanHeaders"), args = "--preclean") diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index 32f219dc3..8a4e23cb0 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -52,3 +52,9 @@ git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" git_head <- function() "6aeef54d87f66376dffcac152b392d671bb4b44c" +git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" +git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" +git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" +git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" +git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" +git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" From a4d69686af083d762a13bc04351d6791037c58e1 Mon Sep 17 00:00:00 2001 From: Ryan Giordano Date: Tue, 20 Nov 2018 15:53:40 -0800 Subject: [PATCH 11/11] Revert to develop versions of unneeded files --- install_StanHeaders.R | 34 ++++++------------------- rstan/R_Makevars | 1 - rstan/makefile | 46 ++++++++++++++++++---------------- rstan/rstan/R/git_rstan_head.R | 29 --------------------- 4 files changed, 31 insertions(+), 79 deletions(-) diff --git a/install_StanHeaders.R b/install_StanHeaders.R index 7fc8fd3c3..17120a10b 100644 --- a/install_StanHeaders.R +++ b/install_StanHeaders.R @@ -3,31 +3,11 @@ # including its submodules. This requires the git2r and devtools packages path_rstan <- tempfile(pattern = "git2r-") -print(path_rstan) -path_stan_dev <- "/home/rgiordan/Documents/git_repos/stan-dev/" - -# git2r::clone(file.path(path_stan_dev, "rstan"), -# path_rstan, branch = "add_hessian4") -# -# git2r::clone(file.path(path_stan_dev, "stan"), -# file.path(path_rstan, "StanHeaders", "inst", "include", "upstream"), -# branch = "add_hessians3") - - git2r::clone(file.path(path_stan_dev, "rstan"), - path_rstan, branch = "develop") - - git2r::clone(file.path(path_stan_dev, "stan"), - file.path(path_rstan, "StanHeaders", "inst", "include", "upstream"), - branch = "develop") - -# Currently broken: -# git2r::clone(file.path(path_stan_dev, "math"), -# file.path(path_rstan, "StanHeaders", "inst", "include", "mathlib"), -# branch = "develop") - -#Currently works: - git2r::clone(file.path(path_stan_dev, "math"), - file.path(path_rstan, "StanHeaders", "inst", "include", "mathlib"), - branch = "master") - +git2r::clone("http://github.com/stan-dev/rstan", path_rstan, branch = "develop") +git2r::clone("http://github.com/stan-dev/stan", + file.path(path_rstan, "StanHeaders", "inst", "include", "upstream"), + branch = "develop") # may want to change this branch +git2r::clone("http://github.com/stan-dev/math", + file.path(path_rstan, "StanHeaders", "inst", "include", "mathlib"), + branch = "develop") # may want to change this branch devtools::install(file.path(path_rstan, "StanHeaders"), args = "--preclean") diff --git a/rstan/R_Makevars b/rstan/R_Makevars index a16a3b1e0..889fafb5f 100644 --- a/rstan/R_Makevars +++ b/rstan/R_Makevars @@ -6,4 +6,3 @@ MAKEFLAGS=-j6 # CXXFLAGS = -O3 $(LTO) -CXXFLAGS+=-std=c++14 \ No newline at end of file diff --git a/rstan/makefile b/rstan/makefile index 6ab30f728..2df2a628c 100644 --- a/rstan/makefile +++ b/rstan/makefile @@ -16,33 +16,33 @@ ifeq ($(uname_S), Linux) endif ifdef R_HOME - R = $(R_HOME)/bin/R - RSCRIPT = $(R_HOME)/bin/Rscript + R = $(R_HOME)/bin/R + RSCRIPT = $(R_HOME)/bin/Rscript else # R = R --arch x86_64 - R = R - RSCRIPT = Rscript + R = R + RSCRIPT = Rscript endif R_MAKEVARS_USER = $(CURDIR)/R_Makevars export R_MAKEVARS_USER -BUILD_OPTS ?= +BUILD_OPTS ?= RSTANVER := $(shell $(RSCRIPT) -e "cat(read.dcf(file = './rstan/DESCRIPTION')[1, deparse(quote(Version))])") -STANPKG := rstan_$(RSTANVER).tar.gz +STANPKG := rstan_$(RSTANVER).tar.gz GIT = $(strip $(shell $(WH) git)) -# The default target of this makefile +# The default target of this makefile help: .PHONY: help -help: +help: @echo '--------------------------------------------------------------------------------' @echo 'main targets:' @echo ' build: to build rstan source package using R CMD build' - @echo ' install: to install rstan using R CMD INSTALL' + @echo ' install: to install rstan using R CMD INSTALL' @echo ' uninstall: to install rstan using R CMD REMOVE rstan' @echo ' clean: to remove the gzip package file' @echo ' check: use R CMD check on rstan package' @@ -55,7 +55,7 @@ help: .PHONY: build check install clean clean-all install_pre_rpkg test-R example_csv install_pre_rpkg: - @R -q -e "options(repos=structure(c(CRAN = 'http://cran.rstudio.com'))); for (pkg in c('inline', 'Rcpp', 'RcppEigen', 'RUnit', 'BH', 'StanHeaders', 'RInside')) if (!require(pkg, character.only = TRUE)) install.packages(pkg, dep = TRUE); sessionInfo()" + @R -q -e "options(repos=structure(c(CRAN = 'http://cran.rstudio.com'))); for (pkg in c('inline', 'Rcpp', 'RcppEigen', 'RUnit', 'BH', 'StanHeaders', 'RInside')) if (!require(pkg, character.only = TRUE)) install.packages(pkg, dep = TRUE); sessionInfo()" build $(STANPKG): ./rstan/DESCRIPTION ifeq ($(GIT),) @@ -65,26 +65,26 @@ else endif $(R) CMD build rstan --md5 $(BUILD_OPTS) --no-build-vignettes --no-manual -example_csv: +example_csv: ifneq ($(OS), win) @if test -f "example/rstan_doc_ex_1.csv"; then cp example/rstan_doc_ex_*.csv ./rstan/inst/misc && echo "Copy done."; \ else echo "Not found example/rstan_doc_ex_1.csv, run example/examplemodel.R first."; fi endif -check: $(STANPKG) - $(R) CMD check --as-cran $(STANPKG) +check: $(STANPKG) + $(R) CMD check --as-cran $(STANPKG) -install: $(STANPKG) +install: $(STANPKG) $(R) CMD INSTALL $(STANPKG) uninstall: $(R) CMD REMOVE rstan test-R: - cd tests; R -q -f runRunitTests.R + cd tests; R -q -f runRunitTests.R -clean: - rm -f $(STANPKG) +clean: + rm -f $(STANPKG) rm -rf $(shell find tests/cpp -type f -name '*.d') rm -rf $(shell find tests/cpp -type f -name '*.d.*') @@ -92,9 +92,9 @@ clean-all: clean rm -f $(shell find . -name 'rstan_*.tar.gz') rm -rf rstan.Rcheck rm -rf tests/cpp/lib $(GTEST_LIB) - rm -rf $(patsubst %.cpp,%$(EXE),$(CPPTESTS)) + rm -rf $(patsubst %.cpp,%$(EXE),$(CPPTESTS)) -# buildbin: # build a binary version +# buildbin: # build a binary version # R CMD INSTALL -l ./tmp --build rstan @@ -108,7 +108,7 @@ STANLIB := $(STAN_MATH_SUBMODULE)lib STANSRC := $(STAN_SUBMODULE)src GTEST_LIB := tests/cpp/gtest.o -LIBPTHREAD = +LIBPTHREAD = ifeq ($(OS), linux) LIBPTHREAD += -lpthread endif @@ -116,11 +116,11 @@ endif CPPTESTS := $(shell find tests/cpp -type f -name '*_test.cpp') STAN_INSTANTIATION_FILES := $(patsubst rstan/src/%.cpp,tests/cpp/lib/src/%.o,$(wildcard rstan/src/*.cpp)) -## comment this out if you need a different version of R, +## comment this out if you need a different version of R, ## and set set R_HOME accordingly as an environment variable R_HOME := $(shell R RHOME) -## include headers and libraries for R +## include headers and libraries for R RCPPFLAGS := $(shell $(R_HOME)/bin/R CMD config --cppflags) RLDFLAGS := $(shell $(R_HOME)/bin/R CMD config --ldflags) @@ -132,6 +132,7 @@ RCPPLIBS := $(shell echo 'Rcpp:::LdFlags()' | $(R_HOME)/bin/R --slave) RINSIDEINCL := $(shell $(R_HOME)/bin/Rscript -e 'RInside:::CxxFlags()') RINSIDELIBS := $(shell $(R_HOME)/bin/Rscript -e 'RInside:::LdFlags()') + $(STAN_INSTANTIATION_FILES) $(GTEST_LIB) $(patsubst %.cpp,%$(EXE),$(CPPTESTS)) $(patsubst %.cpp,%.d,$(CPPTESTS)): CPPFLAGS += -Wall $(RCPPFLAGS) $(STAN_INSTANTIATION_FILES) $(GTEST_LIB) $(patsubst %.cpp,%$(EXE),$(CPPTESTS)) $(patsubst %.cpp,%.d,$(CPPTESTS)): CXXFLAGS += $(RCPPFLAGS) $(RCPPINCL) $(RINSIDEINCL) $(shell $(R_HOME)/bin/R CMD config CXXFLAGS) -I $(STAN_MATH_SUBMODULE) -I rstan/inst/include -I $(STANSRC) $(addprefix -isystem ,$(wildcard $(STANLIB)/*)) -isystem $(GTESTPATH)/include -isystem $(GTESTPATH) -ftemplate-depth=256 -DBOOST_RESULT_OF_USE_TR1 -DBOOST_NO_DECLTYPE -DBOOST_DISABLE_ASSERTS @@ -173,3 +174,4 @@ endif ifneq (,$(filter $(patsubst %.cpp,%$(EXE),$(CPPTESTS)),$(MAKECMDGOALS))) -include $(patsubst %$(EXE),%.d,$(filter $(patsubst %.cpp,%$(EXE),$(CPPTESTS)),$(MAKECMDGOALS))) endif + diff --git a/rstan/rstan/R/git_rstan_head.R b/rstan/rstan/R/git_rstan_head.R index 8a4e23cb0..eb486edd1 100644 --- a/rstan/rstan/R/git_rstan_head.R +++ b/rstan/rstan/R/git_rstan_head.R @@ -29,32 +29,3 @@ git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" git_head <- function() "2e1f913d3ca3678128f159b3d17d3d1f9b82704e" -git_head <- function() "50a10cca0de66ee51cf3832d1aa6ef29b3590659" -git_head <- function() "57444a10c088dbc52105fb394aadb996cd0c9df0" -git_head <- function() "2165604262c9085f0e47408d60bb04c76e5485c3" -git_head <- function() "2165604262c9085f0e47408d60bb04c76e5485c3" -git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" -git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" -git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" -git_head <- function() "be52c5b610c9e8a7c31819c86213cc94e89c2b94" -git_head <- function() "232f85fb57547338339acfd5e950f4de64d7c9d9" -git_head <- function() "232f85fb57547338339acfd5e950f4de64d7c9d9" -git_head <- function() "232f85fb57547338339acfd5e950f4de64d7c9d9" -git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" -git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" -git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" -git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" -git_head <- function() "2d7c20b4be33936c67b88215245a96f4255b3748" -git_head <- function() "e56060c723e942519335c25d4b41741ec9296af2" -git_head <- function() "9a2bbd719e4eabb86453ffaa8d61712a2c727b7c" -git_head <- function() "9a2bbd719e4eabb86453ffaa8d61712a2c727b7c" -git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" -git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" -git_head <- function() "84a8ac9462e22eefa1422620d35e4309d0eefaf5" -git_head <- function() "6aeef54d87f66376dffcac152b392d671bb4b44c" -git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" -git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" -git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" -git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" -git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1" -git_head <- function() "217c5b6f42ff3216540f9aad1082e7d78e59f0b1"