Merge pull request #128 from neurolib-dev/feature/concat_noise

Go wild with stimulation and/or noise

Author: caglorithm

.github/workflows/ci.yml

@@ -28,10 +28,8 @@ jobs:
       - name: Install dependencies
         run: |
           python -m pip install --upgrade pip
-          pip install treon flake8 codecov pytest-cov wheel setuptools jupyterlab matplotlib seaborn mne
+          pip install flake8 codecov pytest-cov wheel setuptools matplotlib seaborn mne
           if [ -f requirements.txt ]; then pip install -r requirements.txt; fi
-          #python setup.py sdist bdist_wheel
-          #pip install --upgrade dist/*
           pip install .
       - name: Lint with flake8
         run: |

.github/workflows/test-notebooks.yml

@@ -0,0 +1,33 @@
+# This workflow will install Python dependencies, run tests and lint with a variety of Python versions
+# For more information see: https://help.github.com/actions/language-and-framework-guides/using-python-with-github-actions
+
+name: notebooks
+
+on:
+  push:
+  pull_request:
+
+jobs:
+  build:
+    runs-on: ${{ matrix.os }}
+    strategy:
+      fail-fast: false
+      matrix:
+        os: [ubuntu-latest, macos-latest]
+        python-version: [3.7]
+
+    steps:
+      - uses: actions/checkout@v2
+      - name: Set up Python ${{ matrix.python-version }}
+        uses: actions/setup-python@v2
+        with:
+          python-version: ${{ matrix.python-version }}
+      - name: Install dependencies
+        run: |
+          python -m pip install --upgrade pip
+          pip install treon wheel setuptools jupyterlab matplotlib seaborn
+          if [ -f requirements.txt ]; then pip install -r requirements.txt; fi
+          pip install .
+      - name: Test notebooks with treon
+        run: |
+          treon examples/

examples/example-0.5-aln-external-stimulus.ipynb

@@ -11,7 +11,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 3,
+   "execution_count": 1,
    "metadata": {},
    "outputs": [],
    "source": [
@@ -27,18 +27,9 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 4,
+   "execution_count": 3,
    "metadata": {},
-   "outputs": [
-    {
-     "name": "stderr",
-     "output_type": "stream",
-     "text": [
-      "/Users/caglar/anaconda/lib/python3.7/site-packages/xarray/core/formatting_html.py:6: UserWarning: Module neurolib was already imported from /Users/caglar/anaconda/lib/python3.7/site-packages/neurolib/__init__.py, but /Users/caglar/Documents/PhD/projects/neurolib is being added to sys.path\n",
-      "  import pkg_resources\n"
-     ]
-    }
-   ],
+   "outputs": [],
    "source": [
     "try:\n",
     "    import matplotlib.pyplot as plt\n",
@@ -54,7 +45,8 @@
     "from neurolib.models.aln import ALNModel\n",
     "\n",
     "# Some useful functions are provided here\n",
-    "import neurolib.utils.functions as func"
+    "import neurolib.utils.functions as func\n",
+    "from neurolib.utils.stimulus import construct_stimulus"
    ]
   },
   {
@@ -79,7 +71,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "We can construct a simple stimulus using the function `func.construct_stimulus`:"
+    "We can construct a simple stimulus using the function `stimulus.construct_stimulus`:"
    ]
   },
   {
@@ -88,7 +80,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "stimulus = func.construct_stimulus(\"rect\", duration=model.params.duration, dt=model.params.dt, stim_amp=1.0, stim_freq=1)"
+    "stimulus = construct_stimulus(\"rect\", duration=model.params.duration, dt=model.params.dt, stim_amp=1.0, stim_freq=1)"
    ]
   },
   {
@@ -97,7 +89,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "stimulus = func.construct_stimulus(\"ac\", duration=model.params.duration, dt=model.params.dt, stim_amp=1.0, stim_freq=1)"
+    "stimulus = construct_stimulus(\"ac\", duration=model.params.duration, dt=model.params.dt, stim_amp=1.0, stim_freq=1)"
    ]
   },
   {
@@ -272,7 +264,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "`neurolib` helps you to create a few basic stimuli out of the box using the function `functions.construct_stimulus()`."
+    "`neurolib` helps you to create a few basic stimuli out of the box using the function `stimulus.construct_stimulus()`."
    ]
   },
   {
@@ -282,7 +274,7 @@
    "outputs": [],
    "source": [
     "# construct a stimulus\n",
-    "ac_stimulus = func.construct_stimulus(stim=\"ac\", stim_freq = 25, duration=model.params.duration, dt=model.params.dt)\n",
+    "ac_stimulus = construct_stimulus(stim=\"ac\", stim_freq = 25, duration=model.params.duration, dt=model.params.dt)\n",
     "\n",
     "# this stimulus is 1-dimensional. neurolib will threfore automatically apply it to *all nodes*.\n",
     "model.params['ext_exc_current'] = ac_stimulus * 5.0 "
@@ -335,7 +327,7 @@
    "outputs": [],
    "source": [
     "# this stimulus is 1-dimensional\n",
-    "ac_stimulus = func.construct_stimulus(stim=\"ac\", stim_freq = 25, duration=model.params.duration, dt=model.params.dt)\n",
+    "ac_stimulus = construct_stimulus(stim=\"ac\", stim_freq = 25, duration=model.params.duration, dt=model.params.dt)\n",
     "\n",
     "# let's make a N-dimensional stimulus vector out of it, by copying and pasting each entry N times\n",
     "ac_stimulus = np.tile(ac_stimulus.T, (model.params.N, 1))\n",

examples/example-1.3-aln-bifurcation-diagram.ipynb

@@ -57,6 +57,7 @@
     "from neurolib.utils.parameterSpace import ParameterSpace\n",
     "from neurolib.optimize.exploration import BoxSearch\n",
     "import neurolib.utils.functions as func\n",
+    "import neurolib.utils.stimulus as stim\n",
     "import neurolib.optimize.exploration.explorationUtils as eu\n",
     "import neurolib.utils.devutils as du\n",
     "from neurolib.utils.loadData import Dataset"
@@ -118,7 +119,7 @@
     "model.params['mui_ext_mean'] = 2.5\n",
     "\n",
     "# construct a stimulus\n",
-    "rect_stimulus = func.construct_stimulus(stim=\"rect\", duration=model.params.duration, dt=model.params.dt)\n",
+    "rect_stimulus = stim.construct_stimulus(stim=\"rect\", duration=model.params.duration, dt=model.params.dt)\n",
     "model.params['ext_exc_current'] = rect_stimulus * 5.0    \n",
     "    \n",
     "model.run()"
@@ -195,7 +196,7 @@
     "    # ---------------------------------------    \n",
     "    model.params['duration'] = defaultDuration\n",
     "    \n",
-    "    rect_stimulus = func.construct_stimulus(stim=\"rect\", duration=model.params.duration, dt=model.params.dt)\n",
+    "    rect_stimulus = stim.construct_stimulus(stim=\"rect\", duration=model.params.duration, dt=model.params.dt)\n",
     "    model.params['ext_exc_current'] = rect_stimulus * 5.0    \n",
     "    \n",
     "    model.run()\n",

neurolib/models/multimodel/builder/aln.py

@@ -8,10 +8,10 @@
 from h5py import File
 from jitcdde import input as system_input
 
+from ....utils.stimulus import OrnsteinUhlenbeckProcess
 from ..builder.base.constants import EXC, INH, LAMBDA_SPEED
 from ..builder.base.network import Network, SingleCouplingExcitatoryInhibitoryNode
 from ..builder.base.neural_mass import NeuralMass
-from .model_input import OrnsteinUhlenbeckProcess
 
 DEFAULT_QUANTITIES_CASCADE_FILENAME = "quantities_cascade.h5"
 
@@ -422,7 +422,7 @@ class ExcitatoryALNMass(ALNMass):
         "lambda",
     ]
 
-    noise_input = [OrnsteinUhlenbeckProcess(mu=0.4, sigma=0.0, tau=5.0)]
+    _noise_input = [OrnsteinUhlenbeckProcess(mu=0.4, sigma=0.0, tau=5.0)]
 
     @staticmethod
     def _rescale_strengths(params):
@@ -619,7 +619,7 @@ class InhibitoryALNMass(ALNMass):
         "lambda",
     ]
 
-    noise_input = [OrnsteinUhlenbeckProcess(mu=0.3, sigma=0.0, tau=5.0)]
+    _noise_input = [OrnsteinUhlenbeckProcess(mu=0.3, sigma=0.0, tau=5.0)]
 
     @staticmethod
     def _rescale_strengths(params):

neurolib/models/multimodel/builder/base/network.py

@@ -1,5 +1,5 @@
 import logging
-from itertools import chain
+from itertools import chain, islice
 
 import numpy as np
 import symengine as se
@@ -47,7 +47,6 @@ def __init__(self, neural_masses):
         # for all masses
         self.num_state_variables = sum([mass.num_state_variables for mass in self])
         self.num_noise_variables = sum([mass.num_noise_variables for mass in self])
-        self.noise_input = sum([mass.noise_input for mass in self], [])
         assert len(self.noise_input) == self.num_noise_variables
         self.idx_state_var = None
         self.initialised = False
@@ -109,6 +108,19 @@ def state_variable_names(self):
     def max_delay(self):
         return 0.0
 
+    @property
+    def noise_input(self):
+        return sum([mass.noise_input for mass in self], [])
+
+    @noise_input.setter
+    def noise_input(self, new_noise):
+        assert len(new_noise) == self.num_noise_variables
+        masses_noise_length = [mass.num_noise_variables for mass in self.masses]
+        new_noise = iter(new_noise)
+        for noise_chunk, mass in zip([list(islice(new_noise, 0, i)) for i in masses_noise_length], self.masses):
+            assert len(noise_chunk) == mass.num_noise_variables
+            mass.noise_input = noise_chunk
+
     def get_nested_params(self):
         """
         Return nested dictionary with parameters from all masses within this
@@ -466,7 +478,6 @@ def __init__(self, nodes, connectivity_matrix, delay_matrix=None):
         self.nodes = nodes
         self.num_state_variables = sum([node.num_state_variables for node in self])
         self.num_noise_variables = sum([node.num_noise_variables for node in self])
-        self.noise_input = sum([mass.noise_input for mass in self], [])
         assert len(self.noise_input) == self.num_noise_variables
         assert connectivity_matrix.shape[0] == self.num_nodes
         if delay_matrix is None:
@@ -556,6 +567,19 @@ def initial_state(self, initial_state):
         assert initial_state.shape[0] == self.num_state_variables
         self._initial_state = initial_state
 
+    @property
+    def noise_input(self):
+        return sum([node.noise_input for node in self], [])
+
+    @noise_input.setter
+    def noise_input(self, new_noise):
+        assert len(new_noise) == self.num_noise_variables
+        nodes_noise_length = [node.num_noise_variables for node in self.nodes]
+        new_noise = iter(new_noise)
+        for noise_chunk, node in zip([list(islice(new_noise, 0, i)) for i in nodes_noise_length], self.nodes):
+            assert len(noise_chunk) == node.num_noise_variables
+            node.noise_input = noise_chunk
+
     @staticmethod
     def _strip_index(symbol_name):
         """

neurolib/models/multimodel/builder/base/neural_mass.py

@@ -4,7 +4,7 @@
 import symengine as se
 from jitcdde import y as state_vector
 
-from ..model_input import ModelInput
+from .....utils.stimulus import ModelInput
 
 
 class NeuralMass:
@@ -43,7 +43,7 @@ class NeuralMass:
     num_noise_variables = 0
 
     # inputs to the mass - usually noise
-    noise_input = []
+    _noise_input = []
 
     # names for the state variables to link them with results
     state_variable_names = []
@@ -147,6 +147,17 @@ def _validate_callbacks(self, callback_list):
         )
         assert all(callable(callback[1]) for callback in callback_list)
 
+    @property
+    def noise_input(self):
+        return self._noise_input
+
+    @noise_input.setter
+    def noise_input(self, new_noise):
+        assert len(new_noise) == self.num_noise_variables
+        self._noise_input = new_noise
+        # update noise parameters
+        self._get_params_from_noise()
+
     def init_mass(self, start_idx_for_noise=None):
         """
         Initialise neural mass. Usually just initialise the state vector,

neurolib/models/multimodel/builder/fitzhugh_nagumo.py

@@ -1,9 +1,9 @@
 import numpy as np
 from jitcdde import input as system_input
 
+from ....utils.stimulus import OrnsteinUhlenbeckProcess
 from ..builder.base.network import Network, Node
 from ..builder.base.neural_mass import NeuralMass
-from .model_input import OrnsteinUhlenbeckProcess
 
 FHN_DEFAULT_PARAMS = {
     "alpha": 3.0,
@@ -40,7 +40,7 @@ class FitzHughNagumoMass(NeuralMass):
         "y_ext",
     ]
     required_couplings = ["network_x", "network_y"]
-    noise_input = [
+    _noise_input = [
         OrnsteinUhlenbeckProcess(mu=0.0, sigma=0.0, tau=5.0),
         OrnsteinUhlenbeckProcess(mu=0.0, sigma=0.0, tau=5.0),
     ]

neurolib/models/multimodel/builder/hopf.py

@@ -1,9 +1,9 @@
 import numpy as np
 from jitcdde import input as system_input
 
+from ....utils.stimulus import OrnsteinUhlenbeckProcess
 from ..builder.base.network import Network, Node
 from ..builder.base.neural_mass import NeuralMass
-from .model_input import OrnsteinUhlenbeckProcess
 
 HOPF_DEFAULT_PARAMS = {
     "a": 0.25,
@@ -35,7 +35,7 @@ class HopfMass(NeuralMass):
     state_variable_names = ["x", "y"]
     required_params = ["a", "w", "x_ext", "y_ext"]
     required_couplings = ["network_x", "network_y"]
-    noise_input = [
+    _noise_input = [
         OrnsteinUhlenbeckProcess(mu=0.0, sigma=0.0, tau=5.0),
         OrnsteinUhlenbeckProcess(mu=0.0, sigma=0.0, tau=5.0),
     ]

neurolib/models/multimodel/builder/thalamus.py

@@ -2,10 +2,10 @@
 from jitcdde import input as system_input
 from symengine import exp
 
+from ....utils.stimulus import OrnsteinUhlenbeckProcess, ZeroInput
 from ..builder.base.constants import EXC, INH, LAMBDA_SPEED
 from ..builder.base.network import SingleCouplingExcitatoryInhibitoryNode
 from ..builder.base.neural_mass import NeuralMass
-from .model_input import OrnsteinUhlenbeckProcess, ZeroInput
 
 TCR_DEFAULT_PARAMS = {
     "tau": 20.0,  # ms
@@ -167,7 +167,7 @@ class ThalamocorticalMass(ThalamicMass):
         "ext_current",
         "lambda",
     ]
-    noise_input = [OrnsteinUhlenbeckProcess(mu=0.0, sigma=0.0, tau=5.0)]
+    _noise_input = [OrnsteinUhlenbeckProcess(mu=0.0, sigma=0.0, tau=5.0)]
 
     def __init__(self, params=None):
         super().__init__(params=params or TCR_DEFAULT_PARAMS)
@@ -332,7 +332,7 @@ class ThalamicReticularMass(ThalamicMass):
         "ext_current",
         "lambda",
     ]
-    noise_input = [ZeroInput()]
+    _noise_input = [ZeroInput()]
 
     def __init__(self, params=None):
         super().__init__(params=params or TRN_DEFAULT_PARAMS)