rethinking::dashboard() doesn't support CmdStanMCMC_CSV

[milesalanmoore]

rethinking::dashboard() doesn't seem to support CmdStanMCMC_CSV objects read into the environment via cmdstanr::as_cmdstan_fit(). I've really enjoyed using the dashboard function (thanks for sharing this!) on my own models. I recognize this function is not documented or actively supported for all object types. Regardless, below I provide a reprex and my own local work around (see next comment) for consideration. This is a common aspect of my work flow where I run threaded models on a server and analyze the samples locally.

I hope that this at least helps someone else who uses CmdStanR output CSVs! Thank you for your hard work.

Reprex

library(rethinking)
library(cmdstanr)

# fit a simple model with cmdstanr
code <- "
data {
  int<lower=0> N;
  vector[N] y;
}
parameters {
  real mu;
  real<lower=0> sigma;
}
model {
  mu ~ normal(0, 10);
  sigma ~ exponential(1);
  y ~ normal(mu, sigma);
}
"

stan_data <- list(
  N = 100,
  y = rnorm(100, mean = 5, sd = 2)
)

model <- cmdstan_model(write_stan_file(code))
model$sample(
  data = stan_data,
  chains = 4,
  output_dir = getwd()
) # not saving to env to demonstrate the following:

fit <- cmdstanr::as_cmdstan_fit(list.files(pattern = "-\\w{6}.csv"))
fit$diagnostic_summary()

# This works fine, i think
rethinking::precis(fit)

# These fail with CmdStanMCMC_CSV objects

rethinking::divergent(fit)
# This is summing the stepsize__ variable from the cmdstanMCMC_CSV object, not
# divergent transitions (in the sum(x[,,2]) call)

rethinking::dashboard(fit)
# This one has returned two different errors to me at different times:

# Error in attr(fit, "cstanfit")$sampler_diagnostics(inc_warmup = warmup) : 
#   $ operator is invalid for atomic vectors

# OR

# Error in trankplot(fit, pars = "lp__", lp = TRUE, add = TRUE) : 
#   requires map2stan, ulam or stanfit/cmdstan object

I have some modifications I've made locally to dashboard and trankplot that seem to prevent all of the above errors. I'll add those changes as a comment below. Fair warning, I am not really sure what the Error in attr(fit, "cstanfit")$sampler_diagnostics(inc_warmup = warmup) : $ operator is invalid for atomic vectors error is about... it seems to only occur for some objects (and I'm certain they are all CmdStanMCMC_CSV class). Therefore my work around might be nothing more than duct tape here...

I am happy to make a fork and submit as a PR if you'd like.

Cheers!

Session Info

> packageVersion("cmdstanr")
[1] ‘0.6.1.9000’
> cmdstanr::cmdstan_version()
[1] "2.37.0"


> sessionInfo()
R version 4.5.2 (2025-10-31)
Platform: x86_64-pc-linux-gnu
Running under: Ubuntu 22.04.5 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0 
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0  LAPACK version 3.10.0

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C               LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8    LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C             LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

time zone: Etc/UTC
tzcode source: system (glibc)

attached base packages:
[1] parallel  stats     graphics  grDevices datasets  utils     methods   base     

other attached packages:
 [1] rstan_2.32.6        StanHeaders_2.32.10 jsonlite_1.9.1      lubridate_1.9.4     forcats_1.0.0       stringr_1.5.0      
 [7] dplyr_1.1.4         purrr_1.0.4         readr_2.1.5         tidyr_1.3.1         tibble_3.3.0        tidyverse_2.0.0    
[13] ggplot2_3.5.1       rethinking_2.42     bayesplot_1.11.1    posterior_1.6.1     cmdstanr_0.6.1.9000

loaded via a namespace (and not attached):
 [1] gtable_0.3.6         shape_1.4.6.1        tensorA_0.36.2.1     xfun_0.51            QuickJSR_1.6.0      
 [6] processx_3.8.6       inline_0.3.21        lattice_0.22-6       bspm_0.5.7           tzdb_0.4.0          
[11] vctrs_0.6.5          tools_4.5.2          ps_1.9.0             generics_0.1.4       stats4_4.5.2        
[16] curl_6.2.1           pkgconfig_2.0.3      data.table_1.17.0    checkmate_2.3.3      RColorBrewer_1.1-3  
[21] distributional_0.5.0 RcppParallel_5.1.10  lifecycle_1.0.3      compiler_4.5.2       farver_2.1.2        
[26] codetools_0.2-20     pillar_1.11.1        MASS_7.3-65          abind_1.4-8          tidyselect_1.2.1    
[31] mvtnorm_1.3-3        stringi_1.8.4        reshape2_1.4.4       labeling_0.4.3       grid_4.5.2          
[36] colorspace_2.1-1     cli_3.6.5            magrittr_2.0.4       loo_2.8.0            pkgbuild_1.4.6      
[41] withr_3.0.2          scales_1.4.0         backports_1.5.0      timechange_0.3.0     matrixStats_1.5.0   
[46] gridExtra_2.3        hms_1.1.3            coda_0.19-4.1        evaluate_1.0.3       knitr_1.49          
[51] V8_6.0.1             viridisLite_0.4.1    rlang_1.1.6          Rcpp_1.0.14          glue_1.8.0          
[56] pkgload_1.4.0        rstudioapi_0.17.1    R6_2.6.1             plyr_1.8.8          

[milesalanmoore]

Quick patch for the divergent function

# extracts n_divergent from stan fit
divergent <- function( fit , warmup=FALSE ) {
  # Check if it's a cmdstan_fit object
  if ( inherits(fit, "CmdStanMCMC_CSV") ) {
    x <- fit$sampler_diagnostics(inc_warmup=warmup)
    return(sum(x[,,"divergent__"]))
  }
  
  x <- attr(fit,"cstanfit")$sampler_diagnostics(inc_warmup=warmup)
  sum(x[,,2])
}

For dashboard()

dashboard <- function( fit , warmup=FALSE , plot=TRUE , trank=TRUE ) {
  # <<<<< This section includes my changes 1 of 2 >>>>
  # Check if it's a cmdstan_fit object
  if ( inherits(fit, "CmdStanMCMC_CSV") ) {
    x <- fit$sampler_diagnostics(inc_warmup=warmup)
    n_chains <- dim(x)[2]
    y <- precis(fit, 3) # get n_eff and Rhat
    n_samples <- dim(x)[1] * dim(x)[2]
    
    # Adjust indexing for cmdstan format (uses named dimensions)
    divergent_idx <- "divergent__"
    treedepth_idx <- "treedepth__"
    energy_idx <- "energy__"
  } else {
    # Original code for other formats
    x <- attr(fit,"cstanfit")$sampler_diagnostics(inc_warmup=warmup)
    n_chains <- dim(x)[2]
    y <- precis(fit, 3)
    n_samples <- dim(x)[1] * dim(x)[2]
    
    # Use numeric indexing for other formats
    divergent_idx <- 2
    treedepth_idx <- 3
    energy_idx <- 6
  }
  # <<<<< End Changes 1 of 2  >>>>

  # display
  if ( plot==TRUE ) {
    set_nice_margins()
    par(cex.axis = 1 , cex.lab=1.2 )
    par(mfrow=c(2,2))
    # n_eff distribution
    Rhat_vals <- as.numeric( round( y[,'rhat'], 2 ) )
    plot( y[,'ess_bulk'] , (Rhat_vals) , xlab="number of effective samples" , ylab="Rhat" , ylim=c( 1 , max(1.1,(Rhat_vals),na.rm=TRUE) ) , log='y' )
    abline( v=0.1*n_samples , lty=1 , col="red" )
    abline( v=n_samples , lty=1 , col=grau() )
    abline( h=1 , lty=2 )
    
    # energy plot
    if ( inherits(fit, "CmdStanMCMC") ) {
      dens( x[,,energy_idx] , adj=0.1 , xlab="HMC energy" )
      mu <- mean(x[,,energy_idx])
      sig <- sd(x[,,energy_idx])
    } else {
      dens( x[,,3] , adj=0.1 , xlab="HMC energy" )
      mu <- mean(x[,,3])
      sig <- sd(x[,,3])
    }
    curve( dnorm( x , mu , sig ) , add=TRUE , col=rangi2 , lwd=1.5 )
    
    # count of divergent iters
    if ( inherits(fit, "CmdStanMCMC") ) {
      n_divergent <- sum(x[,,divergent_idx])
    } else {
      n_divergent <- sum(x[,,2])
    }
    
    plot( NULL , bty="n" , xlab=" " , ylab="" , xlim=c(0,1) , ylim=c(0,1) , xaxt="n" , yaxt="n" )
    text( 0.5 , 0.8 , n_divergent , cex=6 )
    text( 0.5 , 0.5 , "Divergent transitions" )
    if ( n_divergent==0 ) text( 0.5 , 0.2 , "Outlook good" , cex=1.5 )
    if ( n_divergent>0 & n_divergent<10 ) text( 0.5 , 0.2 , "Roll again" , cex=1.5 )
    if ( n_divergent>=10 ) text( 0.5 , 0.2 , "Check yourself before\nyou wreck yourself" , cex=1.5 )
    
    # three trace plots with lowest n_eff
  # <<<<< This section includes my changes 2 of 2 >>>>
    # for CmdStanMCMC objects, we need to attach a cstanfit attribute for trank/traceplot
    if ( inherits(fit, "CmdStanMCMC_CSV") ) {
      fit_with_attr <- fit
      attr(fit_with_attr, "cstanfit") <- fit
      if ( trank==TRUE ) {
        trankplot( fit_with_attr , pars="lp__" , lp=TRUE , add=TRUE )
      } else {
        traceplot( fit_with_attr , pars="lp__" , lp=TRUE , add=TRUE )
      }
    } else {
      if ( trank==TRUE ) {
        trankplot( fit , pars="lp__" , lp=TRUE , add=TRUE )
      } else {
        traceplot( fit , pars="lp__" , lp=TRUE , add=TRUE )
      }
    }
  # <<<<< End Changes 2 of 2  >>>>
  }
  # invisible result
  invisible(x)
}

One change to trankplot() is required to accept CmdStanMCMC_CSV's I think too. Namely I modified it to check if any of the labels returned by class(object) are in this list:

  if ( !(any(class(object) %in% c("map2stan","ulam","stanfit","CmdStanMCMC")) )) stop( "requires map2stan, ulam or stanfit/cmdstan object" )

instead of the current:

rethinking/R/trankplot.R

Line 14 in ac1b3b2

if ( !(class(object)[1] %in% c("map2stan","ulam","stanfit","CmdStanMCMC")) ) stop( "requires map2stan, ulam or stanfit/cmdstan object" )