NEST implementation of the astrocyte_lr_1994 model

The purpose of this notebook is to provide a reference solution for the the astrocyte_lr_1994 model in NEST. The model is based on Li, Y. X., & Rinzel, J. (1994), De Young, G. W., & Keizer, J. (1992), and Nadkarni, S., & Jung, P. (2003).

This notebook demonstrates how the dynamics of astrocyte_lr_1994 is implemented, and generates a recording of the dynamics (test_astrocyte.dat). This recording serves as a reference for the verification of astrocyte_lr_1994 using test_astrocyte.py in the PyNEST tests.

The problem

The equations governing the evolution of the astrocyte model are

\[\begin{split}\left\lbrace\begin{array}{rcl} \frac{d[\mathrm{Ca^{2+}}](t)}{dt} &=& J_{\mathrm{channel}}(t) - J_{\mathrm{pump}}(t) + J_{\mathrm{leak}}(t) \\ J_{\mathrm{channel}}(t) &=& \mathrm{r_{ER,cyt}} \cdot \mathrm{v_{IP_3R}} \cdot m_{\infty}(t)^3 \cdot n_{\infty}(t)^3\cdot h_{\mathrm{\mathrm{IP_3R}}}(t)^3 \cdot \left([\mathrm{Ca^{2+}}]_{\mathrm{ER}} - [\mathrm{Ca^{2+}}](t)\right) \\ J_{\mathrm{pump},k}(t) &=& \frac{\mathrm{v_{SERCA}} \cdot [\mathrm{Ca^{2+}}](t)^2}{K_{\mathrm{m,SERCA}}^2+[\mathrm{Ca^{2+}}](t)^2} \\ J_{\mathrm{leak}}(t) &=& \mathrm{r_{ER,cyt}} \cdot \mathrm{v_L} \cdot \left( [\mathrm{Ca^{2+}}]_{\mathrm{ER}}(t) - [\mathrm{Ca^{2+}}](t) \right) \\ m_{\infty}(t) &=& \frac{[\mathrm{IP_3}](t)}{[\mathrm{IP_3}](t) + \mathrm{K_{d,IP_3,1}}} \\ n_{\infty}(t) &=& \frac{[\mathrm{Ca^{2+}}](t)}{[\mathrm{Ca^{2+}}](t) + \mathrm{K_{d,act}}} \\ [\mathrm{Ca^{2+}}]_{\mathrm{ER}}(t) &=& \frac{[\mathrm{Ca^{2+}}]_{\mathrm{tot}}-[\mathrm{Ca^{2+}}](t)}{\mathrm{r_{ER,cyt}}} \\ \frac{dh_{\mathrm{IP_3}}(t)}{dt} &=& \alpha_{h_{\mathrm{IP_3}}}(t) \cdot (1-h_{\mathrm{IP_3}}(t)) - \beta_{h_{\mathrm{IP_3}}}(t) \cdot h_{\mathrm{IP_3}} (t) \\ \alpha_{h_{\mathrm{IP_3}}}(t) &=& \mathrm{k_{IP_3R}} \cdot \mathrm{K_{d,inh}} \cdot \frac{[\mathrm{IP_3}](t) + \mathrm{K_{d,IP_3,1}}}{[\mathrm{IP_3}](t)+\mathrm{K_{d,IP_3,2}}} \\ \beta_{h_{\mathrm{IP_3}}}(t) &=& \mathrm{k_{IP_3R}} \cdot [\mathrm{Ca^{2+}}](t) \\ \frac{d[\mathrm{IP_3}](t)}{dt} &=& \frac{[\mathrm{IP_3}]^{*} - [\mathrm{IP_3}](t)}{\tau_\mathrm{IP_3}} + \Delta_\mathrm{IP3} \cdot J_\mathrm{syn}(t) \end{array}\right.\end{split}\]

where \([\mathrm{IP_3}]\), \([\mathrm{Ca^{2+}}]\), and \(h_{\mathrm{IP_3}}\) are the state variables of interest.

Technical details and requirements

Implementation of the functions

  • The reference solution is implemented through scipy.

Requirements

To run this notebook, you need:

[1]:
import numpy as np
from scipy.integrate import odeint

Reference solution

Right hand side function

[2]:
def rhs(y, _, p):
    """
    Implementation of astrocyte dynamics.

    Parameters
    ----------
    y : list
        Vector containing the state variables [IP3, Ca_astro, h_IP3R]
    _ : unused var
    p : Params instance
        Object containing the astrocyte parameters.

    Returns
    -------
    dCa_astro : double
        Derivative of Ca_astro
    dIP3 : double
        Derivative of IP3
    dh_IP3R : double
        Derivative of h_IP3R
    """
    IP3 = y[0]
    Ca_astro = y[1]
    h_IP3R = y[2]

    Ca_astro = max(0, min(Ca_astro, p.Ca_tot * (1 + p.ratio_ER_cyt)))
    alpha = p.k_IP3R * p.Kd_inh * (IP3 + p.Kd_IP3_1) / (IP3 + p.Kd_IP3_2)
    beta = p.k_IP3R * Ca_astro
    Ca_ER = (p.Ca_tot - Ca_astro) / p.ratio_ER_cyt
    m_inf = IP3 / (IP3 + p.Kd_IP3_1)
    n_inf = Ca_astro / (Ca_astro + p.Kd_act)
    J_ch = p.ratio_ER_cyt * p.rate_IP3R * ((m_inf * n_inf * h_IP3R) ** 3) * (Ca_ER - Ca_astro)
    J_pump = p.rate_SERCA * (Ca_astro**2) / (p.Km_SERCA**2 + Ca_astro**2)
    J_leak = p.ratio_ER_cyt * p.rate_L * (Ca_ER - Ca_astro)

    dCa_astro = J_ch - J_pump + J_leak
    dIP3 = (p.IP3_0 - IP3) / p.tau_IP3
    dh_IP3R = alpha * (1 - h_IP3R) - beta * h_IP3R

    return dIP3, dCa_astro, dh_IP3R

Complete model

[3]:
def scipy_astrocyte(p, f, simtime, dt, spk_ts, spk_ws):
    """
    Complete astrocyte model using scipy `odeint` solver.

    Parameters
    ----------
    p : Params instance
        Object containing the astrocyte parameters.
    f : function
        Right-hand side function
    simtime : double
        Duration of the simulation (will run between
        0 and tmax)
    dt : double
        Time increment.
    spk_ts, spk_ws : list
        Times and weights of spike input to the astrocyte

    Returns
    -------
    t : list
        Times at which the astrocyte state was evaluated.
    y : list
        State values associated to the times in `t`
    fos : list
        List of dictionaries containing additional output
        information from `odeint`
    """
    t = np.arange(0, simtime, dt)  # time axis
    n = len(t)
    y = np.zeros((n, 3))  # state variables: IP3, Ca_astro, h_IP3R
    # for the state variables, assign the same initial values as in test_astrocyte.py
    y[0, 0] = 1.0
    y[0, 1] = 1.0
    y[0, 2] = 1.0
    fos = []  # full output dict from odeint()
    delta_ip3 = 5.0  # parameter determining the increase in IP3 induced by synaptic input

    # update time-step by time-step
    for k in range(1, n):
        # solve ODE from t_k-1 to t_k
        d, fo = odeint(f, y[k - 1, :], t[k - 1 : k + 1], (p,), full_output=True)
        y[k, :] = d[1, :]

        # apply synaptic inputs (spikes)
        if t[k] in spk_ts:
            y[k, 0] += delta_ip3 * spk_ws[spk_ts.index(t[k])]

        fos.append(fo)

    return t, y, fos

Parameters

[4]:
astro_param = {
    "Ca_tot": 2.0,
    "IP3_0": 0.16,
    "Kd_act": 0.08234,
    "Kd_inh": 1.049,
    "Kd_IP3_1": 0.13,
    "Kd_IP3_2": 0.9434,
    "Km_SERCA": 0.1,
    "ratio_ER_cyt": 0.185,
    "delta_IP3": 5.0,
    "k_IP3R": 0.0002,
    "rate_L": 0.00011,
    "tau_IP3": 7142.0,
    "rate_IP3R": 0.006,
    "rate_SERCA": 0.0009,
}


class Params:
    """
    Class giving access to the astrocyte parameters.
    """

    def __init__(self):
        self.params = astro_param
        self.Ca_tot = astro_param["Ca_tot"]
        self.IP3_0 = astro_param["IP3_0"]
        self.Kd_act = astro_param["Kd_act"]
        self.Kd_inh = astro_param["Kd_inh"]
        self.Kd_IP3_1 = astro_param["Kd_IP3_1"]
        self.Kd_IP3_2 = astro_param["Kd_IP3_2"]
        self.Km_SERCA = astro_param["Km_SERCA"]
        self.ratio_ER_cyt = astro_param["ratio_ER_cyt"]
        self.delta_IP3 = astro_param["delta_IP3"]
        self.k_IP3R = astro_param["k_IP3R"]
        self.rate_L = astro_param["rate_L"]
        self.tau_IP3 = astro_param["tau_IP3"]
        self.rate_IP3R = astro_param["rate_IP3R"]
        self.rate_SERCA = astro_param["rate_SERCA"]


p = Params()

Simulate the implementation

[5]:
# Parameters for the simulation
simtime = 100.0
resolution = 0.1
spike_times = [10.0]
spike_weights = [1.0]

# Simulate and get the recording
t, y, fos = scipy_astrocyte(p, rhs, simtime, resolution, spike_times, spike_weights)
data = np.concatenate((np.array([t]).T, y), axis=1)

# Save the recording (excluding the initial values)
np.savetxt(
    "test_astrocyte.dat",
    data[1:, :],
    header="""\ntest_astrocyte.dat\n
This file is part of NEST.\n
This .dat file contains the recordings of the state variables of an
astrocyte, simulated using the ODEINT solver of SciPy, with the
implementation detailed in
``doc/htmldoc/model_details/astrocyte_model_implementation.ipynb``.
This data is used as reference for the tests in ``test_astrocyte.py``.\n
       Times                     IP3                    Ca_astro                  h_IP3R""",
)

License

This file is part of NEST. Copyright (C) 2004 The NEST Initiative

NEST is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 2 of the License, or (at your option) any later version.

NEST is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.