Estimated non-decision times larger than reaction times with wiener function in RStan

Hello Stan community,

I am trying to fit a drift-diffusion model (DDM) on some reaction time (RT) data. Participants had to choose between a 20 cents (lower boundary) and 60 cents (upper boundary) bet.
During data preprocessing, I removed all trials faster than 150 ms (a frequently used lower-bound for RTs in the literature).
I used a hierarchical model with the basic Wiener function:

data {
  int <lower=1> nSubs;
  int <lower=1> maxTrials;
  int <lower=1> subTrials[nSubs];
  real <lower=-9> betsize[nSubs, maxTrials];
  real <lower=-9> betRT[nSubs, maxTrials];
}


parameters {
  // Group level
  //// Non-decision time
  real <lower=0.01, upper=2> mu_tau;
  real <lower=0.01, upper=2> sd_tau;
  //// Boundary separation
  real <lower=0.1> mu_alpha;
  real <lower=0.01, upper=5> sd_alpha;
  //// Response bias
  real <lower=-4, upper=4> mu_beta;
  real <lower=0> sd_beta;
  //// Drift rate
  real mu_delta;
  real <lower=0> sd_delta;

  // Subject level
  //// Non-decision time
  vector <lower=0.01, upper=2> [nSubs] tau;
  //// Boundary separation
  vector <lower=0.1> [nSubs] alpha;
  //// Response bias
  vector <lower=-4, upper=4> [nSubs] beta;
  //// Drift rate
  vector [nSubs] delta;
}


model {
  // Define values
  real betaval[nSubs];
  
  // Sample parameter values
  //// Group level
  mu_tau ~ uniform(0.01, 2);
  sd_tau ~ uniform(0.01, 2);
  mu_alpha ~ uniform(0.1, 5);
  sd_alpha ~ uniform(0.01, 5);
  mu_beta ~ uniform(-4, 4);
  sd_beta ~ uniform(0, 1);
  mu_delta ~ normal(0, 10);
  sd_delta ~ uniform(0, 10);
  
  //// Subject level
  for (subi in 1:nSubs) {
    tau[subi] ~ normal(mu_tau, sd_tau)T[0.01, ];
    alpha[subi] ~ normal(mu_alpha, sd_alpha)T[0.1, ];
    beta[subi] ~ normal(mu_beta, sd_beta);
    delta[subi] ~ normal(mu_delta, sd_delta);
  }
  
  
  for (subi in 1:nSubs) {
    // Parameter values
    betaval[subi] = Phi(beta[subi]);
    
    for (triali in 1:subTrials[subi]) {
      if (betsize[subi, triali] > 0) {
        // Sample RT
        if (betsize[subi, triali] == 60) {
          betRT[subi, triali] ~ wiener(
            alpha[subi], tau[subi], betaval[subi], delta[subi]
          );
        } else {
          betRT[subi, triali] ~ wiener(
            alpha[subi], tau[subi], 1 - betaval[subi], -delta[subi]
          );
        }
      }
    }
  }
}


generated quantities {
  // Define values
  real betaval[nSubs];
  real log_lik[nSubs];
  
  for (subi in 1:nSubs) {
    betaval[subi] = Phi(beta[subi]);
    log_lik[subi] = 0;
    
    for (triali in 1:subTrials[subi]) {
      if (betsize[subi, triali] > 0) {
        // Log likelihood
        if (betsize[subi, triali] == 60) {
          log_lik[subi] = log_lik[subi] + wiener_lpdf(
            betRT[subi, triali] | alpha[subi], tau[subi], betaval[subi],
            delta[subi]
          );
        } else {
          log_lik[subi] = log_lik[subi] + wiener_lpdf(
            betRT[subi, triali] | alpha[subi], tau[subi], 1 - betaval[subi],
            -delta[subi]
          );
        }
      }
    }
  }
}

When I fit this model and do some posterior checks, I find that several participants have RTs faster than the mean estimated non-decision time (tau) at the subject level. This does not make sense. It should not be possible to make a decision faster than the non-decision time.

I have tried many different boundaries for the priors of tau (e.g., lower=0.001, upper=1). The only time the model estimated tau < min(RT) is when I kept removed RTs < 250 ms. However, this is not great because then I am removing a lot of data and I have no theoretical justification for a cutoff of 250 ms.

I have several questions about this and I hope you can help me.

  • Why does this happen?
  • Why is Stan not “Rejecting initial value”?
    • It should be unable to estimate an RT faster than tau
  • How can I fix this?

Any help is appreciated.
Thank you in advance.

Hi, and welcome to the Stan Discourse.

It’s not my field, but I think there are at least two different ways of seeing this: either it’s a well-accepted assumption of people in the field, or it’s actually impossible. In the former case, you may or may not want to challenge that assumption, in the latter, there’s still the possibility that the non-decision time may vary among subgroups or individuals, so what you are estimating is an average that may end up on average lower than the value you compute.

In addition to the two above, a third, simpler, and less philosophical possibility is that due to the variance in the data, there’s nothing statistically preventing the decision to be less than the non-decision time, it’s simply what is most compatible with the data you have, regardless of the reason why it shouldn’t be. If it’s not a big deal, you can chalk it up to noise, but if it’s something unacceptable to your peers in an analysis you may need change the model to constrain decision times to be larger than non-decision times (e.g. formulate it as non-decision plus additional time, or constrain the variable to a hard-coded maximum and minimum value so the ranges don’t overlap – I don’t think you can use a dynamic constraint in Stan).

Stan will reject the initial value if the log posterior is zero, but it will not reject it if it’s incompatible with an assumption that is not statistically explicit, even if the estimate is incompatible with reality (I always use the example of a demographic parameter I tried to estimate in a model that was related to human lifespan with results that implied average lifetime of 120 years in some estimates, biologically it was wrong, but statistically it was ok)

As mentioned above, it depends on how big a problem it is. If it’s just annoying, but it’s just a nuisance parameter in your model, you don’t have to fix it, just explain what it means. If the relationship between the parameters is important, than you need to formulate the model in a way that precludes the incompatibility.

2 Likes

Maybe the new added full DDM function helps the problem with the initial value of non-decision time t0. But still, we need to make sure the t0 is smaller than the minRT. One easy way to do it is adding another link function just like what you did on the starting point bias. something like: t0 = Phi(tau[subi]) * minRT

Thank you for your answer!
In my field it is unfortunately not really possible that reaction times are slower than the time it takes to make a reaction (without thinking about which option to select, which is the interpretation of non-decision time). Similar to your average lifetime of 120 years, I could not just leave it.

This worked perfectly, I feel like I should have thought about it. Thank you!