Initialization failed: nondecision time must be nonnegative

Dear community,

in my R data simulation for a drift diffusion model using cmdStanr and the wiener_full_lpdf function, I receive negative initial values for my nondecision time. Also, for some iterations the random variable is smaller than the nondecision time.

Here is the relevant code:


data {                                  
  int<lower=1> M; // Number of participants
  ...
  int<lower=1> N; // Tials per participant
  ...
}

parameters {
  // Population level
  ...
  real<lower=0> mu_t0;
  real<lower=0> sigma_t0;
  real<lower=0> mu_xi;
  real<lower=0> sigma_xi;
  ...
  
  // Individual level
  ...
  vector<lower=0>[M] t0_j;
  vector<lower=0>[M] xi_j;
  ...
  
  // Trial level
  matrix<lower=0>[M, N] t0_ij;
}

model {
  // Population level
  ...
  mu_t0 ~ normal(0,.1);
  sigma_t0 ~ normal(0,.1);
  mu_xi ~ normal(0,.1);
  sigma_xi ~ normal(0,.1);
  ...

  // Individual level
  ...
  t0_j ~ normal(mu_t0,sigma_t0);
  xi_j ~ normal(mu_xi,sigma_xi);
  ...
  
  // Trial level
  for (j in 1:M) {
    for (i in 1:N) {
      t0_ij[j, i] ~ normal(t0_j[j],xi_j[j]);
    }
  }
  
  // Likelihood
  for (j in 1:M) {
    for (i in 1:N) {
      if (response[j, i] == 1) {
        target += wiener_full_lpdf(rt[j, i] | a_ij[j, i], t0_ij[j, i], 0.5, v_ij[j, i], sv_j[j], 0, st0_j[j]);
      } else {
        target += wiener_full_lpdf(rt[j, i] | a_ij[j, i], t0_ij[j, i], 0.5, -v_ij[j, i], sv_j[j], 0, st0_j[j]);
      }
    }
  }
}

Here is an excerpt from the output:

Chain 1 Rejecting initial value:
Chain 1   Error evaluating the log probability at the initial value.
Chain 1 Exception: wiener_full_lpdf: Nondecision time is -0.144704, but must be nonnegative! (in ...)
...
Chain 1 Exception: wiener_full_lpdf: Random variable  = 1.02586, but must be greater than nondecision time = 1.14669 (in ...)
...
Chain 1  Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.
Chain 1 Initialization failed.
Warning: Chain 1 finished unexpectedly!

I use Windows 11 and package version 0.6.1

Looks like the nondecision time is t0_ij[j, i]. It’s declared with <lower=0> in the parameters so I don’t know how it could be negative.
The “random variable” problems means rt is sometimes smaller than nondecision time. The code you have posted does not show where rt comes from but I assume it’s in the data section. If it’s a M-by-N matrix, you can use it as an upper bound for the nondecision time parameter and that should fix the problem.

  matrix<lower=0,upper=rt>[M, N] t0_ij;

Thanks for your reply and suggestion!

Indeed, t0_ij[j, i] is the nondecision time.

My data block looks as follows:

data {                                  
  int<lower=1> M; // Number of participants
  int<lower=1, upper=6> K; // Lotteries per participant
  int<lower=1> N; // Trials per participant
  array[M, N] int<lower=1, upper=K> option1; // Option 1 in each trial for each participant
  array[M, N] int<lower=1, upper=K> option2; // Option 2 in each trial for each participant
  array[M, N] int<lower=0, upper=1> feedback; // Feedback for each trial for each participant
  array[M, N] int<lower=0, upper=1> response; // Response for each trial for each participant
  array[M, N] real rt; // Response time for each trial for each participant
}

I tried to constrain t0_ij[j, i] with the upper bound that you suggested but still receive values greater than rt (as well as nonnegative values).

That data block declares rt as an 2d array, which is similar to a matrix but technically a different type. Stan is quite strict when it comes to types. A model like

data {
  int<lower=1> M; // Number of participants
  int<lower=1> N; // Trials per participant
  array[M, N] real rt; // Response time for each trial for each participant
}
parameters {
  matrix<lower=0,upper=rt>[M,N] t0_ij;
}

would fail to compile with an error along the lines of

Upper bound must be a scalar or of type matrix. Instead found type array[,] real.

and you’d need to explicitly convert the type

  matrix<lower=0,upper=to_matrix(rt)>[M,N] t0_ij;

Values outside the constraint range should never happen.
Maybe it’s not running the model you think but a previous model without the constraint, I don’t know. Which interface (RStan, PyStan, CmdStan, …) are you using?