Adding data and get_draws() (issue #18)

.Rbuildignore

@@ -20,7 +20,6 @@ benchmarks\.Rmd
 ^man-roxygen$
 ^data-raw$
 ^docker$
-doc
 .*\.Rcheck
 .*\.tar\.gz
 \.github

Makefile

@@ -9,6 +9,8 @@ $(PKGNAME)_$(VERSION).tar.gz: R/*.R inst/NEWS README.md
 	$(MAKE) clean && \
 		R CMD build --no-build-vignettes --no-manual . 
 
+build: $(PKGNAME)_$(VERSION).tar.gz
+
 inst/NEWS: NEWS.md
 	Rscript -e "rmarkdown::pandoc_convert('NEWS.md', 'plain', output='inst/NEWS')"&& \
 	head -n 80 inst/NEWS

NAMESPACE

@@ -1,6 +1,8 @@
 # Generated by roxygen2: do not edit by hand
 
+S3method("$",fmcmc_last_mcmc)
 S3method("[[",fmcmc_kernel_list)
+S3method("[[",fmcmc_last_mcmc)
 S3method(MCMC,default)
 S3method(MCMC,mcmc)
 S3method(MCMC,mcmc.list)
@@ -17,8 +19,6 @@ export(LAST_CONV_CHECK)
 export(LAST_MCMC)
 export(MCMC)
 export(MCMC_INFO)
-export(MCMC_with_conv_checker)
-export(MCMC_without_conv_checker)
 export(append_chains)
 export(check_initial)
 export(convergence_auto)
@@ -33,6 +33,7 @@ export(cov_recursive)
 export(get_)
 export(get_burnin)
 export(get_conv_checker)
+export(get_draws)
 export(get_elapsed)
 export(get_kernel)
 export(get_logpost)

NEWS.md

@@ -18,9 +18,15 @@
 * The new function `get_logpost()` returns the computed values of the objective
   function from the last `MCMC` run.
 
+* The new function `get_draws()` returns the MCMC draws from the kernel's 
+  proposal function.
+
 * The function `append_chains()` was randomly dropping one sample of the final
   set.
 
+* A new artificial dataset `lifeexpect` is shipped with the package. This simulates
+  1,000 observations of `age` at death using US's statistics.
+
 
 # fmcmc 0.4-0
 

R/append_chains.R

@@ -8,6 +8,36 @@
 #' @return If `mcmc.list`, an object of class `mcmc.list`, otherwise,
 #' an object of class `mcmc`.
 #' @export
+#' @examples 
+#' # Appending two chains
+#' data("lifeexpect")
+#' logpost <- function(p) {
+#'   sum(with(lifeexpect, dnorm(
+#'     age - p[1] - smoke * p[2] - female * p[3],
+#'     sd = p[4], log = TRUE)
+#'   ))
+#' }
+#' 
+#' # Using the RAM kernel
+#' kern <- kernel_ram(lb = c(-100, -100, -100, .00001))
+#' 
+#' init <- c(
+#'   avg_age = 70,
+#'   smoke   = 0,
+#'   female  = 0,
+#'   sd      = 1
+#' )
+#' 
+#' ans0 <- MCMC(initial = init, fun = logpost, nsteps = 1000, seed = 22, kernel = kern)
+#' ans1 <- MCMC(initial = ans0, fun = logpost, nsteps = 2000, seed = 55, kernel = kern)
+#' ans2 <- MCMC(initial = ans1, fun = logpost, nsteps = 2000, seed = 1155, kernel = kern)
+#' ans_tot <- append_chains(ans0, ans1, ans2)
+#' 
+#' # Looking at the posterior distributions (see ?lifeexpect for info about
+#' # the model) 
+#' plot(ans_tot)
+#' 
+#' 
 append_chains <- function(...) UseMethod("append_chains")
 
 #' @export

R/data.R

@@ -0,0 +1,37 @@
+#' Life expectancy in the US (2020)
+#' 
+#' A simulated data set based on statistics of life expectancy in the US at birth.
+#' 
+#' @details The data was generated using official statistics from the CDC and a
+#' study of life expectancy of smokers in the US published in the New England
+#' Journal of Medicine (see references).
+#' 
+#' According to the CDC, data from 2020 indicates that the average life expectancy
+#' of females in the US is 80.5 years vs 75.1 years for males (which declined
+#' with respect to 2019 after COVID hit the US). In Jha et al. (2013), evidence
+#' is presented indicating that individuals who smoke have at least ten years
+#' left of life expectancy compared to non-smokers.
+#' 
+#' The parameter estimates for the data generating process where:
+#' 
+#' - An average of life expectancy of 80.1.
+#' - Smokers live 10 years less than non-smokers.
+#' - Females live 5.4 years longer than males.
+#' 
+#' @format 
+#' A 1,000 rows data frame with three columns:
+#' 
+#' - `age`: Years lived.
+#' - `smoke`: 0/1 variable equal to 1 if the individual smokes.
+#' - `female`: 0/1 variable equal to 1 if the individual is a female at birth.
+#' 
+#' @references 
+#' Jha, P., Ramasundarahettige, C., Landsman, V., Rostron, B., Thun, M., Anderson, R. N.,
+#' ... Peto, R. (2013). 21st-Century Hazards of Smoking and Benefits of Cessation
+#' in the United States. New England Journal of Medicine, 368(4), 341–350.
+#' \doi{10.1056/NEJMsa1211128}
+#' 
+#' Arias, E., Tejada-Vera, B., & Ahmad, F. (2021). Provisional life expectancy
+#' estimates for January through June, 2020.
+#' \url{https://www.cdc.gov/nchs/data/vsrr/VSRR10-508.pdf}
+"lifeexpect"

R/kernel_mirror.R

@@ -34,11 +34,11 @@
 #' @references 
 #' Thawornwattana, Y., Dalquen, D., & Yang, Z. (2018). Designing Simple and
 #' Efficient Markov Chain Monte Carlo Proposal Kernels. Bayesian Analysis, 13(4),
-#' 1037–1063. \url{https://doi.org/10.1214/17-BA1084}
+#' 1037–1063. \doi{10.1214/17-BA1084}
 #' 
 #' Yang, Z., & Rodriguez, C. E. (2013). Searching for efficient Markov chain
 #' Monte Carlo proposal kernels. Proceedings of the National Academy of Sciences,
-#' 110(48), 19307–19312. \url{https://doi.org/10.1073/pnas.1311790110}
+#' 110(48), 19307–19312. \doi{10.1073/pnas.1311790110}
 #' @name kernel_mirror
 #' @examples 
 #' # Normal mirror kernel with 5 adaptations and 1000 steps of warmup (burnin)

R/kernel_ram.R

@@ -35,7 +35,7 @@
 #' @references 
 #' Vihola, M. (2012). Robust adaptive Metropolis algorithm with coerced acceptance
 #' rate. Statistics and Computing, 22(5), 997–1008.
-#' \url{https://doi.org/10.1007/s11222-011-9269-5}
+#' \doi{10.1007/s11222-011-9269-5}
 #' @family kernels
 #' @examples 
 #' # Setting the acceptance rate to 30 % and deferring the updates until

R/last.R

@@ -29,7 +29,7 @@ MCMC_INFO$clear <- function(nchains, env) {
   if (missing(env))
     MCMC_INFO$data. <- replicate(
       nchains, list2env({
-        list(logpost = numeric())
+        list(logpost = numeric(), draws = NULL)
       }
       ), simplify = FALSE)
   else {
@@ -70,6 +70,14 @@ MCMC_INFO$c_ <- function(x, val) {
 
 }
 
+#' Combine
+#' @noRd
+MCMC_INFO$rbind_ <- function(x, val) {
+
+  assign(x = x, value = rbind(MCMC_INFO$ptr[[x]], val), envir = MCMC_INFO$ptr)
+
+}
+
 #' @export
 print.fmcmc_info <- function(x, ...) {
 
@@ -176,6 +184,19 @@ get_logpost <- function() {
 
 }
 
+#' @export
+#' @details The function `get_draws()` retrieves the proposed samples from the
+#' kernel function.
+#' @rdname fmcmc-info
+get_draws <- function() {
+
+  if (get_nchains() == 1L)
+    return(get_("draws")[[1L]])
+  else
+    return(get_("draws"))
+
+}
+
 
 #' @export
 #' @rdname fmcmc-info
@@ -214,13 +235,24 @@ get_thin <- function() get_("seed")
 
 #' @export
 #' @rdname fmcmc-deprecated
-LAST_MCMC <- MCMC_INFO
-class(LAST_MCMC) <- c("fmcmc_last_mcmc", class(LAST_MCMC))
+LAST_MCMC <- structure(list(), class = "fmcmc_last_mcmc")
+
+#' @export
+`[[.fmcmc_last_mcmc` <- function(i, j, ..., exact=TRUE) {
+  .Deprecated("MCMC_INFO", old = "LAST_MCMC")
+  MCMC_INFO[[i]]
+}
+
+#' @export
+`$.fmcmc_last_mcmc` <- function(x, name) {
+  .Deprecated("MCMC_INFO", old = "LAST_MCMC")
+  MCMC_INFO[[name]]
+}
 
 #' Deprecated methods in fmcmc
 #' 
 #' These functions will no longer be included starting version 0.6-0. Instead,
-#' use the function sin [fmcmc-info].
+#' use the functions in [fmcmc-info].
 #' 
 #' @export
 #' @name fmcmc-deprecated