Using their notation of algorithm 1, I was wondering whether there is a method / object that allows me to access the fitted (log) VI approximation density for the samples, that is q(\theta_s)? Or did you guys manually write a function for the log of the, say, mean field VI density at fitted parameters?
Another related question: Is the PSIS diagnostic applicable for variational approximations when working with right censored data and a corresponding likelihood?
It would be very handy to have code for this in R, too, given a stanfit-class object that gets returned by vb. If I wanted to code this, do I have all the ingredients available/exposed through rstan (e.g. the mean and standard deviation of each Gaussian factor)?