Divergence warnings, decreasing adapt_delta doesn't work

Modeling
,
1

I am trying to fit the following model using rstan. I keep getting divergence warnings but decrease adapt_delta doesn’t fix the problem. Is there a way to reparameterize my code to avoid divergence warnings? Or is there an error in my code that is causing these divergence warnings?

data {
    int<lower=0> k; 
    real M[k]; 
    real var[k];  
}
    
parameters {
    real c; 
    real Z0; 
    real Zm; 
    real<lower=0> sigma;
}

model {
    c ~ uniform(0, 1); 
    Z0 ~ normal(5, 5);
    Zm ~ normal(35, 5);
    sigma ~ exponential(0.1);

    // Likelihood
    for(i in 1:k){
        M[i] ~ normal(c*var[i]*(var[i]-Z0)*sqrt(Zm-var[i]), sigma);
    }
}
2

Can you post the pairs plots?

3

It looks awful!

pairs plot
pairs plot697×551 92.2 KB

1 Like
4

If it helps, here is my R code to execute the stan code (using rstan).

init_fun <- function() { 
  list(c = runif(1, 0, 1),
       Z0 = rnorm(1, 5, 2), 
       Zm = rnorm(1, 35, 2),
       sigma = rexp(1, 0.1))
}  

model.stan = stan(file = "stancode.stan", 
                 data = model_data, 
                 iter = 6000, chains = 6,
                 save_dso = F,
                 init = init_fun,
                 control = list(adapt_delta = 0.999999,
                                max_treedepth = 15))
1 Like
5

Hm, maybe try using the offset/multiplier to account for the somewhat-far-from-inits that the priors imply for Z0 & Zm? I’ve also put a multiplier on sigma, below, but given the exponential prior I don’t think that’s likely to make a huge difference.

Oh!! I wonder if the inconsistent declaration/prior for c might be causing trouble! if you have something that has a restricted range like that, it’s best to tell the sampler in the declaration rather than relying on the prior.

Here’s the version with all my suggestions, but try the updated c declaration first by itself as I’m more confident in the sensibility of that as a fix here.

...
parameters {
    real<lower=0,upper=1> c; 
    real<offset=5,multiplier=5> Z0; 
    real<offset=35,multiplier=5> Zm; 
    real<lower=0,multiplier=.1> sigma;
}
model {
    //c ~ uniform(0, 1); //implied by the declared limits, so can omit but leaving here as commented-out for clarity
    Z0 ~ normal(5, 5);
    Zm ~ normal(35, 5);
    sigma ~ exponential(0.1);

    // Likelihood
    for(i in 1:k){
        M[i] ~ normal(c*var[i]*(var[i]-Z0)*sqrt(Zm-var[i]), sigma);
    }
}
1 Like
6

Thank you!

Updating the declaration for c helped a bit but I still have divergence warnings (see first plot below). I also tried to use and offset and multiplier for Z0 and Zm, but I still have divergence warnings and the plots look at bit worse (see second plot below). Should I update my initial values? If so, do you have recommendations on the best initial values?

pairs plot only c adjusted
pairs plot only c adjusted697×551 100 KB

only c adjusted ^

pairs plot Zm, Z0, c adjusted
pairs plot Zm, Z0, c adjusted697×551 81 KB

c, Zm and Z0 updated ^

7

My sense is that this model has some serious problems with nonidentification. If we do

transformed parameters{
  real zm_new[k] = sqrt(Zm - var);
  real z0_new[k] = var - Z0;
 }

Then the likelihood contains the term c*z0_new*zm_new, where all three of the multiplicands don’t show up anywhere else in the likelihood. Except by putting extremely informative priors on at least two of these quantities, I think there’s no way forward here.

1 Like
8

c wasn’t the only parameter missing constraints. There’s term sqrt(Zm-var[i]) which requires that Zm>var[i]. Try this constraint

parameters {
    real<lower=max(var)> Zm; 
}

(incidentally, your model does not compile for me because Identifier 'var' clashes with reserved keyword.)

1 Like