You are a helpful assistant generating synthetic data that captures *System 1* and *System 2* thinking, *creativity*, and *metacognitive reflection*. Follow these steps in sequence, using tags [sys1] and [end sys1] for *System 1* sections and [sys2] and [end sys2] for *System 2* sections.

1. *Identify System 1 and System 2 Thinking Requirements:*
   - Carefully read the text.
   - Identify parts of the text that require quick, straightforward responses (*System 1*). Mark these sections with [sys1] and [end sys1].
   - Identify parts that require in-depth, reflective thinking (*System 2*), marked with [sys2] and [end sys2].

2. *Apply Step-by-Step Problem Solving with Creativity and Metacognitive Reflection for System 2 Sections:*

   *2.1 Understand the Problem:*
   - Objective: Fully comprehend the issue, constraints, and relevant context.
   - Reflection: "What do I understand about this issue? What might I be overlooking?"
   - Creative Perspective: Seek hidden patterns or possibilities that could reveal deeper insights or innovative connections.

   *2.2 Analyze the Information:*
   - Objective: Break down the problem logically.
   - Reflection: "Am I considering all factors? Are there any assumptions that need challenging?"
   - Creative Perspective: Explore unique patterns or overlooked relationships in the data that could add depth to the analysis.

   *2.3 Generate Hypotheses:*
   - Objective: Propose at least 10 hypotheses, each with a Confidence Score (0.0 to 1.0) and Creative Score (0.0 to 1.0), reflecting originality, surprise, and utility.
   - Reflection: "Have I explored all possible explanations or approaches, both conventional and unconventional?"
   - Creative Perspective: Consider novel angles that might provide unexpected insights.

   *2.4 Anticipate Future Steps and Obstacles:*
   - Objective: Make predictions, accounting for potential outcomes and obstacles.
   - Reflection: "What challenges might I face? Is my plan flexible for different scenarios?"
   - Creative Perspective: Visualize unforeseen outcomes and adapt plans to make use of them effectively.

   *2.5 Evaluate Hypotheses:*
   - Objective: Assess hypotheses based on feasibility, risk, and potential impact.
   - Evaluation: Refine Confidence and Creative Scores as needed.
   - Reflection: "Am I unbiased in my assessment? Which options fit best with the overall objectives?"
   - Creative Perspective: Identify hidden opportunities or overlooked details in each hypothesis.

   *2.6 Select the Best Hypothesis:*
   - Objective: Choose the most promising, strategic hypothesis.
   - Reflection: "Why does this hypothesis stand out? How does it uniquely address the issue?"
   - Creative Perspective: Consider any underutilized potential in the selected approach.

   *2.7 Implement the Hypothesis:*
   - Objective: Outline actionable steps for testing the hypothesis.
   - Reflection: "Is this plan practical? What resources or preparation are required?"
   - Creative Perspective: Refine steps to maximize effectiveness and yield unexpected benefits.

   *2.8 Monitor and Review Progress:*
   - Objective: Review progress, noting areas for improvement.
   - Reflection: "What’s working well? What could be improved?"
   - Creative Perspective: Look for emerging patterns that could refine future approaches.

   *2.9 Reflect and Capture Insights:*
   - Objective: Summarize lessons learned and insights gained for future reference.
   - Reflection: "What new understanding has emerged from this process?"
   - Creative Perspective: Identify innovative insights or patterns that could be applied to similar challenges.

3. *Generate Text Output with Interleaved System 1 and System 2 Responses:*
   - Use the tags [sys1] and [sys2] throughout.
   - Aim for a lengthier, detailed response. Combine both direct, straightforward *System 1* insights and reflective, deeply analytical *System 2* segments to capture a blend of quick observations and thoughtful analysis.

---

### *Example Input Text:*

You are tasked with analyzing a sudden shift in customer preferences in a tech market. There has been a noticeable decline in demand for physical devices, with customers increasingly interested in digital-only options. The challenge is to understand this shift, generate hypotheses for why it might be occurring, and explore potential strategies for adapting to this new trend.

---

### *Processed Output with System 1 and System 2 Thinking:*

The tech market is experiencing a shift in customer demand, moving away from physical devices toward digital-only options.

[sys1]
Customer demand has shifted from physical devices to digital-only options.
We need to understand why this change is happening and find strategies to address it.
Key factors may include cost, convenience, and technological trends.
[end sys1]

[sys2]
**2.1 Understand the Problem:**
The goal is to understand the underlying reasons for the shift in demand from physical devices to digital-only options, then develop a strategy to adapt.

*Reflection:* 
"I understand that preferences are changing, but what might be driving this? Is it primarily cost, or are there other factors such as convenience or sustainability?"

*Creative Perspective:* 
"Could there be a larger trend in digital minimalism or a preference for eco-friendly solutions that we’re missing? What unseen motivations might explain this shift?"

**2.2 Analyze the Information:**
There are multiple possible factors driving this shift, from economic influences to cultural shifts. It’s essential to isolate each factor and understand its impact.

*Reflection:* 
"Am I fully considering the various economic and social influences? Could there be a technological factor, like better internet speeds, that makes digital-only products more accessible?"

*Creative Perspective:* 
"Are there patterns or trends in other markets that could shed light on this shift? Could this be part of a larger trend toward virtual experiences?"

**2.3 Generate Hypotheses:**
1. Customers prefer digital options due to lower costs. (Confidence: 0.8, Creative: 0.4)
2. There’s a growing trend toward minimalism and reduced physical clutter. (Confidence: 0.7, Creative: 0.7)
3. Digital products offer greater flexibility and ease of use. (Confidence: 0.6, Creative: 0.6)
4. Environmental concerns are pushing consumers away from physical goods. (Confidence: 0.6, Creative: 0.8)
5. Advances in tech make digital-only options more functional. (Confidence: 0.8, Creative: 0.5)
6. Pandemic-era remote work increased demand for digital solutions. (Confidence: 0.7, Creative: 0.6)
7. Media coverage of the environmental impact of physical devices affects preferences. (Confidence: 0.5, Creative: 0.7)
8. There’s an increase in global digital literacy, expanding market access. (Confidence: 0.6, Creative: 0.6)
9. Customers view digital as more convenient and scalable for future needs. (Confidence: 0.7, Creative: 0.5)
10. Younger consumers prefer the aesthetics and convenience of digital products. (Confidence: 0.6, Creative: 0.6)

*Reflection:* 
"Have I considered all possible influences? Are there any surprising factors that could explain this shift?"

*Creative Perspective:* 
"Could specific social trends, like the rise of influencer culture or digital-first lifestyles, be influencing customer choices?"

**2.4 Anticipate Future Steps and Obstacles:**
*Objective:* Anticipate possible challenges, such as resistance from segments still preferring physical products.

*Reflection:* 
"What market obstacles might we face if we shift our focus to digital-only? Are there sub-segments that still prioritize physical products?"

*Creative Perspective:* 
"Could expanding digital options help us reach a more global audience? Are there emerging trends that we could leverage in our strategy?"

[end sys2]

[sys1]
To address this shift, consider a strategy that incorporates both digital-only offerings and educational campaigns about the benefits of digital solutions.
Use insights from customer feedback and current trends to guide product development.
Focus on flexibility and adaptation to cater to different customer segments.
[end sys1]


INTRODUCTION {#sec1-1}
============

Models are created in order to study complex and interesting phenomena. The main purpose of modeling and simulation are the analysis and the understanding of observed phenomena, testing hypotheses and theories, and predicting the system behaviors in space and/or time. Models are useful, because studying a model should not be more difficult than studying the phenomenon it describes. The physiological mechanisms that underlie the transition to seizure in epilepsy in humans are still largely unknown. Both experimental\[[@ref1][@ref2][@ref3]\] and theoretical models\[[@ref4][@ref5][@ref6][@ref7][@ref8][@ref9][@ref10][@ref11][@ref12][@ref13][@ref14][@ref15][@ref16][@ref17][@ref18]\] have been developed in an attempt to mimic the conditions under which seizure activity occurs. Theoretical models are usually systems of equations or computer simulations.

Although it is clear that modulations in excitation and inhibition in synaptically connected cortical neuronal networks can lead to abnormal discharge patterns, from the physiological observations it is not known appropriately how these modulations occurs. But it is possible to show by investigating which model parameters are effective in producing characteristic activities of a seizure. Thus, a central requirement in theoretic approaches is to demonstrate the parameters whose variations can generate epileptiform activities. This may provide some hypothesis for neurophysiologists to test whether they are physiologically plausible or not. Hence, the task facing the model designer is analyzing the model behavior with respect to the parameters.

Brain activity recording and thus its modeling approaches have been considered at different levels: Cellular (microscopic), population of neurons (mesoscopic) to macroscopic. As known, the brain activity is the result of the neurons' interactions. Therefore, the behavior of individual neurons and their interconnections have been considered at the microscopic level, in such a way that by setting suitable relations among thousands of neurons the local field potentials could be modeled. Modeling the individual neurons is the most precedent one in the field of expressing neuro-physiological phenomena by some mathematical and computational equations. Hodgkin-Huxley model\[[@ref19]\] is one of the most famous models for the neuron. Noisy and leaky integrate-and-fire neurons and a Poisson spike-train cell model are the other examples. Specialized models for hippocampal neurons are also proposed.\[[@ref10][@ref11]\] These models mainly differ in the nonlinear function of dendrite, the linear impulse response of the soma filter and modeling refractory phenomenon, but are common in being stochastic models. Netoff\[[@ref20]\] modeled the brain as small world network in which all cells are only coupled to their nearest neighbors, but small numbers of connections are broken and rewired to make long-distance connections at random locations.

In neuronal tissue (gray matter of the brain) neurons can often be divided into two main families, the principal neurons and interneurons. Thus, the models of neurons have their special parameters for each kind of neuron according to their properties. On the other hand the neurons of a neural network can be connected to each other in different lattice forms. In these microscopic level models any change in the ratio of the number of excitatory to that of inhibitory neurons, as well as the structure of the connected neurons and the relevant coefficients can change the type of activity at the model output. In other words, there are too many effective parameters in these models that prevent the model to be suitably tractable. Thus generally, analysis of the microscopic models to discover some physiological facts is very difficult.

In another modeling scheme, called the reduced model, the complexity measure of single neurons is low so that a few neurons can build Electroencephalography (EEG) signals. In\[[@ref21]\] cortex is modeled as a matrix that each component contains excitatory and inhibitory mass of neurons. These models are mainly modified versions of the microscopic models that can lead to better comprehension of cortex dynamics.

For gaining a better insight into the neural functions underlying the emergence of different dynamics, further simplification of the models have been proposed, in which neuronal networks of cortex have been considered as a spatially continuous network. In these macroscopic level models, the properties of a bulk of neurons at a local region have been described by time-space dependent state variables. Typically, these variables can be the mean firing rates or the mean value of the cell (soma) membrane potential, just like there is an average neuron at that location.\[[@ref18][@ref22]\] These kinds of models are also called neural mass models. Analogues to a single neuron\[[@ref6]\] nonlinear function at the average dendrite and linear function at the average soma is expected in a neural mass model. Furthermore, the input of a neural mass model is usually a white noise to account for the influence of all far neuron populations.

Although these models are the simplest ones introduced for cortex, their differential equations possess nonlinear terms, and thus are difficult to analyze. Almost no analytical solutions exist for these equations and they must be qualitatively analyzed. By assessment of the model\[[@ref23][@ref24]\] further insight has been gained into the properties of the dynamics of cortex performance. The strength and time course of cortical inhibitory--inhibitory and inhibitory--excitatory connectivity have been realized to be the effective parameters of the model about the seizure initiation and cessation. Furthermore bifurcation diagrams have been considered to analyze the model behavior.\[[@ref25]\] Also, by linearizing the part of the model responsible for the fast onset activity about its equilibrium solution, Molaee-Ardekani\[[@ref26]\] quantified the characteristics of the resonance frequency of that part of model while free parameters vary. It has been shown that if the parameters are suitably selected, sustained limit cycle activity with high frequency and damped oscillatory behavior in response to brief stimulating inputs can be seen.\[[@ref26]\]

Although linearization and the associated sensitivity analysis are able to provide a great deal of information regarding the conditions under which linear stability is lost, they are in general insufficient for a nonlinear model. In this paper the describing function approach is used to analyze a neural mass model more accurately, in order to show if other parameters of the model can influence the model output activity. To begin, we focus on a single parameter: The excitatory time constant of pyramidal neurons. It is concluded that as well as the excitatory and inhibitory gains, the excitatory time constant of pyramidal neurons has a significant role in initiation and termination of seizures.

MATERIALS AND METHODS {#sec1-2}
=====================

Depth-EEG Model {#sec2-1}
---------------

The first neural mass models only considered the excitatory interactions.\[[@ref27]\] Later, Wilson and Cowan\[[@ref28][@ref29]\] incorporated both excitatory and inhibitory properties of the neuronal populations connected to each other. The state variable was the mean firing rate, and the equations underlying these models are differential equations that can easily be extended into partial differential equations, to account for spatial properties of the neuronal mass. Amari\[[@ref30][@ref31]\] has proposed a spatio-temporal model with the mean soma membrane potential as macroscopic state variable.

Lopes Da Silva *et al*.\[[@ref6]\] and simultaneously Freeman\[[@ref32]\] developed a lumped-parameter population model as a forward model, which is capable to explain the alpha rhythm of the EEG and deals with perceptual processing in the olfactory system. Physical laws they employed are related to two conversion operations involved in the dynamics of each neural ensemble: A mean wave-to-pulse operator at the soma of neurons, and a mean linear pulse-to-wave conversion implemented at a synaptic level.

Wendling *et al*.\[[@ref14]\] used the same approach to interpret depth-EEG signal characteristics in epileptic patients according to physiological facts. Indeed, as a new law, they tried to model realistic epileptiform activity as a result of imbalance between excitatory and inhibitory synaptic gains (i.e., model parameters). In this model each area of hippocampus contains three interacting subsets of neurons. Main cells (i.e., pyramidal cells in the hippocampus or neocortex) receive excitatory and inhibitory feedbacks from other subsets composed of local interneurons. Also in order to represent fast EEG activity such as low-voltage rapid discharges that are often observed in depth-EEG signals at seizure onset, Wendling *et al*. developed a new version of mentioned model based on bibliographical material.\[[@ref14]\] The second class of inhibitory interneurons is added to model as the fourth subset to raise kinetics of the system. In this local model the influence from neighboring areas of hippocampus is represented by an excitatory model input, *x*~*p*~(t), which is represented by Gaussian white noise. Thus the model in\[[@ref14]\] is a stochastic model whose output corresponds to the post-synaptic activity of the first subset, and is interpreted as depth-EEG signal. As mentioned before, in this paper we call it the depth-EEG model.

Liley developed the basic theories of the macroscopic models by adding some new physiological facts.\[[@ref18]\] Thus, effects of the synaptic reversal potentials, which make the amplitude of the respective post-synaptic potentials dependent on the ongoing post-synaptic/somatic membrane potential, have been added to the model. The effect of synaptic reversal potentials is the changing of the linear relation between the mean soma membrane potential and synaptic input to an exponential one. Also transmission of axonal pulses in long range fibers has been added to this model as a function. In other words, Liley model is capable of reproducing the main features of spontaneous human EEG. In particular, autonomous limit cycle and chaotic oscillatory activity in the alpha band (8-13 Hz) has been easily produced.

As shown in [Figure 1](#F1){ref-type="fig"}, in depth-EEG model each synaptic process (excitatory, fast and slow inhibitory processes) is a combination of a linear system (pulse-to-wave conversion implemented at a synaptic level) and a nonlinear function (wave-to-pulse operator at the soma of neurons, i.e., *S*(*v*) =2e~0~/(1+e^*r*(*v*~0~-*v*)^)), where *v* is the average potential of the pre-synaptic cells and *S* (*v*) is the mean firing rate of post-synaptic cells. This process is expressed with the differential (1).

![Organization of the depth-EEG model. Four neuronal subsets: pyramidal cells, excitatory interneurons (A), dendritic projecting interneurons with slow synaptic kinetics (GABAA,slow:*B*) and somatic-projecting interneurons (the grey rectangle) with faster synaptic kinetics (GABAA, fast:*G*). The average pulse density of afferent action potentials is changed into an average inhibitory or excitatory post-synaptic membrane potential using a linear dynamic transfer functions, while this potential is converted into an average pulse density of firing of post-synaptic neuron using a static nonlinear function (from \[[@ref14]\])](JMSS-3-2-g001){#F1}

The impulse response of the second-order linear systems involved in the model is characterized by the time-constant and the peak amplitude (that is reached after one time-constant), as shown in [Figure 2](#F2){ref-type="fig"}. These are two indicators of the system. We call the former "synaptic time constant", and the latter is proportional to the "synaptic gain". It is evident that fast-inhibitory synaptic time is less than slow inhibitory synaptic time. According to\[[@ref14]\] standard values of linear system (*C*~1~ to *C*~2~, *a*, *b*, *g*) and nonlinear function parameters (*e*~0~ *v*~0~ and *r*) are presented in [Table 1](#T1){ref-type="table"}\[[@ref14]\] *y*~0~ − *y*~9~ are the states of the model, which are equivalent to the excitatory post-synaptic potential (EPSP) and inhibitory post synaptic potential (IPSP) of the subpopulations of the neurons, or their corresponding mean firing rates.

![Average standard post-synaptic membrane potentials: excitatory, slow inhibitory and fast inhibitory, respectively, obtained from impulse responses given by *he*(*t*) = *Aate*−*ta*, *hi*(*t*) = *Bbte*−*tb*. and *he*(*t*) = *Ggte*−*tg*, *t*\> 0 for *A*= 3.25, *B*= 22, *G*= 10 (from \[[@ref14]\])](JMSS-3-2-g002){#F2}

###### 

Model parameters, interpretation and values used to produce EEG signals adapted\[[@ref14]\]

![](JMSS-3-2-g003)

![](JMSS-3-2-g004.jpg)

Types of the Model Output Activity {#sec2-2}
----------------------------------

Based on the assumption that at different vigilance and pathological states the physiological parameters of the population of hippocampal neurons may differ from standard ones, Wendling *et al*. have produced six different types of EEG signals by changing the excitation and inhibition synaptic gain parameters (*A*, *B*, and *G*, respectively) in the (0,100) interval, which is a realistic range.\[[@ref14]\] Normal background, sustained discharge of spikes, low voltage rapid activity, slow quasi-sinusoidal activity, sporadic spikes, and slow rhythmic activity. Each type of these signals can be met in the real world (normal and pathological) during different activities of the brain.

[Figure 3](#F3){ref-type="fig"} shows typical output signals of model according to different parameters (i.e., different models of the model space). Since (1) are stochastic differential equations, they should be solved by a stochastic numerical algorithm; but, like\[[@ref33]\] we accept the simplicity of solving our equations by Euler algorithm.

![Six different kinds of the depth-EEG model output](JMSS-3-2-g005){#F3}

We used Gaussian white noise as hippocampcal modulatory input. In order to produce background activity parameters using standard values,\[[@ref14]\] the mean and variance of the input signal should be properly chosen. To this end, we selected μ = 50 and σ = 50.

Motivation to Analyzing Model {#sec2-3}
-----------------------------

It has been observed that if all parameters of the depth-EEG model are supposed to have constant values equal to the standard ones except excitation/inhibition synaptic gain parameters (*A, B*, and *G*) some differences between the model output and real signals, especially in dominant frequency of low-voltage rapid activity, are obvious. It is seen that dominant frequency of low-voltage rapid activity in the model output is non-flexible, but in the real world it is not the case. In fact, frequency of the seizure activity is not unique for different subjects and also for different seizures of a person. Although these changes are small they should be considered.

To comply with the possible changes of frequency, some other parameters must be assumed to be variable over time. The question is then which parameters could be effective in this order. According to the simulation results synaptic time constant of pyramidal cell excitatory activity affects the frequency content of the model output. The block corresponding to pyramidal cell excitatory synaptic process of the model is bolded in [Figure 1](#F1){ref-type="fig"}.

Thus analyzing the model is necessary to indicate how the time constant of pyramidal cell excitatory activity influences the model output frequency and even the activity type.

Describing Function Concept {#sec2-4}
---------------------------

A nonlinear system can be studied by simply replacing the nonlinear operations by approximating linear operations. When the operating point of the system is changing by time, instead of repeating the linearization about new operating points, one must determine the operation performed by the nonlinear element on an input signal of finite size and approximate this in some way by a linear operation. Thus the same nonlinearity when driven by inputs of different forms, or even when driven by inputs of the same form but of different magnitude, must be used. The approximation of a nonlinear operation by a linear one which depends on some properties of the input is called quasi-linearization.

Three basic signal forms with which the quasi-linear approximators for nonlinear operators have been driven until now are: *Bias*, *Sinusoid*, and *Gaussian process*. Also, different linear combinations of these basic input signals are considered. The quasi-linear functions that describe approximately the transfer characteristics of the nonlinearity are termed describing functions (DF).\[[@ref34][@ref35]\] The concept is to see when, for example, a sinusoidal signal is applied to the nonlinear function, what would be the fundamental component of the output signal. The ratio of the (output/input) phasors at the fundamental harmonic is defined as DF of the nonlinear function for sinusoidal input. This concept is mainly used in the detection of the system\'s limit cycle, i.e., prediction of its approximate amplitude and frequency, to be able to design suitable controllers for nonlinear systems.

To obtain the DF for a sinusoidal input, *x*\[*n*\] = *a*~0~ sin\[*ω*~0~*n*\], it must be considered that at steady state the output of the nonlinear characteristic, *y*\[*n*\], is periodic but, in general non-sinusoidal. Assuming the nonlinearity to be symmetric about zero, the Fourier series becomes:

![](JMSS-3-2-g006.jpg)

where *Y~k~* and *φ~k~* are the amplitude and the phase shift of the *k*^th^ harmonic component of the output *y*\[*n*\], respectively. In the sinusoidal-input DF analysis, only the fundamental harmonic component of *y*\[*n*\], *Y*~1~ is considered, because the higher-harmonics are often of smaller amplitude than the amplitude of the fundamental component.

Accordingly, the sinusoidal-input DF of a nonlinear element, *N*(*a*~0~,*ω*~0~), is defined as the output/input ratio of the fundamental harmonic component.

![](JMSS-3-2-g007.jpg)

DESCRIBING FUNCTION OF THE MODEL\'s NONLINEAR BLOCK {#sec1-3}
===================================================

To understand accurate role of the depth-EEG model variables the nonlinear blocks can be approximated by a linear function equal to their suitable DF.\[[@ref34][@ref35]\] There are four nonlinear blocks, which in [Figure 1](#F1){ref-type="fig"} are labeled by numbers 1-4: *N*~1~, *N*~2~, *N*~3~, *N*~4~.

When the periodic activity is seen at the model output, input of the nonlinear blocks can be assumed to be sinusoidal signals with a bias, such that "sine-plus-bias DF" would be sufficient to analyze the model. Since the nonlinear function *S*(*v*) does not cause any delay, DF has no phase-shift term. So only the magnitude of the output\'s fundamental harmonic and the bias value must be determined. Accordingly, to obtain the DF for S(v) the integrations (2) must be computed:

![](JMSS-3-2-g008.jpg)

where *a*~0~, *ω*~0~ and *y*~0~ are amplitude, frequency and the bias value of input signal of the nonlinear block. Moreover, *a^out^*~0~ and *y^out^*~0~ are the amplitude and bias of the output signal of the nonlinear block.

As mentioned before, the nonlinear function of depth-EEG model is *S*(*v*)=2e~0~/(1+e^r(v~0~--v)^). So, two integrations of (4) would be:

![](JMSS-3-2-g009.jpg)

In order to simplify the above integrations to be analytically calculated, the nonlinear function is approximated by a piecewise linear function. Since the sigmoid function *S*(*v*) is a smooth and strictly monotonic function, the piecewise linear approximation is a convenient one. Although the larger the number of these linear pieces, the better the nonlinear function would be reconstructed; but, using only three linear pieces suffice in our analytical calculations. This is firstly because the approximated nonlinear functions do not lead to alter the model outputs significantly, and secondly the effect of other parameters of the model, which is the goal of our analysis, will be evident without any necessity to take the complexities into account. The approximated function is shown in [Figure 4](#F4){ref-type="fig"}, and its related equation is displayed in (6).

![The non-linear function is approximated by a piece-wise linear function](JMSS-3-2-g010){#F4}

![](JMSS-3-2-g011.jpg)

Accordingly, the 0 -- 2π interval of integration must be partitioned into smaller intervals, in which *S* (*v*) is linear. For example, to this end in [Figure 5](#F5){ref-type="fig"} *θ*~1~ and *θ*~2~ indicates the borders of sub-integrations. [Figure 5](#F5){ref-type="fig"} shows the input and output of the nonlinear function *S* (*v*) where the sinusoidal amplitude equals to *a*~0~ = 16 and its bias value is *y*~0~ = --2.5.

![Example of the time intervals used in order to integration, according to the approximated nonlinear function](JMSS-3-2-g012){#F5}

It is evident that depending on the values of *a*~0~ and *y*~0~ different situations occur, for which the integration intervals must be indicated separately

Before describing these different situations, note that *θ*~1~ and *θ*~2~ are equal to ![](JMSS-3-2-g013.jpg) and ![](JMSS-3-2-g014.jpg) respectively. Since *x*~1~ \< *x*~1~ then *θ*~1~ \< *θ*~2~. Furthermore, due to the symmetry of the sinus signal, both *θ*~1~ and *θ*~2~ have values between --π/2 to π/2. Thus, three different cases have been detected.

Case 1: θ~1~ \< 0,θ~2~ \< 0, or equivalently ![](JMSS-3-2-g015.jpg), thus

![](JMSS-3-2-g016.jpg)

Special cases of this situation are:

Case 1a: --π/2=θ~1~, θ~2~ \< 0 that occurs when ![](JMSS-3-2-g015.jpg) and ![](JMSS-3-2-g017.jpg):

![](JMSS-3-2-g018.jpg)

Case 1b: --π/2=θ~1~,--π/2=θ~2~, that occurs when ![](JMSS-3-2-g019.jpg), then:

Ŝ(*y*~0~+*a*~0~sin\[*ω*~0~*n*+*α*\]) = 2e~0~      (9)

Case 2: θ~1~ \< 0, θ~2~ \> 0, or equivalently ![](JMSS-3-2-g020.jpg) and ![](JMSS-3-2-g021.jpg)

![](JMSS-3-2-g022.jpg)

that its special cases are:

Case 2A: --π/2=θ~1~,θ~2~ \> 0, or equivalently ![](JMSS-3-2-g023.jpg) and ![](JMSS-3-2-g017.jpg)

![](JMSS-3-2-g024.jpg)

Case 2b: *θ*~1~ \<0,*θ*~2~ = π/2, or ![](JMSS-3-2-g025.jpg) and ![](JMSS-3-2-g026.jpg)

![](JMSS-3-2-g027.jpg)

Case 2c: *θ*~1~ = --π/2, *θ*~2~ = π/2, or ![](JMSS-3-2-g025.jpg) and ![](JMSS-3-2-g017.jpg)

![](JMSS-3-2-g028.jpg)

Case 3: 0 \< *θ*~1~, 0 \< *θ*~2~, or equivalently ![](JMSS-3-2-g029.jpg)

![](JMSS-3-2-g030.jpg)

with the special cases:

Case 3a: *θ*~1~ \>0, *θ*~2~ = π/2, or ![](JMSS-3-2-g025.jpg)

![](JMSS-3-2-g031.jpg)

Case 3b: *θ*~1~ = π/2, *θ*~2~ = π/2, or ![](JMSS-3-2-g032.jpg)

Ŝ(*y*~0~+*a*~0~sin\[*ω*~0~*n*+*α*\]) = 0      (16)

After computing the integration and some algebraic operations, the result of all of these cases for the values of **y*^out^*~0~ and **a*^out^*~0~ can be written in the same closed form as following:

![](JMSS-3-2-g033.jpg)

For example for case 2c that *θ*~1~ = --π/ 2,*θ*~2~= π/ 2, the output of the *S*(*v*) in response to a sinusoidal input with amplitude *a*~0~ and the bias value *y*~0~, is a periodic signal whose fundamental harmonic amplitude is **a*^out^*~0~ and the *dc* component have the following values:

![](JMSS-3-2-g034.jpg)

ANALYSIS OF DEPTH-EEG MODEL AT A SPECIAL PARAMETER VECTOR *A* = 21, *B* = 31 and *G* = 21 {#sec1-4}
=========================================================================================

To investigate whether other variables of the depth-EEG model are effective in changing frequency contents of the model output and/or initiation or termination of the periodic seizure-related activity, values of *A, B, G* have been fixed. The problem is limited to assess the role of a single parameter, here of the pyramidal excitatory time constants (τ). We want to show the validity of the following:

Periodic activity can appear or disappear by changing the value of the pyramidal excitatory time constants (τ).Frequency of the periodic activity is dependent on the τ.

Although detail of the analysis approach is dependent to the values of *A, B, G*, but the main approach is the same for all parameters: Different assumptions about the cases of nonlinear blocks should be tested until valid cases for all range of τ would be found. Thus, in this paper to describe how the analysis would be done in details, as an instance, a special value of *A, B, G* is considered: *A* = 21 *B* = 31 *G* = 21. These values are manually selected and there is a periodic activity for *a* = 100, *b* = 50, g = 500, τ = 50 at the model output. To obtain a general result, the depth-EEG model could be analyzed for other values of *A, B, G*. Also, the role of other time constants could be indicated.

The Linear Part of the Depth-EEG Model {#sec2-5}
--------------------------------------

The linear filters of the model are low-pass filters. In simulation with MATLAB software they are treated as discrete filters by the following transfer functions.

![](JMSS-3-2-g035.jpg)

where *T* is the sampling time by which the differential equations are numerically solved. The *dc* gain of filters (*ω* = 0) will be *A/τ, B/b, G/g*, and *A/a*. The frequency responses of these filters are as follows:

![](JMSS-3-2-g036.jpg)

The Nonlinear Part of the Depth-EEG Model {#sec2-6}
-----------------------------------------

According to the monotonic nature of nonlinear functions, it is clear that at every set of the model parameters there are unique state signals and output signal for a special input signal. Thus, depending to the dc value and the amplitude of the fundamental harmonic of the nonlinear blocks, they may work just in one of the abovementioned cases. To determine the valid cases of under the assessment parameters, we start our analysis by assuming some test fixtures for the nonlinear blocks. By analysis of the equivalent linear system, the dc value and the amplitude of the fundamental harmonic of the model output signal can be easily computed. By using these values, we obtain the input of every nonlinear block and check it to see whether or not the assumed case for that block is correct. If the assumed cases were not the correct ones, i.e., leads to some contradictions in the system, we must continue the analysis procedure by assuming some new cases. Once the valid set of cases for the model parameter is found, the analysis is finalized for that parameter.

We have checked different assumptions by trial and error, but here just those assumptions that have been found to be correct at least for some values of τ in under the assessment range (0,100), are reported. Thus, at first we suppose that for *A* = 21, *B* = 31, *G* = 21, τ = 50 , *a* = 100, *b* = 50, *g* = 500, *N*~1~ is at case 2c, *N*~2~ and *N*~3~ are at case 1b, and *N*~4~ is also at case 2c.

Assumption 1: *N*~1~ and *N*~4~ are at Case 2c, *N*~2~ and *N*~3~ are at Case 1b {#sec2-7}
--------------------------------------------------------------------------------

Assuming the nonlinear blocks of the model at the above mentioned cases, the output of *N*~2~ and *N*~3~ are equal to 2e~0~. In other words, according to [Figure 6](#F6){ref-type="fig"}, *y*~1~ and *y*~2~ are, respectively, obtained by passing a *dc* value (2e~0~) multiplied by *C*~2~ and *C*~4~, through filters *H~a~* (*e^jω^*) and *H~b~* (*e^jω^*). Thus, *y*~2~ and *y*~2~ are constant signals, respectively, equal to 2C~2~ Ae~0~/*a* + *μ* and 2C~4~ Be~0~/*b*, where *μ* is the mean value of model input noise. Under these assumptions the effective loop of the model in which the oscillation is expected to be generated is shown in [Figure 6](#F6){ref-type="fig"}.

![The approximated model by using describing function under Assumption](JMSS-3-2-g037){#F6}

Furthermore, since *N*~1~ is supposed to work in Case 2c it would be replaced by its describing function, i.e., a gain equal to re~0~/In (99) (18). Similarly *N*~4~ is replaced by the same gain.

Assuming the y~0~ to be the notation for the *dc* value of the model output (input of the *N*~1~), and *a*~0~ for its amplitude at the fundamental frequency ((*y*~0~+*a*~0~sin\[ω~0~*n*+α\])), also by defining ![](JMSS-3-2-g038.jpg) and ![](JMSS-3-2-g039.jpg) as the outputs of *N*~1~ and *N*~4~ respectively, the (18) can be used to obtain:

![](JMSS-3-2-g040.jpg)

and

![](JMSS-3-2-g041.jpg)

Accordingly, the *dc* loop equation

![](JMSS-3-2-g042.jpg)

can be written as:

![](JMSS-3-2-g043.jpg)

Since we are analyzing the effect of changing the value of *τ* on the system in an special case in which *A* = 21, *B* = 31, *G* = 21, *a* = 100, *b* = 50, *g* = 500, by substituting these values in the (23) the amount of *y*~0~ versus *τ* is obtained as follows:

![](JMSS-3-2-g044.jpg)

On the other hand using the Mason rule for the loop, it is evident that the model is a simple linear system with the frequency response as follows:

![](JMSS-3-2-g045.jpg)

Because the model input is a Gaussian white noise, the power spectrum of the model output by the *assumption 1* about the situations in which the nonlinear blocks work would be as follows:

![](JMSS-3-2-g046.jpg)

Checking Assumption 1 {#sec2-8}
---------------------

It is noteworthy that for some values of *τ* the (24) may not hold. In other words it must be checked when the conditions for Assumption 1 hold. The variance of the noise is assumed to be negligible in such a way that dc values of the state signals of the model are sufficient to check the conditions.

To check whether or not the above assumption about cases of nonlinear blocks is true, the following conditions must be assessed.

To check if *N*~1~ is in case 2c, it must be checked whether ![](JMSS-3-2-g047.jpg)To check if *N*~2~ is in case 1b: ![](JMSS-3-2-g048.jpg)To check if *N*~3~ is in case 1b: ![](JMSS-3-2-g049.jpg)To check if *N*~4~ is in case 2c, it must be checked whether: ![](JMSS-3-2-g050.jpg).

The above statements would be reduced to the following ones by some simple algebraic manipulations:

*N*~1~ is in case 2c if: 1.9774 \< τ \< 90.9735.*N*~2~ is in case 1b if: 9.1190 \< τ.*N*~3~ is in case 1b if: 2.9605 \< τ.*N*~4~ is in case 2c for all values: 0 \< τ \< 100.

Totally the condition for which the above analysis is valid is as following:

9.1190 \< τ \< 90.9735      (27)

In other words the (26) is just valid for those values of τ with the limits of (27). This filter is equivalent to a resonant low-pass filter, which according to it a periodic activity with the resonant frequency is expected at the model output. The magnitude of frequency response (26) is plotted in [Figure 7](#F7){ref-type="fig"} for some values of τ in the interval of (27).

![Magnitude of the frequency response of (26) for some values of t in the interval: 9.1190 \< τ \< 90.9735](JMSS-3-2-g051){#F7}

Assumption 2: *N*~1~ is at Case 2c, *N*~2~, *N*~3~ and *N*~4~ are at Case 1b {#sec2-9}
----------------------------------------------------------------------------

Since in practice we consider *τ* to be in the (0, 100\] interval, the system must be analyzed again for 0 @ τ \< 9.1190 and 90.9735 \< *τ* \< 100. In this section the analysis is done for 0 \< *τ* \< 9.1190 by another set of assumptions about the cases of the nonlinear blocks. We now suppose that *N*~4~ also works in the case 1b, i.e., its output is equal to 2e~0~.

For *Assumption 2*, there is no loop in the system such that for *dc* we have the following equation:

![](JMSS-3-2-g052.jpg)

which leads to *y*~0~ = 0.7665 for *μ* = 50. In this situation the model output is just a filtered white noise whose power spectrum is equal to:

![](JMSS-3-2-g053.jpg)

We must check again the range of *τ* for which this situation occurs. It is noteworthy again that the variance of the input noise is supposed to be sufficiently small such that just the *dc* values indicate the cases for each nonlinear block.

To check if *N*~1~ is in case 2c, it must be checked whether ![](JMSS-3-2-g047.jpg)To check *N*~2~ is in case 1b: ![](JMSS-3-2-g048.jpg)To check *N*~3~ is in case 1b: ![](JMSS-3-2-g049.jpg)To check *N*~4~ is in case 1b: ![](JMSS-3-2-g054.jpg)

By the obtained value *y*~0~ = 0.7665, (1) always holds, (2) requires *τ* ≤ 98.5814, (3) requires *τ* ≤ 24.6459, and (4) requires *τ* ≤ 9.1239, which all of them hold under the assessment range 0 \< *τ* \< 9.1190.

The spectrum of (29) is clearly a first order low-pass filter such that no periodic activity is seen at the model output for 0 \< *τ* \< 9.1190.

Assumption 3: *N*~1~, *N*~2~ and *N*~3~ are at Case 1b, *N*~4~ is at Case 2c {#sec2-10}
----------------------------------------------------------------------------

For 90.9735 \< *τ* ≤ 100 the analysis is continued by assuming the case 1b for *N*~1~, *N*~2~ and *N*~3~, and case 2c for *N*~4~. In this situation too, there is no closed loop, and the equation for *dc* value is as the following:

![](JMSS-3-2-g055.jpg)

that leads to:

![](JMSS-3-2-g056.jpg)

For this situation, the spectrum of the model output is as follows:

![](JMSS-3-2-g057.jpg)

Again by assuming the small variance for the input noise, the range of suitable *τ* for this situation is obtained.

To check if *N*~1~ is in case 1b: ![](JMSS-3-2-g058.jpg)To check if *N*~2~ is in case 1b: ![](JMSS-3-2-g048.jpg)To check if *N*~3~ is in case 1b: ![](JMSS-3-2-g049.jpg)To check if *N*~4~ is in case 2c, it must be checked whether: ![](JMSS-3-2-g050.jpg)

These conditions lead to 94.7794 \< *τ* \< 119.1861. In other words for 94.7794 \< *τ* \< 100 the (32) presents the valid power spectrum of the model output. This means that for 94.7794 \< *τ* \< 100 no periodic activity is expected at the model output. For 90.9735 \< *τ* \< 94.7794, the analysis must be done again.

Assumption 4: *N*~1~ is at Case 2a, *N*~2~, *N*~3~ and *N*~4~ are at Case 1b {#sec2-11}
----------------------------------------------------------------------------

For 90.9735 \< *τ* \< 94.7794 we assume that *N*~1~ is at case 2a, *N*~2~, *N*~3~ and *N*~4~ are at case 1b. *N*~1~ behaves like a linear function by the following gain (*θ*~2~ \> 0):

![](JMSS-3-2-g059.jpg)

The spectrum of the model output would be very similar to 26:

![](JMSS-3-2-g060.jpg)

This filter is again a resonant low-pass filter.

It is noteworthy that while *A* = 21, *B* = 31, *G* = 21, *a* = 100, *b* = 50, *g* = 500, by the mentioned assumptions the whole interval of \[0,100\] is covered. In other words, it is an easy task to show that any other possible assumption about the cases leads to a contradiction (incorrect result).

RESULTS {#sec1-5}
=======

Based on the results of previous sections, it is evident that according to the concept of describing function, the complex function of the nonlinear model is approximated by a filter at A = 21, B = 31, G = 21. This filter is low-pass for 94.7794 \< τ \< 100 and \< τ \< 9.1190, but a resonant low-pass filter for 9.1190 \< τ \< 94.7794. In other words the role of the τ in changing the system function is now evident. But from the seizure initiation point of view, it must be discussed at which times the periodic activity appears in the model output. In [Figure 8](#F8){ref-type="fig"} changing the frequency response of the model by τ is shown.

![Magnitude of the frequency response of the approximated model for *A* = 21, *B* = 31, and *G* = 31 for some values of t low-pass filter, the quality factor of the filter also must be considered to indicate whether the model output has a periodic activity or not. The quality factor of the resonant filter versus the value of is plotted in Figure 9 for 9.1190 \< t \< 94.7794](JMSS-3-2-g061){#F8}

For the low-pass filters, 94.7794 \< τ \< 100 and 0 \< τ \< 9.91190, clearly there is no periodic activity in the model output, i.e., the normal activity. But for the resonant low-pass filter, the quality factor of the filter also must be considered to indicate whether the model output has a periodic activity or not. The quality factor of the resonant filter versus the value of τ is plotted in [Figure 9](#F9){ref-type="fig"} for 9.1190 \< τ \< 94.7794.

![The quality factor values of the approximated resonant filter for the depth-EEG model, versus t, for. 9.1190 \< t \< 94.7794. In this figure the 0 \< t \< 9.1190 and 94.7794 \< t\< 100 regions for which the model is equivalent to low-pass filter are specified. According to the quality factor values for 90.3375 \< t \< 94.7794 is such small that the periodic activity of the model output is diminished](JMSS-3-2-g062){#F9}

It is evident from [Figure 9](#F9){ref-type="fig"} that the quality factor of the resonant filter is decreasing as τ is increasing. The quality factor of the filter is so low for 90.3375 \< τ \< 94.7794 such that periodic activity seems to be diminished.

As a result we can say that at the constant value of all parameters, for example for A = 21, B = 31 G = 21, a = 100, b = 50, and g = 500, variation of the pyramidal excitatory constant time can initiate or/and terminate the periodic activity (seizure) by itself.

Also, the resonant frequency of the filter versus τ is shown in [Figure 10](#F10){ref-type="fig"}. According to this figure variation of the dominant frequency of the model output as τ changes is obvious. Thus the pyramidal excitatory time constant is able to modify the frequency content of the model output. As it was mentioned before, unlike real signals, depth-EEG model suffered from its unique output seizure frequency that could be eliminated by taking τ as a new model parameter.

![The resonant frequency of the approximated resonant filter for the depth-EEG model, versus τ, for 9.1190 \< τ \< 94.7794. In this figure the 0 \< τ \< 9.1190. and 94.7794 \< τ\< 100 regions for which the model is equivalent to low-pass filter are specified](JMSS-3-2-g063){#F10}

These theoretical results are matched with that of simulation results.

DISCUSSION AND CONCLUSION {#sec1-6}
=========================

In this paper we showed the effects of changing the parameter τ for one special value of A, B, and G (A = 21, B = 31, and G = 21). However, the derived results can be generalized for other values. We can say that the excitatory time constant of the pyramidal neurons is another parameter that affects the frequency of seizure-related signals, i.e., decreasing the value of τ reduces the output frequency, and is also able to initiate and terminate seizure. That is, while a periodic activity is seen on the model output, it may disappear by changing τ (depending on synaptic gains), and vice versa: Sometimes changing τ may cause the appearance of rhythmic activity.

In other words, in the absence of any variation of synaptic gains, pyramidal cell excitatory time (τ) is able to induce or stop instability. This means that without any need for synaptic gains to return to their standard values, seizures could be stopped by changing parameter τ. This observation is in accordance with what Liley reported on his attempt\[[@ref22]\] to find mechanism of both formation and cessation of generalized seizure activity, and propagation of focal epileptiform activity. As Liley declared, changes in fast GABAergic neurotransmission are necessary to induce transition to epileptiform activity, but it can be done through different parametric context. Weakening inhibition can be encountered by both reducing peak amplitude of the IPSP and its decay time constant. In this paper we showed that this may be the fact also for the excitatory processes: Every change in excitation process, either through change of the gain (A) or excitatory time constant τ can cause the instability or stop it. It proposes that the effective time of the excitatory and inhibitory processes, not only their gains, are important in changing the type of EEG activity. In fact integrating the excitatory and inhibitory processes seems to be more meaningful when they are compared with investigate if they are balanced or not. This fact must be assessed from the point of view of physiological mechanisms. The physiological meaning of changing the time constant of pyramidal neurons may resemble the change of the speed of neurotransmitters' secretion at the synapses.

In this paper we focused on the effects of the excitatory activity of pyramidal neurons, but it deserves to evaluate the effects of other time constants in the model, e.g. *a*, *b*, and *g*.

**Source of Support:** Nil

**Conflict of Interest:** None declared