Abstract

Impairments in decision-making are frequently observed in neurodegenerative diseases, but the mechanisms underlying such pathologies remain elusive. In this work, we study, on the basis of novel time-delayed neuronal population model, if the delay in self-inhibition terms can explain those impairments. Analysis of proposed system reveals that there can be up to three positive steady states, with the one having the lowest neuronal activity being always locally stable in nondelayed case. We show, however, that this steady state becomes unstable above a critical delay value for which, in certain parameter ranges, a subcritical Hopf bifurcation occurs. We then apply psychometric function to translate model-predicted ring rates into probabilities that a decision is being made. Using numerical simulations, we demonstrate that for small synaptic delays the decision-making process depends directly on the strength of supplied stimulus and the system correctly identifies to which population the stimulus was applied. However, for delays above the Hopf bifurcation threshold we observe complex impairments in the decision-making process; that is, increasing the strength of the stimulus may lead to the change in the neuronal decision into a wrong one. Furthermore, above critical delay threshold, the system exhibits ambiguity in the decision-making.

1. Introduction

Gamma-Aminobutyric Acid (GABA) is the most prevalent inhibitory neurotransmitter in the human brain [1, 2]. There is a body of evidence that aging has a major influence on the effectiveness of the GABAergic synapses [3]. Moreover, as found in the recent studies with the use of magnetic resonance spectroscopy (MRS, [4, 5]), local GABA concentrations in the frontal areas influence cognitive performance in aging adults. Taken together, aging can cause a reduced release of GABA to the intersynaptic cleft and decrease the quality of the synaptic transmission. This can result in an increase of the synaptic delays in the local inhibition.

On the other hand, in the aging process, the sensory capacity is declining, which affects the cognitive functions [68]. In particular, the working memory—which involves active manipulation of information—is affected, and this effect can further influence the perceptual decision-making [9]. (One important note here is that there is a difference between perceptual decision-making and making the abstract complex choices: the latter was reported not to be impaired in the elderly subjects [10], and some studies even report that elderly subjects are actually more efficient at making complex decisions [11, 12]. Therefore, in the paper we refer to the perceptual decision-making.).

In 1996, Salthouse [13] proposed a processing-speed theory of age-related deficits in cognition, for example, in working memory and perceptual decision-making. According to this theory, reduction in the processing speed can cause impairments in the cognitive functions for two major reasons. Cognitive functions can decline either because necessary logic operations cannot be executed within the time limit, or because the higher-order operations are blocked by slow execution of lower-order operations. This subject-matter was further experimentally investigated from different angles [1416], but the consensus in the field is that, in general, the elderly subjects are slower in perceptual decision-making than youngsters.

In this work, we propose a mechanism linking the two aspects of aging in cortical networks: the neurodegeneration in the local inhibitory synapses and the processing-speed related impairments in perceptual decision-making. This mechanism is based on a neuronal population model of decision-making based on a winner-take-all mechanism. The novelty lies in combining a winner-take-all mechanism well routed in the decision-making neuroscience, with the system of delayed differential equations representing the local inhibition within the two competing populations. With the use of this model, we are able to demonstrate that, for small synaptic delays in the local inhibition within the competing populations, the decision-making process depends directly on the strength of the stimulus, and the network is able to correctly identify the direction the stimulus came from. However, large delays can lead to a subcritical Hopf bifurcation resulting in complex decision-making process impairments. In particular, we demonstrate that, above the Hopf bifurcation point, increasing the strength of the stimulus can confuse the network and cause a wrong decision to be made. Furthermore, for delay values above critical threshold the system exhibits ambiguity in the decision-making. This effect can explain how the experimentally found difficulties in decision-making in elderly adults can be caused by loss in cognitive capacities [14].

The paper is organized in the following way. In Section 2, we introduce the model. In Section 3, qualitative analysis of the model, focusing on the stability analysis and Hopf bifurcation appearance, is made. Firstly, we present the analytic results for a symmetric model with no delay (Section 3.1) and then with positive delay (Section 3.2). In Section 4, we present the simulations undertaken to explore the dynamic repertoire of the model. In Section 5, we critically discuss the results and give recommendations for the future research.

2. Perceptual Decision-Making Model

The famous perceptual experiment on rhesus monkeys [17] by Shadlen and Newsome, involved a binary classification task: the monkeys had to assess whether the majority of the moving dots on the screen moved to the left or to the right. In the literature, one model proposed to model the decision-making in this experiment was the slow reverberation mechanism by Wang [18]. In Wang’s model, two populations of densely interconnected, spiking neurons compete with each other once being supplied by noisy inputs. This model involves simulations of two competing pools of stochastic spiking neurons. Since Wang’s work was published, some rodent [19] and computational [20] models were proposed to study perceptual choices.

In this work, we focus on modeling the most basic perceptual decision-making in neuronal networks. In such conditions, the network needs to disambiguate between two sensory stimuli. We consider changes in the firing rates of two positively interconnected self-inhibiting neuronal populations that receive external inputs and , respectively (Figure 1).

In order to describe the temporal dynamics of considered system (Figure 1), we introduce a system of delay differential equations (DDEs [21]) in the following form:where , denotes a model time scale, denotes the delay in self-inhibition, and is a coefficient describing the maximal capacity of a synapse, related to its anatomy. The function characterizes interactions between populations through synaptic plasticity which undergoes a Hebbian rule [22]; that is, the dynamics of synaptic weights is firing rate-dependent. For in-depth analytical and numerical investigations we have chosen a biologically plausible [23] sigmoid function .

In reality both populations considered in the proposed model are embedded in a larger network; therefore even in the absence of the population-specific sensory stimulus, they receive a constant input. Therefore, in the resting state, this system receives equal constant inputs to both nodes and, due to the system’s symmetry, both populations will be firing with the same rates. The symmetry breaks down if one of the populations receives an additional, external stimulation. The decision-making in this system means decoding which of the two populations received an additional stimulus (without estimation of the stimulus magnitude). The decoding is based on the difference between the two firing rates: the larger the difference , the more likely the decision that the population received the stimulation. The evidence behind each of the two options accumulates over time; therefore the psychometric function for the first population takes the integral form ofwhere is a parameter influencing the slope of the sigmoid function with respect to the cumulative difference . Similarly, we define the psychometric function for the second population asValues and can only asymptotically approach ; therefore we add a condition that if, at a given time point , exceeds the threshold value of (where is a given precision), the decision is made.

3. Qualitative Behavior of the Model

In this section, we provide the analytical results regarding the behavior of solutions of (1) for and constant symmetric input . Notice that the qualitative dynamics of (1) does not depend on .

3.1. Behavior of the Model for and Constant Symmetric Input

In this subsection, we present a detailed analysis of the model dynamics for and , as it is a crucial first step in the analysis of time-delayed models. While looking for steady states of (1), we need to solve the system of equationswhich yields . It is then obvious that both coordinates of any steady state are the same. Let us denote a steady state by . Clearly, , and the number of steady states depends on the shape of the graph of . Notice that the reference case is specific, as for the function is asymptotically linear , while for it tends to as .

Let us consider . Then we have , and therefore is increasing for , achieves its maximum at , and decreases to for . This implies that there are two steady states for and no steady state for , while for there is a bifurcation. The reference value , so there exist two steady states, and .

For , we have , and looking for zeros of we obtain , and we see that this quadratic equation has two positive solutions for , no real solution for , and one positive solution for . This means that for the function is increasing, for it has one maximum and tends either to (for ) or to (for , while for it is first increasing, then decreasing, and eventually increasing linearly to .

Corollary 1. For there is one steady state of (1); for there are two steady states for small values of and no steady state for larger values; for there is one steady state for small and sufficiently large values of , while for intermediate values 3 steady states exist.

Notice that due to its symmetric structure the system described by (1) always has solutions lying within a straight line . Clearly, assuming from both equations of (1) we obtainThe number of steady states determines the dynamics of (5). Let us assume that . Then there is only one steady state and the right-hand side of (5) is positive for and negative for . Therefore, is globally attractive. If , then there are two steady states for small , and the right-hand side is positive for and and negative for . This means that is locally stable, while is unstable. Moreover, solutions for tends to . If there is no steady state, then all solutions tend to , as they are increasing and unbounded. For there can be up to three steady states. Assume that there are three steady states . In this case and are stable, while is unstable. Moreover, for solutions tend to .

It is obvious that the symmetry of (1) implies that the phase space is divided into two symmetric subspaces by the straight line . Moreover, the dynamics of (5) is crucial in determining the whole model dynamics. Notice that and for we have . Similarly, for we have . Moreover, these straight lines are asymptotes for null-clines of (1). Therefore, asymptotic dynamics for large values of , is related to the dynamics of (5) for large , as it determines the direction of the vector-field in the region between the null-clines. Three generic types of the model dynamics when at least one steady state exists are presented in Figure 2.

3.2. Model Behavior for and Constant Symmetric Input

Now, we aim to study the influence of the delay on the dynamics of the model. Our main goal is to show that there exists such time delay for which the steady state loses stability and a Hopf bifurcation occurs. Moreover, we would like to analyze the type of this bifurcation.

Before starting the analysis, we first state some general results that will be useful in this section. Let us consider a general DDEwith a smooth function . Let denote a characteristic function for (6) at . Assume thatwhere and are polynomials, . Together with (6) we considerfor which is a characteristic function.

Lemma 2. Assume that and have no common imaginary root and has a pair of purely imaginary simple eigenvalues for some critical value , and moreover these eigenvalues satisfy transversality condition . If , then is a pair of simple eigenvalues of for which satisfies transversality condition. Moreover, the eigenvalues of (6) cross imaginary axis in the same direction as the eigenvalues of (8).

Proof. We only need to check transversality condition. The derivative for (6) is obtained from the relationAs while , this implies which is the relation determining transversality condition for (8).

Lemma 2 holds under weaker assumptions, that is, for and being analytic functions, not only polynomials; compare [24] for more details. However, for our analysis, it is sufficient to consider the presented version.

The next lemma is a simple consequence of Proposition  1 from [24].

Lemma 3. Let , . Then there exists a pair of purely imaginary eigenvalues , , for   for which eigenvalues cross imaginary axis from left to right. Moreover, a steady state loses stability for and cannot gain it again for .

Proof. Following [24] we define an auxiliary function . It is obvious that is a simple zero of . Looking for critical delay related to this pair of eigenvalues , we need to solve a system of equationsIt is obvious that we obtain a sequence of critical delays , . However, as the derivative of the auxiliary function is positive, eigenvalues always cross imaginary axis from left to right (independently of ), which means that a switch of stability appears only at .

Now, we turn to the main topic of this subsection that is analysis of (1) for . While studying a Hopf bifurcation, we follow the ideas introduced in [25]. Let us rewrite the system described by (1) in its functional form (cf. [26, 27]); that is,where and for , is a linear part, and is nonlinear part of the system. Hence,where , , and . With the linear operator we are able to associate an operator which is a solution to (1) for initial data . In this way we obtain a strongly continuous semigroup generated by an infinitesimal generator (cf. [26, 27]). If for some critical value the generator has an eigenvalue , then a Hopf bifurcation occurs under the assumptions that the eigenvalues are simple and cross imaginary axis with nonzero speed when crosses .

Let us denote a characteristic matrix by ; that is,Looking for eigenvalues, we need to find zeros of the characteristic functionIt is obvious thatClearly, we can use Lemmas 2 and 3 in the analysis of stability switches for the steady state . In the considered case, the steady state is stable for ; that is, has negative zeros, yielding . Hence, both quasi-polynomials and satisfy assumptions of Lemma 3. This means that there are two sequences of critical delays and associated with and , respectively. However, the switch of stability can occur only for or , depending on the magnitude of those delays. Clearly, according to Lemma 2 eigenvalues cross imaginary axis in the same direction for both and , which means that they cross from left to right, and therefore the steady state loses stability for the smallest critical delay.

As a result, we can state the following theorem.

Theorem 4. Let . The steady state of (1) is locally asymptotically stable for and unstable for , and at a Hopf bifurcation occurs.

Proof. We only need to check that . Notice that . Moreover, the function , is decreasing. As we obtain

For the reference values of parameters for we obtain and . Considering we obtain and . For we have and

Hence, we study the bifurcation at and corresponding . Let us denote to shorten the notation. Therefore, . Notice that from the form of we obtain .

In the following, we base on the ideas presented in [25]. We know that is a purely imaginary eigenvalue of the infinitesimal generator for if there exists a vector such that and then is an eigenvector for at . Moreover, is also an eigenvalue for the adjoint operator for and is an eigenvector, where . In the considered case, we are able to choose such that , where is the derivative with respect to the first variable, here .

When looking for we obtainand thereforeClearly, and we can choose . Hence, is the eigenvector for the eigenvalue .

Next, we need to find a vector that satisfieswhich means that coordinates of satisfy the same relation as for , that is, . Next, calculatingwe obtain , that is, , and eventually

The type of the studied bifurcation is determined by a coefficient of the third term in Taylor expansion of an orbit. This coefficient readswhere denotes a first derivative with respect to the second variable , while , andwhere , denotes th derivative with respect to the first variable, (in fact is a constant function as it does not depend on ), and

If , then the bifurcation is supercritical; that is, periodic solutions appear for and are stable in such a case. If , then the bifurcation is subcritical: periodic solutions exist for , and as the steady state is stable in this case, periodic orbits are necessarily unstable.

First, we calculate the denominator of . We haveand we easily check that the real part of this expression is equal to .

Next, we would like to calculate the numerator of . To this end, we need to calculate the derivatives of the nonlinear part necessary to calculate , and we omit bar in the notation . Denoting by , , and test functions from , we obtainEvaluating the second and third derivative at we obtainwhere , is denoted by , to shorten the notation.

Using the formula above, we calculateTo calculate we need to evaluate . We haveMoreover,and thereforeFinally,To obtain the last parameter we first need to calculate .

Let us denote and calculate the first term of this matrix; that is,Hence,Next we calculate , whereFinally, we evaluateand therefore

Now, we are in a position to check the sign of , where . Clearly,and thereforewhere and .

Notice that, for small , the first term of the expression above is dominating, and therefore in the range of parameters we are interested in the fact that, as , the first term is positive and of the order of units, suggesting that . This implies that the bifurcation we study is subcritical, which therefore yields instability of appearing periodic orbits. However, in general it is necessary to calculate the exact value of , as it is difficult to guess its sign.

Now, we calculate for the reference parameters with . Let us recall that , , , , , and moreover . We obtain and . Consecutive terms in the brackets are equal to , , and , while the last fraction is . Eventually,which means that, in our reference case, the bifurcation is subcritical.

In order to complement the analysis above, we can also show that in the reference case the second bifurcation appearing for is subcritical as well. We again denote by , to shorten the notation. Now, the pair , satisfies yielding , and from this relation we obtain the vector . Hence, is the eigenvector for the eigenvalue . Next, we find such that . As before, coordinates of satisfies the same relation as for , that is, , . Calculating the denominator of we obtainnextand we easily see that the real part of this expression is negative.

Next, we calculateEventually,and hence,For our reference values the expression in the brackets equals .

At the end we sum up the results of the Hopf bifurcation analysis in the following corollary.

Corollary 5. For reference parameter values the system described by (1) undergoes two subsequent subcritical Hopf bifurcations.

4. Numerical Simulations

In order to simulate the decision-making process, we assume a transient change in one of the inputs to the nodes; that is, we start from the resting state of the network and then solve (1) with andwhere is the stimulation duration time. In a properly working decision-making network, we should expect that the bigger the value of is, that is, the stimulus, the faster the decision in favor of population according to the psychometric function (2) will be.

In all the numerical experiments, we have chosen and (for this value the time scale allows to reproduce delays present in real systems). Note that the qualitative dynamics of the system described by (1) does not depend on . Moreover, we choose a reference value to reflect that there is a certain baseline input to the network. Other parameter values fixed across all simulations were , , and  s. We perform numerical simulations of the model with respect to the magnitude of delay and stimulus strength . For all of the simulations we use the built-in MATLAB delayed differential equation solver dde23 [28] with lowered default tolerances (RelTol and AbsTol equal to ). We performed the simulations on the Neuroscience Gateway platform [29].

One important aspect of the considered (1) is that we cannot guarantee the nonnegativity of solutions as the delayed self-inhibition term has a negative sign. It is obvious, however, that the firing rate cannot have a negative value and, in order for the model to be physiologically valid, we need to impose a barrier for the firing rates at zero. Thus, in the simulations, whenever one of the coordinates reaches value of zero, we stop the simulation and solve the reduced system with one of the coordinates equal to zero until its derivative becomes positive again; that is, its delayed value becomes lower than .

We start our numerical experiments from simulations in which stimulus of a fixed magnitude is applied to the network with different delay in the self-inhibition terms. As it can be expected from the analysis, for small delays, that is, below the Hopf bifurcation threshold, the solutions exhibit vanishing oscillations around the stable steady state; see left panel in Figure 3(a). The stimulus in this particular case is too small for any certain decision to be made, and both psychometric functions remain separated from the value of 1; see left panel in Figure 3(b). It is clear, however, that if the stimulus was stronger, then a certain decision would be made, that is, would cross the threshold of at some point.

Interestingly, for delays above the critical value at which first Hopf bifurcation occurs, we observe very nonintuitive behavior; see middle panels in Figures 3(a) and 3(b). Namely, the psychometric function value shows multiple perceptual switches; that is, there are multiple time points in which there is a change from into . For even larger delays we observe that a certain decision has been made for the same stimulus of magnitude that was insufficient to make a certain decision in the case of small delays; see right panels Figures 3(a) and 3(b).

Because of the observed nonintuitive behavior of the psychometric function, we decided to calculate a decision map; that is, to evaluate the psychometric function values on the model solutions for different values of delay and the stimulus strength , see Figure 4. As expected, for delay values below the critical Hopf bifurcation threshold, the certainty of the decision depends directly on the strength of the stimulus, and the networks is always able to correctly identify that the stimulus was applied to population . Interestingly, this is not the case for larger delays and we observe that for the same value of delay above the critical first bifurcation threshold the decision can be opposite depending on the stimulus strength, compare Figure 4(a). Moreover, above the bifurcation threshold we can have more than 10 transient decision switches before any certain decision is made; compare Figure 4(b).

In Figure 5, in order to complete presentation of the numerical experiments, we show the solution for which decision made by the network is in favor of the population instead of the population to which the stimulus was applied.

5. Discussion

In this study, we model perceptual decision-making in elderly individuals. Although there is an extensive evidence that the cognitive impairments in aging relate to the slow processing of information [13, 16, 30, 31], the cortical mechanisms underlying cognitive impairments in aging remain elusive. In principle, the decision-making experiments in humans are mostly performed in young and healthy individuals [32, 33], which is one of the reasons behind the lack of knowledge upon the cortical mechanisms of aging.

In order to bridge this gap, we propose to study the influence of synaptic delays on making perceptual choices with use of a population model of decision-making based on a winner-take-all mechanism. In this simple model, the network needs to make a binary choice between the two options. The inspiration for this model was a spiking neuronal network model with a slow reverberation mechanism by Wang [18], created to explain Shadlen and Newsome’s experiments on perceptual decision-making in rhesus monkeys [17]. Unlike Wang, we do not use simulations of stochastic spiking neuronal networks in our study though, but we take an analytic approach instead. This allows us to study the rich dynamic repertoire of networks with delay and the influence of the delay on the outcome decisions.

We achieve two main results that can contribute to the understanding of the associations between the delayed GABA-signaling and the cognitive impairments during aging. Firstly, the decision-making performance is dependent on the synaptic delays in the local inhibition. For the delays below the critical value, we observe a clear association between the magnitude of the stimulus and the probability of making the correct decision (Figure 4). This result is concordant with classic experiments on perceptual decision-making, as typically, a monotonic psychometric curve links the strength of the stimulus with the decision accuracy [17, 34]. However, for the synaptic delays exceeding the critical threshold, the system falls into a new dynamical regime and is no longer precise in making decisions. Dependent on both the synaptic delays and the signal magnitude, it can even achieve an accuracy of zero, by always choosing the wrong option. As known from experiments on reversal learning in rhesus monkeys [35] and humans [36], aging subjects exhibit habitual behaviors and are reluctant to relearn rules. Therefore in some cases, falling into a behavioral schema in which the wrong decision is being taken in repetitive fashion is possible in elderly subjects.

The second important result from our study is that, for delay values above critical threshold, the system exhibits ambiguity in decision-making, reflected by decision switches. This result can account for the slow reaction times in perceptual decision-making in the elderly [3741]; however the subsequent experimental validation is necessary in order to test this hypothesis.

Our model predicts that impairment in the local inhibition in the cortex can result in the impaired decision-making. Although GABA concentration in prefrontal cortex and perceptual decision-making are both affected by aging, there is a lack of computational models characterizing the causal link between the two. Therefore, the model should be validated in laboratory conditions. The prediction given by the model is hard to test in the human cohorts because recordings from interneurons in the cortex are invasive. However, there are now tools in translational psychiatry that make this validation viable. For instance, the decision-making quality can be evaluated in mice in multiple experimental paradigms [19] and that results can be than correlated with the speed of synaptic transmission evaluated postmortem in an in vitro experiment [42].

One remark to make is that our model is qualitative rather than quantitative, and the same mechanism can be encountered at different time scales of inhibition and for different configurations of inputs and stimuli. There are multiple GABAergic receptors in the cortex, and they have their own characteristic timescales. In example, the fast mode of inhibition is related to the GABA-A receptors [43, 44] which give synaptic delays lasting for several tens of milliseconds (excluding the afterdepolarizations lasting for several tens of milliseconds as well [4549]). On the other hand, the slowest mode of inhibition is related to the metabotropic GABA-B receptors which have a time scale of a few hundred milliseconds [50]. Still little is known about the structure and functions of these receptors [51]. In practice, the local inhibition in the nodes of the cortical network is most probably a combination of the multiple interacting processes at different time scales. Our model is a demonstration of the principle and does not require specification of GABAergic receptors that lead a primary role in the oscillatory mechanisms considered in this work.

As a summary, the model proposed in this work yields new insights into the mechanisms of aging in cortical circuits, mediated by neurodegeneration in the local inhibitory synapses.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Authors’ Contributions

The authors developed the model as a team. Urszula Foryś conducted the derivations in a collaboration with Martyna Płomecka and Katarzyna Piskała. Jan Poleszczuk performed and visualized numerical simulations. Natalia Z. Bielczyk formulated the introduction and the discussion part of the manuscript.

Acknowledgments

The authors would like to thank Maciej Borodzik, Piotr Kozakowski, and Franciszek Rakowski for the discussion and advice. The research leading to these results has received funding from the University of Warsaw, the European Community’s Seventh Framework Programme (FP7/2007–2013) under Grant Agreement no. 305697 (OPTIMISTIC), the European Community’s Seventh Framework Programme (FP7/2007–2013) under Grant Agreement no. 278948 (TACTICS), European Union’s Seventh Framework Programme for Research, Technological Development and Demonstration under Grant Agreement no. 603016 (MATRICS), and National Science Center in Poland under Grant Agreement no. 2011/03/N/ST1/00109.