nl-DDM: a non-linear drift-diffusion 1 model accounting for the dynamics 2 of single-trial perceptual decisions 3

The Drift-Diffusion Model (DDM) is widely accepted for two-alternative forced-choice 12 decision paradigms thanks to its simple formalism, straightforward interpretation, and close ﬁt to 13 behavioral and neurophysiological data. However, this formalism presents strong limitations to 14 capture inter-trial dependency and dynamics at the single-trial level. We propose a novel model, 15 the non-linear Drift-Diffusion Model (nl-DDM), that addresses these issues by allowing the 16 existence of several trajectories to the decision boundary. We show that the ﬁtting accuracy of 17 our model is comparable to the accuracy of the DDM, with the non-linear model performing 18 better than the drift-diffusion model for an equivalent complexity. To give better intuition on the 19 meaning of nl-DDM parameters, we compare the DDM and the nl-DDM through correlation 20 analysis. This paper provides evidence of the functioning of our model as an extension of the 21 DDM. Our model paves the way toward more accurately analyzing single-trial dynamics for 22 perceptual decisions and accounts for pre-and post-stimulus inﬂuences. 23


Introduction 25
Perceptual decision-making has been studied extensively from behavioral (Ratcliff and McKoon,26 2008; Ratcliff and Smith, 2004), neurophysiological (Gold and Shadlen, 2001), and computational 27 (Gold and Shadlen, 2007) perspectives, as it is omnipresent in daily activities. When decisions are 28 timed, evidence accumulation models describe human and animal behavior well. They assume 29 that decisions are made when enough sensory evidence from the external world has been gath- The decision state is represented by a decision variable traveling from a starting point (for example, drawn from a uniform distribution, centered around 0 and of width 2 . It is represented as "SP" on the figure) to a boundary ("Correct boundary" or "Incorrect boundary") under the influence of a drift. Here, the drift depends on the current state of the decision. Depending on the position of 0 relative to , the drift will hence have different shapes. The trajectory is also impacted by white noise so that real trajectories are similar to the thin blue lines. From the stimulus onset, the decision process is delayed by a certain non-decision time ( ). Over an ensemble of decisions, probability density functions of correct and error response times can be created, as displayed here.
In order to provide an intuition for the other parameters, we consider first the potential function 138 derived from the drift term ( Figure 2). It is a function ( ) defined from a drift ( ) such that: In our case, we therefore have: . ( The decision variable can be seen as a ball traveling along the potential function.  From Figure 2, we can see that there are two potential sinks at and − , as well as a source at , 142 which derive directly from the topology of the system. Therefore, ± are the decision boundaries in the wrong direction is easier to correct because a small perturbation in the other direction can 151 counterbalance that effect. This is not as much the case when the wells are deep because then the 152 decision variable is driven rapidly to the stable fixed point, making perturbations less reversible. 153 We can also observe the impact of on the potential function in Figure 2C. Similar to the DDM, 154 fitting of response times can be obtained by solving the Fokker-Planck equation corresponding to 155 the Langevin equation defined above (Shinn et al., 2020b). Then, a non-decision time comes 156 into play in order to shift the resulting distribution to account for biological transmission delays.

157
To better understand the parameters of our model in comparison to the DDM, it can be useful 158 to define a mean drift rate across all trajectories. Since the deterministic trajectories only approach 159 the decision boundary asymptotically, we define an estimate of the mean drift rate. Considering 160 that the maximum drift for each trajectory causes the largest variation in decision value, we can 161 approximate the mean drift of each trajectory by its maximum, and subsequently average over all 162 the trajectories to get an estimate of the mean drift. Put in equations, we obtain: The noise term does not intervene as we assumed a Gaussian white noise. We observe a disconti-164 nuity in , due to the presence of an unstable fixed point at that location. Trajectories determined 165 by 0 = will finish in either well under the influence of noise, and the mean of the noise being zero, the two scenarios are equally likely. Consequently, the mean drift for these trajectories is the 167 average between and , with (respectively ) is the maximum negative (respectively 168 positive) drift rate achievable by the system. The graph of the max drift as a function of starting 169 point is given in Figure 3.  This equivalence is coherent with the interpretation of and as the impact of the stimulus 182 on the decision, and shows that in the absence of a stimulus, the two models follow the same 183 behavior. Because the nl-DDM is not assuming reflection symmetry, the presence of a stimulus 184 impacts the trajectories generated by the two models in different ways.

185
Behavioral results

186
For decision-making analysis, it is helpful to obtain each participant's response times and decision 187 accuracy, particularly for decision model fitting. whether a word presented on screen existed or not. The second one is a dataset not presented be-192 fore, in which participants were shown visual stimuli on screen and had to classify them according 193 to their type (either "face" or "number").

194
To ensure the correctness of both datasets in terms of behavioral measurements, we describe of stimuli is not significant in terms of accuracy (Table 1) nor in terms of response times (Table 2).

208
In the "face is left button" stimulus-response mapping, where participants were instructed 209 to click left upon face stimulus presentation and right when they were presented with a num-  (4), which depends both on and . The darker line represents the variation of the mean drift thus defined as a function of , while the pale blue curve is the variation of the mean drift as a function of . Since is strictly positive, we also represented the absolute value of the mean drift (dotted line). That allows for comparing the magnitude difference of the mean drift rate when or vary. We see that varying changes the mean drift rate more strongly than similar variations of at a given value of .
98.48 ± 1.12% ( = 15), whereas participants who underwent the "face is right button" stimulusresponse mapping, participants ( = 10) responded within 541 ± 30 ms and an accuracy of 98.77 ± 213 0.60%. The effect of the stimulus-response mapping on accuracy and response time was not sig-214 nificant (Tables 3 and 4). We do note however a marginal interaction effect between stimulus-215 response mapping and stimulus type on the accuracy of participants ( = 0.052, Table 1).

216
These results show the uniformity of participant responses across mappings and stimuli. All 217 participants, mappings and stimuli were considered together in the subsequent analyses.  speed condition (12%), and 2∕25 in our dataset (8%)), so 81% of all fitted models were kept.

230
Comparison of loss values 231 The first metric we used to compare the models is the loss value after fitting. Fitting is done by min-232 imizing the negative log-likelihood, which gives information on how close the curve of theoretical 233 response times is to empirical response times histograms. For a measure that takes into considera-234 tion the number of parameters and samples, we also computed the Bayesian Information Criterion 235 (BIC). All the test results on fitting performance are summarized in Table 5.   Speed-accuracy trade-off 247 We computed the behavior prediction of each model type to ensure that the results are consistent 248 with empirical observations. For that, we used a metric described by Roitman and Shadlen (2002), 249 whereby the loss is computed as the sum of the mean squared error on mean response time and 250 the mean squared error on the predicted accuracy over all conditions. We observe that there is no   Table 5, and Figure 7).

254
Although we fit the parameters separately for each stimulus type, we merge all the results to build 255 relations between the parameters of the DDM and the parameters of the nl-DDM. Of the resulting 256 fitted models, we rejected the participants that were rejected in our previous analysis (participants 257 6 and 11). In addition, participant 22 was rejected due to a fitted boundary outside of the other 258 models' range. Hence, 44 models were taken into account.

259
Given the mathematical formalism described above, we expect to find a negative correlation boundaries ± . In other words, increasing will decrease the drift, hence the negative correlation.

270
The correlation matrix of the nl-DDM and DDM parameters across all models is given in Figure   271 8. Note that we only took into account our dataset for parameter comparison, and while fitting the 272 models to this set, we assumed that the starting points followed a uniform distribution spanning creases, the other should decrease. In our data, since the accuracy is similar across all participants 282 and conditions, and the noise term is kept constant, these two terms are strongly correlated. Note 283 that the effect of each parameter is still different, as shown on Figure 2B and C. While increasing 284 deepens both wells, increasing will not only deepen the wells but also pull them apart. Effectively, 285 the relation between and is not linear, as seen on Figure 9. 286 We also observed the known relations within DDM parameters in the correlation matrix (

293
Concerning the cross-model comparison, we observe that is negatively correlated to the DDM 294 boundary and positively to 0 . The link with the decision boundary is expected, as in the nl-

295
DDM regulates the depth of the decision wells, that is, the time necessary to reach each decision.

296
More specifically, increasing makes the wells more attractive and hence results in fast decisions.

297
Conversely, increasing the DDM boundary will result in longer response times as more time will be 298 needed for the decision variable to reach the boundary, given the linear drift. 299 We also observe a significant negative correlation between the parameter of the nl-DDM and associated with decision bias 0 , the greater the amount of information needed to reach a decision. In this case, PC3 would reflect this decision caution: the greater the uncertainty associated 329 with decision bias 0 , the more cautious one must be, and in turn greater amounts of evidence one 330 needs to process. We have presented in this paper a non-linear model of decision-making. This model is a form of However, it must be noted that these parameters are not entirely equivalent as we did not find a 350 perfect mapping between them, meaning that the nl-DDM is conceptually different from the DDM. populations could produce the dynamics we have described. However, the main assumption of 356 the reduction was that the network was invariant through reflection. We argue that the mecha-357 nisms described by the nl-DDM are in fact similar to these of the DWM, but offer a broader range 358 of application beyond the case of symmetrical models.

359
A question that arises from our analyses is the different assumptions made on the starting point. 360 We took in both cases a uniform distribution, but while fitting our dataset, we assumed that this to show that with fewer degrees of freedom, our model could fit behavioral data better than the in fast and confident observations, although little evidence has been accumulated (we would be located at a plateau in our model), that is, even if the stimulus was not well perceived.

416
The dynamics that we propose here are not the sole product of mathematical formalism and 417 constraints, but have deep roots in empirical observations made in neurophysiological studies.

418
More specifically, three phases can be identified in the decision trajectories: an initial inertia stage, such that [ 0 − , 0 + ] ⊆] − , [ (Laming, 1968), or from other parametric distributions (Ratcliff and Rouder, 1998). We will consider uniformly distributed starting points in our fitting to provide Figure 11. Timeline of a single trial. Each trial is preceded by a rest period, followed by a baseline period (necessary for EEG processing, not reported here), each lasting 1.5 seconds. A noise clip consisting of 9 random-dot frames of 100 ms each indicates the arrival of the stimulus in a non-stimulus-specific fashion. The stimulus then appears on a noisy visual background for 100 ms. The same noisy background frame then lasts until the participant's response and times out after 2 seconds otherwise.
Pre-existing dataset from Wagenmakers et al. (2008) 556 To discard the possibility of better performances emerging from the fitting algorithm or data acqui-557 sition, we also lead our analyses on a pre-existing dataset taken from Wagenmakers et al. (2008).

558
17 human participants performed a classification task, as they were randomly presented with real 559 or invented words. The invented words were generated from real words by changing a vowel, and were discarded to avoid anticipatory responses. More details can be found in Wagenmakers et al.

567
(2008), and the dataset can be accessed from here.

568
Behavioral analyses 569 We are interested in comparing model parameters between the DDM and the nl-DDM. It is impor-570 tant to check whether participants' performance across stimulus-response mappings and stimuli is 571 coherent in terms of response times and accuracy. Indeed, the experimental paradigm we defined 572 entails two types of stimuli and two motor commands for the choices. In addition, we have cre-573 ated two experimental groups, which were instructed to respond with opposite motor commands.

574
First, we computed the percentage of stimuli in each class to verify that the stimuli were globally 575 equiprobable for each participant. Since we designed the experiment to display each stimulus with 576 the same probability at each trial, we expect this number to be close to 50%. Otherwise, participants 577 could opt for a strategy that prioritizes one response against the other. Then, we performed two 578 mixed-model ANOVAs, testing response times and accuracy respectively. The stimulus-response 579 mapping was considered a between-subject factor and the stimulus type a within-subject factor. The classical way of fitting evidence-accumulation models is by fitting one drift for each stimulus cat-582 egory separately. In that case, the positive and negative boundaries still correspond to correct and 583 incorrect responses respectively, and the starting points are taken from the same distribution re-584 gardless of the stimulus. Consequently, one pair of boundaries ± , the middle of the starting point 585 distribution 0 and its half-width , and two drifts 0 and 1 (corresponding respectively to "face" 586 and "number+sound" trials) have to be fitted in the DDM. Similarly, one pair of stable fixed points 587 (attractors, also corresponding to the decision boundaries) ± , one time scale and two unstable 588 fixed points 0 and 1 (repellers, that will tune the drift in the "face" and "number+sound" stimuli 589 respectively) are needed for the nl-DDM. In both cases we fix the noise parameter to = 0.3. As explained by Ratcliff (1978), since the speed-accuracy trade-off is determined by the boundary sep-591 aration, fitting two parameters among drift, boundary, and noise is constraining enough. Hence,  Hence, we compare each loss pairwise, using three repeated-measure one-sided paired-sample 641 -tests. Indeed, we want to test whether the nl-DDM is better than the DDM with these three met-642 rics, hence testing the hypothesis loss nl-DDM < loss DDM . Since we are comparing 3 losses, we set the 643 threshold for significance to = 0.017, corresponding to the Bonferroni-corrected 5% threshold. First, we computed the correlation matrix between all the parameters of both models. This allows for a first look into first-order interactions between model parameters, within and across model types. The correlation coefficients were computed using Pearson's , defined as: , = cov( , ) .
Next, to compare parameters of the DDM to parameters of the nl-DDM more quantitatively, 655 we performed principal component analysis on the correlation matrix of DDM and nl-DDM fitted 656 parameters. The goal is indeed to find how parameters relate to each other. This becomes possible 657 by observing the coefficients of the decomposition matrix.