Help Extracting Fitted Values from CmdStanR

Hi everyone,

I’m a beginner and I need some help extracting fitted values from a CmdStanR model. Here is my model code:

data {
  int n; //number of observations in the data
  vector[n] mpg; //vector of length n for the car's MPG
  vector[n] weight_c; //vector of length n for the car's weight
  vector[n] cylinders_c; ////vector of length n for the car's cylinders
  vector[n] hp_c; //vector of length n for the car's horsepower
}

parameters {
  real alpha; //the intercept parameter
  real beta_w; //slope parameter for weight
  real beta_cyl; //slope parameter for cylinder
  real beta_hp; //slope parameter for horsepower
  real<lower = 0> sigma; //variance parameter and restrict it to positive values
}

model {
  //linear predictor mu
  vector[n] mu;
  
  //write the linear equation
  mu = alpha + beta_w * weight_c + beta_cyl * cylinders_c + beta_hp * hp_c;
  
  //likelihood function
  mpg ~ normal(mu, sigma);
}

generated quantities {
  //replications for the posterior predictive distribution
  real y_rep[n] = normal_rng(alpha + beta_w * weight_c + beta_cyl * 
  cylinders_c + beta_hp * hp_c, sigma);
}

With Rstan, I would just run:

y.rep <- extract(linear.fit.3)[["y_rep"]]

Then I can use bayesplot to plot the posterior predictive densities:

ppc_dens_overlay(y = cars.data$mpg, yrep = y.rep[1:100, ])

Where I’m stuck is that I’m not sure what the equivalent of extract is for CmdStanR.

I appreciate any help and your patience!

One way to do it would create a stanfit object from the cmdstan csv output files:

stanfit <- rstan::read_stan_csv(linear.fit.3$output_files())

then you could proceed as you would using RStan.

Another way would be to use the posterior package to get a data frame that you can play with:

library(posterior)

# an array of draws from posterior
draws_array <- linear.fit.3$draws() 

# convert to data frame 
draws_df <- as_draws_df(draws_array) #  or as_draws_matrix() for matrix

The posterior package also can be used in conjunction with bayesplot:

mcmc_hist(linear.fit.3$draws("y_rep"))

See the getting started with cmdstanr documentation for more.

2 Likes