How to run ADVI from R

Interfaces
,
#1

Hi,
I am fairly new to Bayesian estimation and am looking to use ADVI to estimate my model in hopes of better performance.

This is my model

data {
  int<lower=0> N;//Number of observations
  int<lower=1> J;//Number of predictors with random slope
  int<lower=1> K;//Number of predictors with non-random slope
  int<lower=1> L;//Number of customers/groups
  int<lower=0,upper=1> y[N];//Binary response variable
  int<lower=1,upper=L> ll[N];//Number of observations in groups
  matrix[N,K] x1;
  matrix[N,J] x2;
}
parameters {
  vector[J] rbeta_mu; //mean of distribution of beta parameters
  vector<lower=0>[J] rbeta_sigma; //variance of distribution of beta parameters
  vector[J] beta_raw[L]; //group-specific parameters beta
  vector[K] beta;
}
transformed parameters {
  vector[J] rbeta[L];
  for (l in 1:L)
    rbeta[l] = rbeta_mu + rbeta_sigma .* beta_raw[l]; // coefficients on x
}
model {
  rbeta_mu ~ normal(0,5);
  rbeta_sigma ~ gamma(1,1);
  beta~normal(0,5);
  for (l in 1:L)
    beta_raw[l] ~ std_normal();

  for(n in 1:N)
    y[n]~bernoulli_logit(x1[n] * beta + x2[n] * rbeta[ll[n]]);
}

I am unfamiliar on how the interface works from R for advi. Any suggestions or resources are appreciated.

Also if the model code needs to be modified please tell me how…

Thank you.

#2

If you can sample your model with the stan command in rstan, you should just be able to use the vb command in rstan to do variational Bayes on your model. Whether it works well or not is another question, but computationally it should just run.

#3

From R there are two ways to run advi, one using the RStan interface and one using the new CmdStanR interface.

If you’re using the RStan interface:

# https://mc-stan.org/rstan/reference/stanmodel-method-vb.html
model <- stan_model("file.stan") # replace "fiie.stan" with whatever your stan file is called
fit <- vb(model, ...) # replace ... with other arguments, e.g. data

If you’re using the CmdStanR interface:

# https://mc-stan.org/cmdstanr/reference/model-method-variational.html
model <- cmdstan_model("file.stan") # replace "fiie.stan" with whatever your stan file is called
fit <- model$variational(...) # replace ... with other arguments, e.g. data
#4

When it works well ADVI is much faster than MCMC so you will indeed get better performance. However, there are many many cases when ADVI doesn’t work well and so any increase in speed doesn’t really matter.

But there may also be things you can do to your Stan program to improve performance without having to change the algorithm. For example, see if you can vectorize any calls to probability distributions, e.g., when you use bernoulli_logit(). That can speed up a model quite a bit.

#5

I did try to vectorize the operations such tha the model block was

model {
  vector[N] p;
  rbeta_mu ~ normal(0,3);
  rbeta_sigma ~ gamma(1,1);
  beta~normal(0,2);
  for (l in 1:L)
    beta_raw[l] ~ std_normal();

  p = x1 * beta + (x2 .* rbeta[ll]) * ones; 
  y~bernoulli_logit(p);
}

However this did not provide any noiticeable speed up.