You may also check out my GoogleScholarProfile .

Book

Papers

Move over the publication to see a short/description or the Abstract or follow the (DOI) link to the article at the journal page.
  1. Pathological responses to single pulse electrical stimuli in epilepsy: the role of feedforward inhibition. G.J. Hebbink, G.J.M. Huiskamp, S.A. van Gils, F.S.S. Leijten and H.G.E. Meijer, European Journal of Neuroscience, 2019, In Press Abstract
    Description:

    Delineation of epileptogenic cortex in focal epilepsy patients may profit from single pulse electrical stimulation during intracranial EEG recordings. Single pulse electrical stimulation evokes early and delayed responses. Early responses represent connectivity. Delayed responses are a biomarker for epileptogenic cortex, but up till now, the precise mechanism generating delayed responses remains elusive. We used a data-driven modelling approach to study early and delayed responses. We hypothesized that delayed responses represent indirect responses triggered by early response activity and investigated this for 11 patients. Using two coupled neural masses we modelled early and delayed responses by combining simulations and bifurcation analysis. An important feature of the model is the inclusion of feedforward inhibitory connections. The waveform of early responses can be explained by feedforward inhibition. Delayed responses can be viewed as second‐order responses in the early response network which appear when input to a neural mass falls below a threshold forcing it temporarily to a spiking state. The combination of the threshold with noisy background input explains the typical stochastic appearance of delayed responses. The intrinsic excitability of a neural mass and the strength of its input influence the probability at which delayed responses to occur. Our work gives a theoretical basis for the use of delayed responses as a biomarker for the epileptogenic zone, confirming earlier clinical observations. The combination of early responses revealing effective connectivity, and delayed responses showing intrinsic excitability, makes single pulse electrical stimulation an interesting tool to obtain data for computational models of epilepsy surgery.

    DOI (open access)
  2. Homoclinic saddle to saddle-focus transitions in 4D systems. M. Kalia, Yu.A. Kuznetsov and H.G.E. Meijer, Nonlinearity, June 2019, Vol. 32, No. 6, pp 2024-2054 Abstract
    Description:

    A saddle to saddle-focus homoclinic transition when the stable leading eigenspace is 3-dimensional (called the 3DL-bifurcation) is analyzed. Here a pair of complex eigenvalues and a real eigenvalue exchange their position relative to the imaginary axis, giving rise to a 3-dimensional stable leading eigenspace at the critical parameter values. This transition is different from the standard Belyakov bifurcation, where a double real eigenvalue splits either into a pair of complex-conjugate eigenvalues or two distinct real eigenvalues. In the wild case, we obtain sets of codimension 1 and 2 bifurcation curves and points that asymptotically approach the 3DL-bifurcation point and have a structure that differs from that of the standard Belyakov case. We give an example of this bifurcation in a perturbed Lorenz-Stenflo 4D ODE model.

    DOI
  3. A Comparison of Evoked and Non-evoked Functional Networks, G.J. Hebbink, D. van Blooijs, G.J.M. Huiskamp, F.S.S. Leijten, S.A. van Gils and H.G.E. Meijer, Brain Topography, May 2019, Vol. 32, No. 3, pp 405-417 Abstract
    Description:

    The growing interest in brain networks to study the brain's function in cognition and diseases has produced an increase in methods to extract these networks. Typically, each method yields a different network. Therefore, one may ask what the resulting networks represent. To address this issue we consider electrocorticography (ECoG) data where we compare three methods. We derive networks from on-going ECoG data using two traditional methods: cross-correlation (CC) and Granger causality (GC). Next, connectivity is probed actively using single pulse electrical stimulation (SPES). We compare the overlap in connectivity between these three methods as well as their ability to reveal well-known anatomical connections in the language circuit. We find that strong connections in the CC network form more or less a subset of the SPES network. GC and SPES are related more weakly, although GC connections coincide more frequently with SPES connections compared to non-existing SPES connections. Connectivity between the two major hubs in the language circuit, Broca's and Wernicke's area, is only found in SPES networks. Our results are of interest for the use of patient-specific networks obtained from ECoG. In epilepsy research, such networks form the basis for methods that predict the effect of epilepsy surgery. For this application SPES networks are interesting as they disclose more physiological connections compared to CC and GC networks.

    DOI
  4. Evoked directional network characteristics of epileptogenic tissue derived from Single Pulse Electrical Stimulation. D. van Blooijs, F.S.S. Leijten, P.C. van Rijen, H.G.E. Meijer and G.J.M. Huiskamp, Human Brain Mapping. 2018, Vol. 39, No. 11, pp. 4611-4622. Abstract
    Description:

    Objective: We investigated effective networks constructed from Single Pulse Electrical Stimulation (SPES) in epilepsy patients who underwent intracranial electrocorticography (ECoG). Using graph analysis, we compared network characteristics of tissue within and outside the epileptogenic area.
    Methods: In 21 patients with subdural electrode grids (1 cm interelectrode distance), we constructed a binary, directional network derived from SPES early responses (<100 ms). We calculated in-degree, out-degree, betweenness centrality, the percentage of bidirectional, receiving and activating connections, and the percentage of connections towards the (non-)epileptogenic tissue for each node in the network. We analyzed whether these network measures were significantly different in seizure onset zone (SOZ)-electrodes compared to non-SOZ electrodes, in resected area (RA)-electrodes compared to non-RA electrodes, and in seizure free compared to not-seizure free patients.
    Results: Electrodes in the SOZ/RA showed significantly higher values for in-degree and out- degree, both at group level, and at patient level, and more so in seizure free patients. These differences were not observed for betweenness centrality. There were also more bidirectional and fewer receiving connections in the SOZ/RA in seizure free patients. It appears that the SOZ/RA is densely connected with itself, with only little input arriving from non-SOZ/non-RA electrodes.
    Conclusion: These results suggest that meso-scale effective network measures are different in epileptogenic compared to normal brain tissue. Local connections within the SOZ/RA are increased and the SOZ/RA is relatively isolated from the surrounding cortex.
    Significance: This offers the prospect of enhanced prediction of epilepsy-prone brain areas using SPES.

    DOI
  5. Numerical Bifurcation Analysis of Double $+1$ Multiplier in Z_3-Symmetric Maps. R. Mazrooei-Sebdani, Z. Eskandari-Shahraki, H.G.E. Meijer, J. Diff. Eqs. and Appl. 2018, Vol. 24, No. 9, pp. 1402-1413. Abstract
    Description:

    We consider a resonance 1:1 bifurcation with Z3-symmetry in a discrete-time dynamical system. We employ standard normalization and center manifold reduction techniques. With this, we obtain the normal form up to cubic order and also explicit formulas for the critical normal form coefficients. We provide an implementation of our algorithm that can be used in the numerical bifurcation toolbox matcontm. We illustrate our analysis with a Cournot triopoly model from economics.

    DOI
  6. Numerical bifurcation analysis of homoclinic orbits embedded in one-dimensional manifolds of maps. N. Neirynck, W. Govaerts, Yu.A. Kuznetsov, H.G.E. Meijer, ACM Trans. Math. Softw. 2018, Vol. 44, No. 3, Article 25 (19 pages) Abstract
    Description:

    We describe new methods for initializing the computation of homoclinic orbits for maps in a state space with arbitrary dimension and for detecting their bifurcations. The initialization methods build on known and improved methods for computing one-dimensional stable and unstable manifolds. The methods are implemented in MatcontM, a freely available toolbox in Matlab for numerical analysis of bifurcations of fixed points, periodic orbits and connecting orbits of smooth nonlinear maps. The bifurcation analysis of homoclinic connections under variation of one parameter is based on continuation methods and allows to detect all known codimension 1 and 2 bifurcations in 3D maps, including tangencies and generalized tangencies. MatcontM provides a graphical user interface, enabling interactive control for all computations. As the prime new feature we discuss an algorithm for initializing connecting orbits in the important special case where either the stable or unstable manifold is one-dimensional, allowing to compute all homoclinic orbits to saddle points in three-dimensional maps. We illustrate this algorithm in the study of the adaptive control map, a 3D map introduced in 1991 by Frouzakis, Adomaitis and Kevrekidis, to obtain a complete bifurcation diagram of the resonance horn in a 1:5 Neimark-Sacker bifurcation point, revealing new features.

    DOI
  7. Phenomenological network models: lessons for epilepsy surgery. G.J. Hebbink, H.G.E. Meijer, G. Huiskamp, S.A. van Gils, F.S.S. Leijten, Epilepsia, 2017, Vol 58(10), pp. e147-151. Abstract
    Description:

    The current opinion in epilepsy surgery is that successful surgery is about removing pathological cortex in the anatomic sense. This contrasts with recent developments in epilepsy research where epilepsy is seen as a network disease. Computational models offer a framework to investigate the influence of networks, as well as local tissue properties, and to explore alternative resection strategies. Here we study, using such a model, the influence of connections on seizures and how this might change our traditional views of epilepsy surgery. We use a simple network model consisting of four interconnected neuronal populations. One of these populations can be made hyperexcitable, modeling a pathological region of cortex. Using model simulations the effect of surgery on the seizure rate is studied. We find that removal of the hyperexcitable population is, in most cases, not the best approach to reduce the seizure rate. Removal of normal populations located at a crucial spot in the network, the ''driver'', is typically more effective in reducing seizure rate. This work strengthens the idea that network structure and connections may be more important than localizing the pathological node. This can explain why lesionectomy may not always be sufficient.

    DOI PDF
  8. Synaptic damage underlies EEG abnormalities in postanoxic encephalopathy: a computational study. B.J. Ruijter, J. Hofmeijer, H.G.E. Meijer, M.J.A.M. van Putten, Clinical Neurophysiology, Vol. 218(9) pp. 1682-1695, 2017 Abstract
    Description:

    Objective: In postanoxic coma, EEG patterns indicate the severity of encephalopathy and typically evolve in time. We aim to improve the understanding of pathophysiological mechanisms underlying these EEG abnormalities.
    Methods: We used a mean field model comprising excitatory and inhibitory neurons, local synaptic connections, and input from thalamic afferents. Anoxic damage is modeled as aggravated short-term synaptic depression, with gradual recovery over many hours. Additionally, excitatory neurotransmission is potentiated, scaling with the severity of anoxic encephalopathy. Simulations were compared with continuous EEG recordings of 155 comatose patients after cardiac arrest.
    Results: The simulations agree well with six common categories of EEG rhythms in postanoxic encephalopathy, including typical transitions in time. Plausible results were only obtained if excitatory synapses were more severely affected by short-term synaptic depression than inhibitory synapses.
    Conclusions: In postanoxic encephalopathy, the evolution of EEG patterns presumably results from gradual improvement of complete synaptic failure, where excitatory synapses are more severely affected than inhibitory synapses. The range of EEG patterns depends on the excitation-inhibition imbalance, probably resulting from long-term potentiation of excitatory neurotransmission.

    DOI
  9. Estimation and Identifiability of Model Parameters in Human Nociceptive Processing Using Yes-No Detection Responses to Electrocutaneous Stimulation. Huan Yang, H.G.E. Meijer, J.R. Buitenweg, S.A. van Gils, Frontiers in Psychology 2016, Vol. 7, pp. 1884 Abstract
    Description:

    Healthy or pathological states of nociceptive subsystems determine different stimulus-response relations measured from quantitative sensory testing. In turn, stimulus-response measurements may be used to assess these states. In a recently developed computational model, six model parameters characterize activation of nerve endings and spinal neurons. However, both model nonlinearity and limited information in yes-no detection responses to electrocutaneous stimuli challenge to estimate model parameters. Here, we address the question whether and how one can overcome these difficulties for reliable parameter estimation. First, we fit the computational model to experimental stimulus-response pairs by maximizing the likelihood. To evaluate the balance between model fit and complexity, i.e., the number of model parameters, we evaluate the Bayesian Information Criterion. We find that the computational model is better than a conventional logistic model regarding the balance. Second, our theoretical analysis suggests to vary the pulse width among applied stimuli as a necessary condition to prevent structural non-identifiability. In addition, the numerically implemented profile likelihood approach reveals structural and practical non-identifiability. Our model-based approach with integration of psychophysical measurements can be useful for a reliable assessment of states of the nociceptive system.

    DOI
  10. Initialization of homoclinic solutions near Bogdanov-Takens points: Lindstedt-Poincare compared with regular perturbation method. B. Al-Hdaibat, W. Govaerts, Yu.A. Kuznetsov, H.G.E. Meijer, SIAM Journal on Applied Dynamical Systems 2016, Vol. 15, No. 2, pp. 952-980 Abstract
    Description:

    To continue a branch of homoclinic solutions starting from a Bogdanov-Takens (BT) point in parameter and state space, one needs a predictor based on asymptotics for the bifurcation parameter values and the corresponding small homoclinic orbits in the phase space. We derive two explicit asymptotics for the homoclinic orbits near a generic BT point. A recent generalization of the Lindstedt-Poincare (L-P) method is applied to approximate a homoclinic solution of a strongly nonlinear autonomous system that results from blowing up the BT normal form. This solution allows us to derive an accurate second-order homoclinic predictor to the homoclinic branch rooted at a generic BT point of an n-dimensional ordinary differential equation (ODE). We prove that the method leads to the same homoclinicity conditions as the classical Melnikov technique, the branching method, and the regular perturbation (R-P) method. However, it is known that the R-P method leads to a `parasitic turn' near the saddle point. The new asymptotics based on the L-P method do not have this turn, making them more suitable for numerical implementation. We show how to use these asymptotics to calculate the initial data to continue homoclinic orbits in two free parameters. The new homoclinic predictors are implemented in the MATLAB continuation package MatCont to initialize the continuation of homoclinic orbits from a BT point. Two examples with multidimensional state spaces are included.

    DOI
  11. Accurate Approximation of Homoclinic Solutions in Gray-Scott Kinetic Model, Yu.A. Kuznetsov, H.G.E. Meijer, B. Al-Hdaibat and W. Govaerts, Internat. J. Bifur. Chaos Appl. Sci. Engrg. 2015, Vol. 25, No. 9: 1550125 Abstract
    Description:

    The second-order predictor for the homoclinic orbit is applied to the Gray-Scott model. The problem is used to illustrate the approximation of the homoclinic orbits near a generic Bogdanov-Takens bifurcation in n-dimensional systems of differential equations. In the process we show that it is necessary to take (usually ignored) cubic terms in the Bogdanov-Takens normal form into account to derive a correct second-order prediction for the homoclinic bifurcation curve. The analytic solutions are compared with those obtained by numerical continuation.

    DOI
  12. Computational modeling of Adelta-fiber-mediated nociceptive detection of electrocutaneous stimulation. H. Yang, H.G.E. Meijer, R.J. Doll, J.R. Buitenweg, S.A. van Gils, Biological Cybernetics 2015, Vol. 109, No. 4: 479-491 Abstract
    Description:

    Sensitization is an example of malfunctioning of the nociceptive pathway in either the peripheral or central nervous system. Using quantitative sensory testing, one can only infer sensitization, but not determine the defective subsystem. The states of the subsystems may be characterized using computational modeling together with experimental data. Here, we develop a neurophysiologically plausible model replicating experimental observations from a psychophysical human subject study. We study the effects of single temporal stimulus parameters on detection thresholds corresponding to a 0.5 detection probability. To model peripheral activation and central processing, we adapt a stochastic drift-diffusion model and a probabilistic hazard model to our experimental setting without reaction times. We retain six lumped parameters in both models characterizing peripheral and central mechanisms. Both models have similar psychophysical functions, but the hazard model is computationally more efficient. The model-based effects of temporal stimulus parameters on detection thresholds are consistent with those from human subject data.

    DOI (open access)
  13. Modeling Focal Epileptic Activity in the Wilson-Cowan model with Depolarization Block. H.G.E. Meijer, T.L. Eissa, B. Kiewiet, J.F. Neuman, C.A. Schevon, R.G. Emerson, R.R. Goodman, G.M. McKhann Jr, C.J. Marcuccilli, A.K. Tryba, J.D. Cowan, S.A. van Gils and W. van Drongelen, Journal of Mathematical Neuroscience 2015, 5:7 Abstract
    Description:

    Measurements of neuronal signals during human seizure activity and evoked epileptic activity in experimental models suggest that, in these pathological states, the individual nerve cells experience an activity driven depolarization block, i.e. they saturate. We examined the effect of such a saturation in the Wilson-Cowan formalism by adapting the nonlinear activation function; we substituted the commonly applied sigmoid for a Gaussian function. We discuss experimental recordings during a seizure that support this substitution. Next we perform a bifurcation analysis on the Wilson-Cowan model with a Gaussian activation function. The main effect is an additional stable equilibrium with high excitatory and low inhibitory activity. Analysis of coupled local networks then shows that such high activity can stay localized or spread. Specifically, in a spatial continuum we show a wavefront with inhibition leading followed by excitatory activity. We relate our model simulations to observations of spreading activity during seizures.

    DOI PDF
  14. Tracking pattern evolution through extended center manifold reduction and singular perturbations. L. Sewalt, A. Doelman, H.G.E. Meijer, V. Rottschfer, A. Zagaris, Physica D 2015, Vol. 298-299, pp. 48-67 Abstract
    Description:

    In this paper we develop an extended center manifold reduction method: a methodology to analyze the formation and bifurcations of small-amplitude patterns in certain classes of multi-component, singularly perturbed systems of partial differential equations. We specifically consider systems with a spatially homogeneous state whose stability spectrum partitions into eigenvalue groups with distinct asymptotic properties. One group of successive eigenvalues in the bifurcating group are widely interspaced, while the eigenvalues in the other are stable and cluster asymptotically close to the origin along the stable semi-axis. The classical center manifold reduction provides a rigorous framework to analyze destabilizations of the trivial state, as long as there is a spectral gap of sufficient width. When the bifurcating eigenvalue becomes commensurate to the stable eigenvalues clustering close to the origin, the center manifold reduction breaks down. Moreover, it cannot capture subsequent bifurcations of the bifurcating pattern. Through our methodology, we formally derive expressions for low-dimensional manifolds exponentially attracting the full flow for parameter combinations that go beyond those allowed for the (classical) center manifold reduction, i.e. to cases in which the spectral gap condition no longer can be satisfied. Our method also includes an explicit description of the flow on these manifolds and thus provides an analytical tool to study subsequent bifurcations. Our analysis centers around primary bifurcations of transcritical type -- that can be either of co-dimension 1 or 2 -- in two- and three-component PDE systems. We employ our method to study bifurcation scenarios of small-amplitude patterns and the possible appearance of low-dimensional spatio-temporal chaos. We also exemplify our analysis by a number of characteristic reaction-diffusion systems with disparate diffusivities.

    DOI PDF
  15. Exploiting pallidal plasticity for stimulation in Parkinson's disease, M.A.J. Lourens, B. Schwab, J.A. Nirody, H.G.E. Meijer and S.A. van Gils, Journal of Neural Engineering, 2015, 12: 026005 Abstract
    Description:

    Objective. Continuous application of high-frequency deep brain stimulation (DBS) often effectively reduces motor symptoms of Parkinson's disease patients. While there is a growing need for more effective and less traumatic stimulation, the exact mechanism of DBS is still unknown. Here, we present a methodology to exploit the plasticity of GABAergic synapses inside the external globus pallidus (GPe) for the optimization of DBS. Approach. Assuming the existence of spike-timing-dependent plasticity (STDP) at GABAergic GPe-GPe synapses, we simulate neural activity in a network model of the subthalamic nucleus and GPe. In particular, we test different DBS protocols in our model and quantify their influence on neural synchrony. Main results. In an exemplary set of biologically plausible model parameters, we show that STDP in the GPe has a direct influence on neural activity and especially the stability of firing patterns. STDP stabilizes both uncorrelated firing in the healthy state and correlated firing in the parkinsonian state. Alternative stimulation protocols such as coordinated reset stimulation can clearly profit from the stabilizing effect of STDP. These results are widely independent of the STDP learning rule. Significance. Once the model settings, e.g., connection architectures, have been described experimentally, our model can be adjusted and directly applied in the development of novel stimulation protocols. More efficient stimulation leads to both minimization of side effects and savings in battery power.

    DOI
  16. Bifurcation analysis of a model for atherosclerotic plaque evolution. M.A.K. Bulelzai, J.L.A. Dubbeldam, H.G.E. Meijer, Physica D 2014, Vol. 278-9, pp. 31-43 Abstract
    Description:

    We analyze two ordinary differential equation (ODE) models for atherosclerosis. The ODE models describe long time evolution of plaques in arteries. We show how the dynamics of the first atherosclerosis model (model A) can be understood using codimension-two bifurcation analysis. The Low-Density Lipoprotein (LDL) intake parameter (d) is the first control parameter and the second control parameter is either taken to be the conversion rate of macrophages (b) or the wall shear stress (σ). Our analysis reveals that in both cases a Bogdanov-Takens (BT) point acts as an organizing center. The bifurcation diagrams are calculated partly analytically and to a large extent numerically using AUTO07 and MATCONT. The bifurcation curves show that the concentration of LDL in the plaque as well as the monocyte and the macrophage concentration exhibit oscillations for a certain range of values of the control parameters. Moreover, we find that there are threshold values for both the cholesterol intake rate dcritdcrit and the conversion rate of the macrophages bcritbcrit, which depend on the values of other parameters, above which the plaque volume increases with time. It is found that larger conversion rates of macrophages lower the threshold value of cholesterol intake and vice versa. We further argue that the dynamics for model A can still be discerned in the second model (model B) in which the slow evolution of the radius of the artery is coupled self-consistently to changes in the plaque volume. The very slow evolution of the radius of the artery compared to the other processes makes it possible to use a slow manifold approximation to study the dynamics in this case. We find that in this case the model predicts that the concentrations of the plaque constituents may go through a period of oscillations before the radius of the artery will start to decrease. These oscillations hence act as a precursor for the reduction of the artery radius by plaque growth.

    DOI PDF
  17. Improved homoclinic predictor for Bogdanov-Takens bifurcation. Yu.A. Kuznetsov, H.G.E. Meijer, B. Al Hdaibat and W. Govaerts, Internat. J. Bifur. Chaos Appl. Sci. Engrg., 2014, Vol. 24, No. 4: 1450057 Abstract
    Description:

    An improved homoclinic predictor at a generic codim 2 Bogdanov-Takens (BT) bifucation is derived. We use the classical "blow-up" technique to reduce the canonical smooth normal form near a generic BT bifurcation to a perturbed Hamiltonian system. With a simple perturbation method, we not only rederive the known approximation for the bifurcation parameter values but also obtain an explicit first-order correction of the unperturbed homoclinic orbit. To obtain the normal form on the center manifold, we apply the standard parameter-dependent center manifold reduction combined with the normalization, that is based on the Fredholm solvability of the homological equation. However, by systematically solving all linear systems appearing from the homological equation, we remove an ambiguity in the parameter transformation existing in the literature. The actual implementation of the improved predictor in MatCont and numerical examples illustrating its effciency are discussed.

    DOI
  18. Travelling waves in models of neural tissue: from localised structures to periodic waves, H.G.E. Meijer and S. Coombes, EPJ Nonlinear Biomedical Physics 2014, 2:3 Abstract
    Description:

    We consider travelling waves (fronts, pulses and periodics) in spatially extended one dimensional neural field models. We demonstrate for an excitatory field with linear adaptation that, in addition to an expected stable pulse solution, a stable anti-pulse can exist. Varying the adaptation strength we unravel the organizing centers of the bifurcation diagram for fronts and pulses, with a mixture of exact analysis for a Heaviside firing rate function and novel numerical schemes otherwise. These schemes, for non-local models with space-dependent delays, further allow for the construction and continuation of periodic waves. We use them to construct the dispersion curve -- wave speed as a function of period -- and find that they can be oscillatory and multi-valued, suggesting bistability of periodic waves. A kinematic theory predicts the onset of wave instabilities at stationary points in the dispersion curve, leading to period doubling behaviour, and is confirmed with direct numerical simulations. We end with a discussion of how the construction of dispersion curves may allow a useful classification scheme of neural field models for epileptic waves.

    DOI (open access) Poster version (CNS, 2013 Paris) see DOI . Follow the link for documented simulation and continuation Code . PDF
  19. Travelling waves in a neural field model with refractoriness, H.G.E. Meijer and S. Coombes, Journal of Mathematical Biology 2014, Vol. 68, No. 5, pp. 1249-1268 Abstract
    Description:

    At one level of abstraction neural tissue can be regarded as a medium for turning local synaptic activity into output signals that propagate over large distances via axons to generate further synaptic activity that can cause reverberant activity in networks that possess a mixture of excitatory and inhibitory connections. This output is often taken to be a firing rate, and the mathematical form for the evolution equation of activity depends upon a spatial convolution of this rate with a fixed anatomical connectivity pattern. Such formulations often neglect the metabolic processes that would ultimately limit synaptic activity. Here we reinstate such a process, in the spirit of an original prescription by Wilson and Cowan [1972, ``Excitatory and inhibitory interactions in localised populations of model neurons'', Biophysical Journal 12, 1-24.], using a term that multiplies the usual spatial convolution with a moving time average of local activity over some refractory time-scale. This modulation can substantially affect network behaviour, and in particular give rise to periodic travelling waves in a purely excitatory network (with exponentially decaying anatomical connectivity), which in the absence of refractoriness would only support travelling fronts. We construct these solutions numerically as stationary periodic solutions in a co-moving frame (of both an equivalent delay differential model as well as the original delay integro-differential model). Continuation methods are used to obtain the dispersion curve for periodic travelling waves (speed as a function of period), and found to be reminiscent of those for spatially extended models of excitable tissue. A kinematic analysis (based on the dispersion curve) predicts the onset of wave instabilities, which are confirmed numerically.

    DOI (open access) PDF
  20. Tracking of nociceptive thresholds using adaptive psychophysical methods, R.J. Doll, J.R. Buitenweg, H.G.E. Meijer and P. Veltink, Behavior Research Methods 2014, Vol. 46, No. 1, pp. 55-66 Abstract
    Description:
    Tracking psychophysical thresholds over time is beneficial for the investigation of dynamic behavior of underlying mechanisms. For example, noxious events can activate endogenous analgesic mechanisms reflected in an increased nociceptive threshold. Psychophysical techniques for tracking these thresholds should include efficient methods for stimulus selection and threshold estimation. This study compares the performance of different adaptive stimulus selection procedures and estimation methods in simulation and human psychophysical experiments.
    Monte Carlo simulations were performed to compare the bias and precision of threshold estimates between three stimulus selection procedures (simple staircase, random staircase and a minimum entropy procedure) and two estimation methods (logistic regression and Bayesian estimation). Logistic regression and Bayesian estimations resulted in similar precision only when prior probability distributions were appropriately chosen. The minimum entropy (with an appropriate prior distribution) and simple staircase procedures showed the highest precision but the random staircase procedure was less sensitive to different procedure specific settings than the other procedures.
    The simple staircase and random staircase stimulus selection procedures were compared in a human subject study (n=30), using logistic regression for threshold estimation. Electrocutaneous stimulation was used to track the nociceptive perception threshold before, during and after a cold pressor task as a conditioning stimulus. The two stimulus selection procedures were able to detect habituation as well as changes induced by the conditioning stimulus, but the random staircase procedure achievedthe highest precision.
    Based on these results, we recommend the random staircase procedure in combination with logistic regression for nonstationary threshold tracking experiments.

    DOI
  21. Analysis of bifurcations of limit cycles with Lyapunov exponents and numerical normal forms, V. De Witte, W. Govaerts, Yu.A. Kuznetsov, H.G.E. Meijer, Physica D 2014, Vol. 269, pp. 126-141 Abstract
    Description:

    Explicit computational formulas for coefficients of the periodic normal forms of the three most complex codim 2 bifurcations of limit cycles with dimension of the center manifold equal to 4 or to 5 in generic autonomous ODEs are derived. The resulting formulas are independent of the dimension of the phase space and involve solutions of certain boundary-value problems as well as multilinear functions from the Taylor expansion of the ODE right-hand side near the cycle. The formulas allow one to distinguish between the complicated bifurcation scenarios which can happen near these codim 2 bifurcations of limit cycles, where 3-tori and 4-tori can be present. We apply our techniques to the study of a known laser model, a novel model from population biology, and one for mechanical vibrations; these models exhibit Limit Point--Neimark-Sacker, Period-Doubling--Neimark-Sacker and double Neimark-Sacker bifurcations. Lyapunov exponents are computed to numerically confirm the results of the normal form analysis, in particular with respect to the existence of stable invariant tori of various dimensions and chaos.

    DOI PDF
  22. Epileptic spikes are phase-locked to rhythmic delta activity in non-convulsive status epilepticus. R. Hindriks, H.G.E. Meijer, S.A. van Gils and M.J.A.M van Putten, Frontiers in Systems Neuroscience 2013, Vol. 7: 111 Abstract
    Description:

    The EEG of patients in non-convulsive status epilepticus (NCSE) often displays delta oscillations or generalized spike-wave discharges. In some patients, these delta oscillations coexist with intermittent epileptic spikes. In this study we verify the prediction of a computational model of the thalamo-cortical system that these spikes are phase-locked to the delta oscillations. We subsequently describe the physiological mechanism underlying this observation as suggested by the model. It is suggested that the spikes reflect inhibitory stochastic fluctuations in the input to thalamo-cortical relay neurons and phase-locking is a consequence of differential excitability of relay neurons over the delta cycle. Further analysis shows that the observed phase-locking can be regarded as a stochastic precursor of generalized spike-wave discharges. This study thus provides an explanation of intermittent spikes during delta oscillations in NCSE and might be generalized to other encephathologies in which delta activity can be observed.

    DOI
  23. Functional neuronal activity and connectivity within the subthalamic nucleus in Parkinson's disease, M.A.J. Lourens, H.G.E. Meijer, M.F. Contarino, P. van den Munckhof, P.R. Schuurman, S.A. van Gils and L.J. Bour, Clinical Neurophysiology 2013, Vol. 124, No. 5, pp 967-981 Abstract
    Description:

    Objective: Characterization of the functional neuronal activity and connectivity within the subthalamic nucleus (STN) in patients with Parkinson's disease (PD).
    Methods: Single units were extracted from micro-electrode recording (MER) of 18 PD patients who underwent STN deep brain stimulation (DBS) surgery. The firing rate and pattern of simultaneously recorded spike trains and their coherence were analyzed. To provide a precise functional assignment of position to the observed activities, for each patient we mapped its classified multichannel STN MERs to a generic atlas representation with a sensorimotor part and a remaining part.
    Results: Within the sensorimotor part we found significantly higher mean firing rate (P < 0.05) and significantly more burst-like activity (P < 0.05) than within the remaining part. The proportion of significant coherence in the beta band (13-30 Hz) is significantly higher in the sensorimotor part of the STN than elsewhere (P=0.015).
    Conclusions: The STN sensorimotor part distinguishes itself from the remaining part with respect to beta coherence, firing rate and burst-like activity and postoperatively was found as the preferred target area. Significance: Our firing behavior analysis may help to discriminate the STN sensorimotor part for the placement of the DBS electrode.

    DOI
  24. Analysis of Stability and Bifurcations of Fixed Points and Periodic Solutions of a Lumped Model of Neocortex with Two Delays. S. Visser S, H.G.E. Meijer, M.J.A.M. van Putten and S.A. van Gils, The Journal of Mathematical Neuroscience 2012, 2:8 Abstract
    Description:

    A lumped model of neural activity in neocortex is studied to identify regions of multi-stability of both steady states and periodic solutions. Presence of both steady states and periodic solutions is considered to correspond with epileptogenesis. The model, which consists of two delay differential equations with two fixed time lags, is mainly studied for its dependency on varying connection strength between populations. Equilibria are identified, and using linear stability analysis, all transitions are determined under which both trivial and non-trivial fixed points lose stability. Periodic solutions arising at some of these bifurcations are numerically studied with a two-parameter bifurcation analysis.

    DOI
  25. The Hopf-van der Pol System: Failure of a Homotopy Method, H.G.E. Meijer and T. Kalmar-Nagy, Differential Equations and Dynamical Systems: Vol. 20, No. 3, pp 323-328, 2012 Abstract
    Description:

    The purpose of this article to provide an explicit example where continuation based on the homotopy method fails. The example is a one-parameter homotopy for periodic orbits between two well-known nonlinear systems, the normal form of the Hopf bifurcation and the van der Pol system. Our analysis shows that various types of obstructions can make approximation over the whole range of the homotopy parameter impossible. The Hopf-van der Pol system demonstrates that homotopy methods may fail even for seemingly innocent systems.

    DOI
  26. From parkinsonian thalamic activity to restoring thalamic relay using Deep Brain Stimulation: new insights from computational modeling, H.G.E. Meijer, M. Krupa, H. Cagnan, M.A.J. Lourens, T. Heida, H.C.F. Martens, L.J. Bour and S.A. van Gils, Journal of Neural Engineering, 2011, 8: 066005 Abstract
    Description:

    We present a computational model of a thalamocortical relay neuron for exploring basal ganglia thalamocortical loop behavior in relation to Parkinson's disease and deep brain stimulation (DBS). Previous microelectrode, single-unit recording studies demonstrated that oscillatory interaction within and between basal ganglia nuclei is very often accompanied by synchronization at parkinsonian rest tremor frequencies (3-10 Hz). These oscillations have a profound influence on thalamic projections and impair the thalamic relaying of cortical input by generating rebound action potentials. Our model describes convergent inhibitory input received from basal ganglia by the thalamocortical cells based on characteristics of normal activity, and/or low-frequency oscillations (activity associated with Parkinson's disease). In addition to simulated input, we also used microelectrode recordings as inputs for the model. In the resting state, and without additional sensorimotor input, pathological rebound activity is generated for even mild parkinsonian input. We have found a specific stimulation window of amplitudes and frequencies for periodic input, which corresponds to highfrequency DBS, and which also suppresses rebound activity for mild and even more prominent parkinsonian input. When low-frequency pathological rebound activity disables the thalamocortical cell's ability to relay excitatory cortical input, a stimulation signal with parameter settings corresponding to our stimulation window can restore the thalamocortical cell's relay functionality.

    DOI PDF
  27. The pedunculopontine nucleus as an additional target for deep brain stimulation. M.A.J. Lourens, H.G.E. Meijer, T. Heida, E. Marani and S.A. van Gils, Neural networks, 2011, Vol. 24, No. 6, pp 617-630 Abstract
    Description:

    The pedunculopontine nucleus has been suggested as a target for DBS. In this paper we propose a single compartment computational model for a PPN Type I cell and compare its dynamic behavior with experimental data. The model shows bursts after a period of hyperpolarization and spontaneous firing at 8 Hz. Bifurcation analysis of the single PPN cell shows bistability of fast and slow spiking solutions for a range of applied currents. A network model for STN, GPe and GPi produces basal ganglia output that is used as input for the PPN cell. The conductances for projections from the STN and the GPi to the PPN are determined from experimental data. The resulting behavior of the PPN cell is studied under normal and Parkinsonian conditions of the basal ganglia network. The effect of high frequency stimulation of the STN is considered as well as the effect of combined high frequency stimulation of the STN and the PPN at various frequencies. The relay properties of the PPN cell demonstrate that the combined high frequency stimulation of STN and low frequency (10 Hz, 25 Hz, 40 Hz) stimulation of PPN hardly improves the effect of exclusive STN stimulation. Moreover, PPN-DBS at low stimulation amplitude has a better effect than at higher stimulation amplitude. The effect of PPN output on the basal ganglia is investigated, in particular the effect of STN-DBS and/or PPN-DBS on the pathological firing pattern of STN and GPe cells. PPN-DBS eliminates the pathological firing pattern of STN and GPe cells, whereas STN-DBS and combined STN-DBS and PPN-DBS eliminate the pathological firing pattern only from STN cells.

    DOI
  28. Comparing Epileptiform Behavior of Mesoscale Detailed Models and Population Models of Neocortex. S. Visser, H.G.E. Meijer, H.C. Lee, W. van Drongelen, M.J.A.M. van Putten and S.A. van Gils, Journal of Clinical Neurophysiology, 2010, Vol. 27, No. 6, pp 471-478 Abstract
    Description:

    Two models of the neocortex are developed to study normal and pathologic neuronal activity. One model contains a detailed description of a neocortical microcolumn represented by 656 neurons, including superficial and deep pyramidal cells, four types of inhibitory neurons, and realistic synaptic contacts. Simulations show that neurons of a given type exhibit similar, synchronized behavior in this detailed model. This observation is captured by a population model that describes the activity of large neuronal populations with two differential equations with two delays. Both models appear to have similar sensitivity to variations of total network excitation. Analysis of the population model reveals the presence of multistability, which was also observed in various simulations of the detailed model.

    DOI
  29. Frequency-selectivity of a thalamocortical relay neuron during Parkinson's disease and deep brain stimulation: a computational study. H. Cagnan, H.G.E. Meijer, S.A. van Gils, M. Krupa, T. Heida, M. Rudolph, W.J. Wadman, and H.C.F. Martens, European Journal of Neuroscience, 2009, Vol. 30, No. 7, pp 1306-1317 Abstract
    Description:

    In this computational study, we investigated (i) the functional importance of correlated basal ganglia (BG) activity associated with Parkinson's disease (PD) motor symptoms by analysing the effects of globus pallidus internum (GPi) bursting frequency and synchrony on a thalamocortical (TC) relay neuron, which received GABAergic projections from this nucleus; (ii) the effects of subthalamic nucleus (STN) deep brain stimulation (DBS) on the response of the TC relay neuron to synchronized GPi oscillations; and (iii) the functional basis of the inverse relationship that has been reported between DBS frequency and stimulus amplitude, required to alleviate PD motor symptoms [A. L. Benabid et al. (1991)Lancet, 337, 403-406]. The TC relay neuron selectively responded to and relayed synchronized GPi inputs bursting at a frequency located in the range 2-25 Hz. Input selectivity of the TC relay neuron is dictated by low-threshold calcium current dynamics and passive membrane properties of the neuron. STN-DBS prevented the TC relay neuron from relaying synchronized GPi oscillations to cortex. Our model indicates that DBS alters BG output and input selectivity of the TC relay neuron, providing an explanation for the clinically observed inverse relationship between DBS frequency and stimulus amplitude.

    DOI
  30. Emergence and bifurcations of slow manifolds in nonlinear wave equations. T. Bakri, H.G.E. Meijer and F. Verhulst, Journal of Nonlinear Science, 2009, Vol. 19, No. 5, pp 571-596 Abstract
    Description:

    Persistence and bifurcations of Lyapunov manifolds can be studied by a combination of averaging-normalization and numerical bifurcation methods. This can be extended to infinite-dimensional cases when using suitable averaging theorems. The theory is applied to the case of a parametrically excited wave equation. We find fast dynamics in a finite, resonant part of the spectrum and slow dynamics elsewhere. The resonant part corresponds with an almost-invariant manifold and displays bifurcations into a wide variety of phenomena among which are 2- and 3-tori.

    DOI (open access) PDF
  31. Numerical Bifurcation Analysis. H.G.E. Meijer, F. Dercole and B. Oldeman, in Meyers, Robert (Ed.) Encyclopedia of Complexity and Systems Science, Part 14, Pages 6329-6352, Springer New York, 2009. Abstract
    Description:

    The theory of dynamical systems studies the behavior of solutions of systems, like nonlinear ordinary differential equations (ODEs), depending upon parameters. Using qualitative methods of bifurcation theory, the behavior of the system is characterized for various parameter combinations. In particular, the catalog of system behaviors showing qualitative differences can be identified, together with the regions in parameter space where the different behaviors occur. Bifurcations delimit such regions. Symbolic and analytical approaches are in general infeasible, but numerical bifurcation analysis is a powerful tool that aids in the understanding of a nonlinear system. When computing power became widely available, algorithms for this type of analysis matured and the first codes were developed. With the development of suitable algorithms, the advancement in the qualitative theory has found its way into several software projects evolving over time. The availability of software packages allows scientists to study and adjust their models and to draw conclusions about their dynamics.

    DOI PDF
  32. Numerical continuation of connecting orbits of maps. W. Govaerts, R. Khoshsiar Ghaziani, Yu.A. Kuznetsov, and H.G.E. Meijer, Journal of Difference Equations and Applications, 2009, Vol. 15, No. 8, pp 849-875 Abstract
    Description:

    We present new or improved methods to continue heteroclinic and homoclinic orbits to fixed points in iterated maps and to compute their fold bifurcation curves, corresponding to the tangency of the invariant manifolds. The proposed methods are applicable to general n-dimensional maps and are implemented in matlab. They are based on the continuation of invariant subspaces (CIS) algorithm, which is presented in a novel way. The systems of defining equations include the Riccati equations appearing in CIS for bases of the generalized stable and unstable eigenspaces. We use the bordering techniques to continue the folds, and provide full algorithmic details on how to treat the Jacobian matrix of the defining system as a sparse matrix in matlab. For a special -- but important in applications -- case n=2 we describe the first matlab implementation of known algorithms to grow one-dimensional stable and unstable manifolds of the fixed points of noninvertible maps. The methods are applied to study heteroclinic and homoclinic connections in the generalized Hnon map.

    DOI
  33. Switching to nonhyperbolic cycles from codim 2 bifurcations of equilibria in ODEs. Yu.A. Kuznetsov, H.G.E. Meijer, W. Govaerts and B. Sautois, Physica D 2008, Vol. 237, pp 3061-3068. Abstract
    Description:

    The paper provides full algorithmic details on switching to the continuation of all possible codim 1 cycle bifurcations from generic codim 2 equilibrium bifurcation points in n-dimensional ODEs. We discuss the implementation and the performance of the algorithm in several examples, including an extended Lorenz-84 model and a laser system.

    DOI
  34. New features of the software MatCont for bifurcation analysis of dynamical systems. A. Dhooge, and W. Govaerts, Yu.A. Kuznetsov, H.G.E. Meijer and B. Sautois, Mathematical and Computer Modelling of Dynamical Systems 2008, Vol. 14, No. 2, pp 147-175 Abstract
    Description:

    Bifurcation software is an essential tool in the study of dynamical systems. From the beginning (the first packages were written in the 1970's) it was also used in the modelling process, in particular to determine the values of critical parameters. More recently, it is used in a systematic way in the design of dynamical models and to determine which parameters are relevant. MatCont and Cl_MatCont are freely available matlab numerical continuation packages for the interactive study of dynamical systems and bifurcations. MatCont is the GUI-version, Cl_MatCont is the command-line version. The work started in 2000 and the first publications appeared in 2003. Since that time many new functionalities were added. Some of these are fairly simple but were never before implemented in continuation codes, e.g. Poincar maps. Others were only available as toolboxes that can be used by experts, e.g. continuation of homoclinic orbits. Several others were never implemented at all, such as periodic normal forms for codimension 1 bifurcations of limit cycles, normal forms for codimension 2 bifurcations of equilibria, detection of codimension 2 bifurcations of limit cycles, automatic computation of phase response curves and their derivatives, continuation of branch points of equilibria and limit cycles. New numerical algorithms for these computations have been published or will appear elsewhere; here we restrict to their software implementation. We further discuss software issues that are in practice important for many users, e.g. how to define a new system starting from an existing one, how to import and export data, system descriptions, and computed results.

    DOI / PDF .
  35. Numerical Methods for Two-Parameter Local Bifurcation Analysis of Maps. W. Govaerts, R. Khoshsiar Ghaziani, Yu.A. Kuznetsov and H.G.E. Meijer, SIAM Journal on Scientific Computing 2007, Vol. 29, No. 6, pp 2644-2667 Abstract
    Description:

    We discuss new and improved algorithms for the bifurcation analysis of fixed points and periodic orbits (cycles) of maps and their implementation in matcont, a MATLAB toolbox for continuation and bifurcation analysis of dynamical systems. This includes the numerical continuation of fixed points of iterates of the map with one control parameter, detecting and locating their bifurcation points (i.e., limit point, period-doubling, and Neimark-Sacker) and their continuation in two control parameters, as well as detection and location of all codimension 2 bifurcation points on the corresponding curves. For all bifurcations of codim 1 and 2, the critical normal form coefficients are computed, both numerically with finite directional differences and using symbolic derivatives of the original map. Using a parameter-dependent center manifold reduction, explicit asymptotics are derived for bifurcation curves of double and quadruple period cycles rooted at codim 2 points of cycles with arbitrary period. These asymptotics are implemented into the software and allow one to switch at codim 2 points to the continuation of the double and quadruple period bifurcations. We provide two examples illustrating the developed techniques: a generalized Henon map and a juvenile/adult competition model from mathematical biology.

    DOI
  36. Remarks on interacting Neimark-Sacker bifurcations. Yu.A. Kuznetsov and H.G.E. Meijer, Journal of Difference Equations and Applications 2006, Vol. 12, No. 10, pp 1009-1035, 2006 Abstract
    Description:

    We study codimension-2 bifurcations of fixed points of dissipative diffeomorphisms with a pair of complex eigenvalues together with either an eigenvalue $-1$ or another such a pair. In the previous studies only cubic normal forms were considered. However, in some cases the unfolding requires higher order terms and these are investigated here. We (re)derive the normal forms and reduce them to a single amplitude map. This map is similar to the amplitude system for the double-Hopf bifurcation for vector fields. We show how the critical normal form coefficients determine the general bifurcation picture for this amplitude map. Representative nonsymmetric perturbations of the normal forms are studied numerically. Our case studies show a detailed picture near various bifurcation curves, which is richer than known theoretical predictions. For arbitrary maps with these bifurcations we give explicit formulas for critical normal form coefficients on center manifolds. We apply them to an example from robotics where we are able to demonstrate the existence of a bubble-structure, which was only observed in perturbed normal forms before.

    DOI
  37. Generalized Hnon map and bifurcations of homoclinic tangencies. V.S. Gonchenko, Yu.A. Kuznetsov and H.G.E. Meijer, SIAM Journal on Applied Dynamical Systems 2005, Vol. 4, No. 2, pp 407-436 Abstract
    Description:

    The Generalized Hnon Map(GHM) appears in codim 2 bifurcation of homoclinic orbits for maps, for example a neutral saddle with a homoclinic tangency (NHT). The standard Hnon map with the standard parameters itself is interesting for the strange attractor, the homoclinic orbit and the resulting chaos. However, from the bifurcation point of view it turns out that many appearing bifurcations of the standard Hnon map are degenerate, for example, the strong resonances. This is not the case for GHM. One interesting feature is that the GHM has the NHT-bifurcation, so we show that the existence of infinitely many of such bifurcations.

    DOI
  38. Numerical normal forms for codim 2 bifurcations of fixed points with at most two critical eigenvalues. Yu.A. Kuznetsov and H.G.E Meijer, SIAM Journal on Scientific Computing 2005, Vol. 26, No. 6, pp 1932-1954 Abstract
    Description:

    We derive the formulas for the center manifold reduction and the critical normal form coefficients for the codimension two bifurcations for diffeomorphisms of the plane. Some formulas involve third or fifth order derivatives. This poses no problems for maps, but for limit cycles of ODE's it does. A RungeKuttaFehlberg-78 integrator adapted with automatic differentiation(ADOL-C 1.8) is used to obtain the required data.

    DOI
  39. The Fold-Flip bifurcation. Yu.A. Kuznetsov, H.G.E Meijer and L. van Veen, Internat. J. Bifur. Chaos Appl. Sci. Engrg. 2004, Vol. 14, No. 7, pp 2253-2282 Abstract
    Description:

    The fold-flip bifurcation occurs if a map has a fixed point with multipliers +1 and -1 simultaneously. In this paper the normal form of this singularity is calculated explicitly. Both local and global bifurcations of the unfolding are analysed by exploring a close relationship between the derived normal form and the truncated amplitude system for the fold-Hopf bifurcation of ODEs. Two examples are presented, the generalized H\'enon map and an extension of the Lorenz-84 model. In the latter example the first, second and third order derivatives of the Poincar\'e map are computed to find the normal form coefficients.

    DOI

Preprints/Conference Proceedings/Posters

Ph.D Thesis

Title: Codimension 2 Bifurcations of Iterated Maps
Promotor: prof.dr. F. Verhulst
Co-Promotor: dr. Yu.A Kuznetsov
An electronic version is available here . Or mail me for a printed copy.

Master Thesis