Advice on slab width choice for regularized horseshoe

Hi all,

I’m using some code that was given by Michael Betancourt for the regularized horseshoe. The code runs fine, but I was wondering if there were any rules-of-thumb or general ways of thinking through how to choose the values in the \hl{transformed parameters} block.

Thanks,

David

modelString = "
data {
  int <lower=1> n;          // number of observations
  int <lower=1> p;          // number of predictors
  real readscore[n];        // outcome
  matrix[n,p] X;            // inputs

}


transformed data {
  real p0 = 5;
  real slab_scale = 2.5;
  real slab_scale2 = square(slab_scale);
  real slab_df = 4;
  real half_slab_df = 0.5 * slab_df;
  
}

parameters {
  vector[p] beta_tilde;
  vector<lower=0>[p] lambda;
  real<lower=0> c2_tilde;
  real<lower=0> tau_tilde;
  real alpha;
  real<lower=0> sigma;

}


transformed parameters {
  vector [p] beta ; // regression coefficients
  real tau0 = (p0 / (p - p0)) * (sigma / sqrt(1.0 * n));
  real tau = tau0 * tau_tilde; // tau ~ cauchy(0, tau0)
    
    // c2 ~ inv_gamma(half_slab_df, half_slab_df * slab_scale2)
    // Implies that marginally beta ~ student_t(slab_df, 0, slab_scale)
    
    real c2 = slab_scale2 * c2_tilde;
    
    vector[p] lambda_tilde =
      sqrt( c2 * square(lambda) ./ (c2 + square(tau) * square(lambda)) );
    
    // beta ~ normal(0, tau * lambda_tilde)
    beta = tau * lambda_tilde .* beta_tilde;
    
    }
 
model {
  beta_tilde ~ normal(0, 1);
  lambda ~ cauchy(0, 1);
  tau_tilde ~cauchy(0,1);
  c2_tilde ~ inv_gamma(half_slab_df, half_slab_df);
  
  alpha ~ normal(0, 2);
  sigma ~ cauchy(0, 1);
  
  readscore ~ normal(X * beta + alpha, sigma);

}

// For posterior predictive checking and loo cross-validation
generated quantities {
  vector[n] readscore_rep;
  vector[n] log_lik;
  for (i in 1:n) {
  readscore_rep[i] = normal_rng(alpha + X[i,:] * beta, sigma);
  log_lik[i] = normal_lpdf(readscore[i] | alpha + X[i,:] * beta, sigma);
}
  }


"

Slab corresponds to prior for big coefficients. So you can think of what would be the prior if you would have only the big coefficients in the model. It’s specifically useful e.g. for logistic regression with (almost) separable data.