Conflict between prior and data

I use Stan only occasionally so it is possible that what follows is a really trivial question. Apologies.

I am using three priors based on distributions with finite support, such as the beta distribution. Draws from each prior distribution are transformed into the parameters of a custom distribution that I have written.

The problem is that some combinations of draws from the three priors yield parameters for the custom distribution that are inconsistent with the data and an error is thrown.

What I want to happen in such cases is that the log likelihood is simply increased by zero and we move on to the next set of draws.

In a bit more detail I have three parameters:

parameters {
  real<lower = 0, upper = 1> q1;
  real<lower = 0, upper = 1>  q2;
  real<lower = 0, upper = 1>  q3; 
}
...
 model {
   // priors
   q1 ~ beta(shape11, shape21); 
   q2 ~ beta(shape12, shape22); 
   q3 ~ beta(shape13, shape23);  
 ...
}

where the shape11, ..., shape23 parameters are data.
The q1,...,q3 parameters are transformed to produce a parameter mu of my custom distribution the code for which includes the quantity

log1p(y - mu)

where y is the vector of observations.

The trouble is that depending on the observations some elements of 1 + y - mu can be negative which throws an error and terminates the run.
I think that what I need to know is how to amend my custom function so that when a negative element of 1 + y - mu occurs the value returned is 0.

In R I could do this using tryCatch() but I don’t know of an equivalent feature in Stan.

type or paste code here

Can we see the maths/code for how mu is computed? It may be possible to constrain the qs so that 1 + y - mu is always positive.

Here is the code:

transformed parameters{
  real qtilde1 = (h1- l1) * q1 + l1;
  real qtilde2 = (h2- l2) * q1 + l2;
  real qtilde3 = (h3- l3) * q1 + l3;
  real mu;
  real<lower = 0> nu = qtilde3 / qtilde2;
  real xi = log10(nu); 
  real<lower = 0> sigma;
  
  
  
  sigma  = qtilde2 * xi / (nu * (nu - 1));
  mu = qtilde1 - qtilde2 / nu;
  }

The supports for the qtilde1,..., qtilde3 are the ranges (l1, h1),... where the ranges are given as data. They represent physical quantities of which an expert is believed to have some intuition. They are then transformed to mu, sigma and xi which are the parameters of an extreme value distribution.

Hard constraints can be implemented in a Stan program using the reject function which prevents the any inconsistent configurations from being accepted as new Markov chain states. That said I would strongly advise against using this it as Stan programs with nontrivial reject functions require careful initializations and break the gradients on which Hamiltonian Monte Carlo relies. See for example the discussion in Section 2.2.4 of General Taylor Models.

@hhau gives the ideal approach. If you can figure out how to constraint the qs so that mu is always consistent then you won’t need the reject function and your sampling should go much more smoothly.