I’m usually used to fitting stan models using rstanarm and getting outputs where the regression coefficients correspond to the variable/predictor name, but when using rstan and getting stanfit objects, I notice that that the betas are listed out by the index number, so what I was wondering is if there’s a way to be able to recover the original names of the coefficients, otherwise it might make it a bit difficult to parse through the outputs when there are a large number of predictors, especially for collaborators who are less familiar with the data.

This is a problem that stems from the fact that the interface you are working in (R in your case) has no idea what user’s Stan program is doing so it cannot assign names to output. The rstanarm, brms, and some other packages work around this because the package creators know what their Stan programs do and can put the right names in the right places from R. But you can do it yourself with something like

stanfit@sim$fnames_oi <- new_names

where new_names is a character vector whose size is equal to the number of posterior quantities (including the value of the log posterior kernel).

I wrote some stuff in ezStan that attempts to help re-apply names when contrast matrices are supplied as input to Stan models. Specifically, see ezStan::stan_summary. Here’s a lecture and code walking through an example hierarchical model, but the variable renaming inezStan::stan_summary should work for any model where you’re using contrast matrices.