From c8594b32358510da6d4d3a27f02e7bbf655a86d8 Mon Sep 17 00:00:00 2001 From: jgabry Date: Fri, 24 Jul 2026 11:16:42 -0600 Subject: [PATCH 1/4] =?UTF-8?q?Add=20plot=C2=A0of=20vector=20parameter=20i?= =?UTF-8?q?n=20posterior=20vignette?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit closes #443 --- vignettes/posterior.Rmd | 70 +++++++++++++++++++++++++++++++++++------ 1 file changed, 61 insertions(+), 9 deletions(-) diff --git a/vignettes/posterior.Rmd b/vignettes/posterior.Rmd index 1d056dbc6..f0d81c1fb 100644 --- a/vignettes/posterior.Rmd +++ b/vignettes/posterior.Rmd @@ -1,5 +1,5 @@ --- -title: "Working with Posteriors" +title: "Working with posteriors" output: rmarkdown::html_vignette: toc: true @@ -15,12 +15,20 @@ vignette: > ```{r child="children/_settings-knitr.Rmd"} ``` +```{r setup, message=FALSE, warning=FALSE} +library(cmdstanr) +library(posterior) +library(ggplot2) +``` + ## Summary statistics -We can easily customize the summary statistics reported by `$summary()` and `$print()`. +We can easily customize the summary statistics reported by `$summary()` and +`$print()`. ```{r} -fit <- cmdstanr::cmdstanr_example("schools", method = "sample") +# this model produces convergence warnings that we'll ignore +fit <- cmdstanr_example("schools", method = "sample") fit$summary() ``` @@ -90,7 +98,8 @@ fit$summary(variables = NULL, "Strictly Positive" = strict_pos) # fit$print(variables = NULL, "Strictly Positive" = strict_pos) ``` -For more information, see `posterior::summarise_draws()`, which is called by `$summary()`. +For more information, see `posterior::summarise_draws()`, which is called +internally by `$summary()`. ## Extracting posterior draws/samples @@ -125,12 +134,13 @@ one array per variable. Setting `with_chains = FALSE` combines the chains, giving the same general structure as the list returned by `rstan::extract()`: ```{r variable-arrays} -draw_arrays <- posterior::extract_list_of_variable_arrays( +draw_arrays <- extract_list_of_variable_arrays( fit$draws(), variables = c("mu", "theta"), with_chains = FALSE ) str(draw_arrays) +dim(draw_arrays$theta) ``` The first dimension of each array indexes draws, and any remaining dimensions @@ -141,7 +151,7 @@ variable as a multidimensional random variable, with its posterior draws handled behind the scenes: ```{r structured-draws} -draws_rvars <- posterior::as_draws_rvars( +draws_rvars <- as_draws_rvars( fit$draws(c("mu", "theta")) ) theta_rvar <- draws_rvars$theta @@ -152,13 +162,13 @@ theta_difference <- theta_rvar[1] - theta_rvar[2] theta_difference hist( - posterior::draws_of(theta_difference), + draws_of(theta_difference), main = "Difference between theta[1] and theta[2]", xlab = "theta[1] - theta[2]" ) -# Direct access to the underlying draws is also available -theta_array <- posterior::draws_of(theta_rvar) +# Direct access to the underlying draws is also available with posterior::draws_of +theta_array <- draws_of(theta_rvar) dim(theta_array) ``` @@ -166,3 +176,45 @@ The object `theta_rvar` behaves like the vector declared in the Stan program. `theta_array` provides direct access to its underlying draws, with the first dimension indexing draws. See the [`rvar` vignette](https://mc-stan.org/posterior/articles/rvar.html) for details. + + +### Plotting the draws of a vector + +Because `theta_array` has draws in the first dimension and the vector index (the +eight schools) in the second, we can reshape it into a long data frame and +overlay the individual draws. + +```{r vector-draws-plot} +theta_plot <- draw_arrays$theta + +theta_df <- data.frame( + .draw = rep(seq_len(nrow(theta_plot)), times = ncol(theta_plot)), + school = rep(seq_len(ncol(theta_plot)), each = nrow(theta_plot)), + theta = c(theta_plot) +) + +ggplot(theta_df, aes(school, theta, group = .draw)) + + ggplot2::geom_line(alpha = 0.01) +``` + +The reshaping above uses only base R. Tidyverse users can produce the same plot +directly from the draws data frame (`format = "df"`) with +`tidyr::pivot_longer()`, extracting the vector index from variable names like +`theta[1]`: + +```{r vector-draws-plot-tidy, eval=FALSE} +fit$draws("theta", format = "df") |> + tidyr::pivot_longer( + cols = dplyr::starts_with("theta"), + names_to = "school", + names_transform = readr::parse_number, + values_to = "theta" + ) |> + ggplot(aes(school, theta, group = .draw)) + + geom_line(alpha = 0.01) +``` + +Here `school` is simply the index into the `theta` vector. In many models the +vector index corresponds to a meaningful covariate, for example the time points +of a time series. In that case you can replace `school` with the associated +covariate values to plot each draw as a function of that covariate. From 45e5c993c696c676f99a544259962ecd00512e49 Mon Sep 17 00:00:00 2001 From: jgabry Date: Fri, 24 Jul 2026 11:23:15 -0600 Subject: [PATCH 2/4] Update posterior.Rmd --- vignettes/posterior.Rmd | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/vignettes/posterior.Rmd b/vignettes/posterior.Rmd index f0d81c1fb..d0f3c07f5 100644 --- a/vignettes/posterior.Rmd +++ b/vignettes/posterior.Rmd @@ -7,7 +7,7 @@ output: params: EVAL: !r identical(Sys.getenv("NOT_CRAN"), "true") vignette: > - %\VignetteIndexEntry{Working with Posteriors} + %\VignetteIndexEntry{Working with posteriors} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- From 0d5473587ca25e4e8c9b0989ba3383e2bc1e5ac5 Mon Sep 17 00:00:00 2001 From: jgabry Date: Fri, 24 Jul 2026 11:29:21 -0600 Subject: [PATCH 3/4] use non-centered parameterization --- vignettes/posterior.Rmd | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/vignettes/posterior.Rmd b/vignettes/posterior.Rmd index d0f3c07f5..fc0f27fe8 100644 --- a/vignettes/posterior.Rmd +++ b/vignettes/posterior.Rmd @@ -27,8 +27,7 @@ We can easily customize the summary statistics reported by `$summary()` and `$print()`. ```{r} -# this model produces convergence warnings that we'll ignore -fit <- cmdstanr_example("schools", method = "sample") +fit <- cmdstanr_example("schools_ncp", method = "sample") fit$summary() ``` @@ -39,7 +38,7 @@ posterior::default_summary_measures() To change the variables summarized, use the `variables` argument: ```{r} -fit$summary(variables = c("mu", "tau")) +fit$summary(variables = c("mu", "tau", "theta")) ``` We can also change which functions are used: @@ -78,7 +77,7 @@ fit$summary(c("mu", "tau"), quantile, .args = list(probs = c(0.025, .05, .95, .9 Each summary function is applied separately to each variable and receives a matrix whose rows are saved iterations and whose columns are chains. ```{r} -fit$summary(variables = NULL, dim, colMeans) +fit$summary(variables = "theta", dim, colMeans) ``` @@ -94,7 +93,7 @@ to the requested number of digits. ```{r} strict_pos <- function(x) if (all(x > 0)) "yes" else "no" -fit$summary(variables = NULL, "Strictly Positive" = strict_pos) +fit$summary(variables = c("mu", "tau", "theta"), "Strictly Positive" = strict_pos) # fit$print(variables = NULL, "Strictly Positive" = strict_pos) ``` From a878117f2664aa14059fd841eadb22d486807385 Mon Sep 17 00:00:00 2001 From: jgabry Date: Tue, 4 Aug 2026 11:27:23 -0600 Subject: [PATCH 4/4] Use default bayesplot theme with white background --- vignettes/posterior.Rmd | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/vignettes/posterior.Rmd b/vignettes/posterior.Rmd index fc0f27fe8..0120ba3ee 100644 --- a/vignettes/posterior.Rmd +++ b/vignettes/posterior.Rmd @@ -19,6 +19,7 @@ vignette: > library(cmdstanr) library(posterior) library(ggplot2) +theme_set(bayesplot::theme_default()) ``` ## Summary statistics @@ -193,7 +194,7 @@ theta_df <- data.frame( ) ggplot(theta_df, aes(school, theta, group = .draw)) + - ggplot2::geom_line(alpha = 0.01) + geom_line(alpha = 0.01) ``` The reshaping above uses only base R. Tidyverse users can produce the same plot