Error message with Betancourt specification of horseshoe prior

Hi all,

I’m trying to run an analysis using the horseshoe specification in Bayes Sparse Regression (2018). The full code is here

``````library(rstan)
library(loo)
library(bayesplot)
library(xtable)

n=nrow(PISA18sampleScale)
X = PISA18sampleScale[,2:11]
p=ncol(X)

modelString = "
data {
int<lower=1> n; // Number of data
int<lower=1> p; // Number of covariates
matrix[n, p] X;
}

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

model {
// tau ~ cauchy(0, sigma)
// beta ~ normal(0, tau * lambda)
vector[p] beta = beta_tilde .* lambda * sigma * tau_tilde;

beta_tilde ~ normal(0, 1);
lambda ~ cauchy(0, 1);
tau_tilde ~ cauchy(0, 1);

alpha ~ normal(0, 1);
sigma ~ cauchy(0, 1);

readscore ~ normal(X' * beta + alpha, sigma);
}
"
# Start estimation

nChains = 4
nIter= 10000
thinSteps = 10
warmupSteps = floor(nIter/2)

myfitHorseshoe = stan(data=data.list,model_code=modelString,chains=nChains,
iter=nIter,warmup=warmupSteps,thin=thinSteps)
``````

I’m getting the following error message

``````SAMPLING FOR MODEL '6a7711900e92e051ef7d7e1268a95faf' NOW (CHAIN 1).
Chain 1: Unrecoverable error evaluating the log probability at the initial value.
Chain 1: Exception: multiply: Columns of m1 (100) and Rows of m2 (10) must match in size  (in 'model98d60f2608e_6a7711900e92e051ef7d7e1268a95faf' at line 29)

[1] "Error in sampler\$call_sampler(args_list[[i]]) : "
[2] "  Exception: multiply: Columns of m1 (100) and Rows of m2 (10) must match in size  (in 'model98d60f2608e_6a7711900e92e051ef7d7e1268a95faf' at line 29)"
error occurred during calling the sampler; sampling not done
>
``````

In addition this error, I’m curious why the specification tau and beta are commented out.

I take it your code is this with `M` renamed to `p`. You have matrix `X` transposed, Betancourt defines `matrix[M,N] X` but your code says `matrix[n,p] X`.
You can leave that and simply remove the transpose operator `'` from the last line

``````  readscore ~ normal(X * beta + alpha, sigma); // not X'
``````

The model is parametrized with `tau_tilde` and `beta_tilde` instead of (more natural?) `tau` and `beta`. (I guess it’s like noncentered parametrization or something?) The model block defines the distributions for these parameter and the comments document the implied distributions for the derived parameters.
In other words, `tau = tau_tilde * sigma` so then `tau_tilde ~ cauchy(0, 1)` is an indirect way to declare `tau ~ cauchy(0, sigma)`.

Thanks very much. That was a silly mistake. It seems to be working but throwing divergent transition errors. I’m playing with the adapt_delta value and hopefully that will help.