Note

Click here to download the full example code

# 2.1. The Jansen-Rit Neural Mass Model¶

Here, we will introduce the Jansen-Rit model, a neural mass model of the dynamic interactions between 3 populations:

pyramidal cells (PCs)

excitatory interneurons (EINs)

inhibitory interneurons (IINs)

Originally, the model has been developed to describe the waxing-and-waning of EEG activity in the alpha frequency range (8-12 Hz) in the visual cortex 1. In the past years, however, it has been used as a generic model to describe the macroscopic electrophysiological activity within a cortical column 2. A graphic representation of such a model, placed inside a brain network, can be found in the figure below.

The structure of the full Jansen-Rit model is depicted in Figure 1 B. As visualized for the pyramidal cell population in :ref:`fig1`C, the model can be decomposed into a number of generic mathematical operators that can be used to build the dynamic equations for the Jansen-Rit model. Essentially, the membrane potential deflections that are caused at the somata of a population by synaptic input, are modeled by a convolution operation with an alpha kernel. This choice has been shown to reflect the dynamic process of polarization propagation from the synapse via the dendritic tree to the soma 3. The convolution operation can be expressed via a second-order differential equation:

where represents the average post-synaptic potential and and are the efficacy and the timescale of the synapse, respectively. As a second operator, the translation of the average membrane potential deflection at the soma to the average firing of the population is given by a sigmoidal function:

In this equation, and represent the average firing rate and membrane potential, respectively, while , and are constants defining the maximum firing rate, firing threshold variance and average firing threshold within the modeled population, respectively.

By using the linearity of the convolution operation, the dynamic interactions between PCs, EINs and IINs can be expressed via 6 coupled ordinary differential equations that are composed of the two operators defined above:

where , , are used to represent the average membrane potential deflection caused by the excitatory synapses at the PC population, the inhibitory synapses at the PC population, and the excitatory synapses at both interneuron populations, respectively.

Below, we will demonstrate how to load this a model into pyrates and perform numerical simulations with it.

## 2.1.1. References¶

- 1
B.H. Jansen & V.G. Rit (1995)

*Electroencephalogram and visual evoked potential generation in a mathematical model of coupled cortical columns.*Biological Cybernetics, 73(4): 357-366.- 2
A. Spiegler, S.J. Kiebel, F.M. Atay, T.R. Knösche (2010)

*Bifurcation analysis of neural mass models: Impact of extrinsic inputs and dendritic time constants.*NeuroImage, 52(3): 1041-1058, https://doi.org/10.1016/j.neuroimage.2009.12.081.- 3
P.A. Robinson, C.J. Rennie, J.J. Wright (1997)

*Propagation and stability of waves of electrical activity in the cerebral cortex.*Physical Review E, 56(826), https://doi.org/10.1103/PhysRevE.56.826.

### 2.1.1.1. Step 1: Importing the frontend class for defining models¶

As a first step, we import the `pyrates.frontend.CircuitTemplate`

class, which allows us to set up a model
definition in PyRates.

```
from pyrates.frontend import CircuitTemplate
```

### 2.1.1.2. Step 2: Loading a model template from the model_templates library¶

In the second step, we load a model template for the Jansen-Rit model that comes with PyRates via the
`from_yaml()`

method of the `CircuitTemplate`

. This method returns a `CircuitTemplate`

instance.
Have a look at the yaml definition of the model that can be found via the path used for the `from_yaml()`

method. You will see that all variables and parameters are already defined there. As an alternative, you can also use
the `circuit_from_yaml()`

function that you can simply import via `from pyrates import circuit_from_yaml`

.
These are the basic steps you perform, if you want to load a model that is defined inside a yaml file.
To check out the different model templates provided by PyRates, have a look at the `PyRates.model_templates`

module.

```
jrc = CircuitTemplate.from_yaml("model_templates.neural_mass_models.jansenrit.JRC")
```

### 2.1.1.3. Step 3: Numerical simulation of a the model behavior in time¶

After loading the model template, numerical simulations can be performed via the `run()`

method.
Calling this function will solve the initial value problem of the above defined differential equations for a time
interval from 0 to the given simulation time.
This solution will be calculated numerically by a differential equation solver in the backend, starting with a defined
step-size. As part of this process, a `pyrates.backend.BaseBackend`

instance is created that contains a graph
representation of the network equations. How exactly these equations will be solved, can be customized via the
arguments of the `run()`

method. You can choose between different backends, numerical solvers and integration
step-sizes, for instance. Here, we choose the scipy backend, which uses a 4(5)th order Runge-Kutta method as default
numerical solver. The default of the `run()`

method is the forward Euler method that is implemented in PyRates
itself.

```
results = jrc.run(simulation_time=2.0,
step_size=1e-4,
sampling_step_size=1e-3,
outputs={'V_pce': 'pc/rpo_e_in/v',
'V_pci': 'pc/rpo_i/v'},
backend='default',
solver='scipy')
```

### 2.1.1.4. Step 4: Visualization of the solution¶

The output of the `run()`

method is a `pandas.Dataframe`

, which comes with a `plot()`

method for
plotting the timeseries it contains.
This timeseries represents the numerical solution of the initial value problem solved in step 4 with respect to the
state variables and of the model.

```
import matplotlib.pyplot as plt
plt.plot(results)
```

To visualize the average membrane potential at the PC somata, simply plot the sum of the excitatory and inhibitory post-synaptic potentials and :

```
v_pc = results['V_pce'] + results['V_pci']
plt.plot(v_pc)
plt.legend(['V_pce', 'V_pci', 'V_pce - V_pci'])
plt.ylabel('V')
plt.xlabel('time')
plt.show()
```

**Total running time of the script:** ( 0 minutes 0.000 seconds)