Skip to content
Snippets Groups Projects
Commit 889aabf4 authored by Didi Hou's avatar Didi Hou Committed by Administrator
Browse files

/

parent 98d0a90b
No related branches found
No related tags found
1 merge request!35Pre-release MAM v1.1.0
...@@ -85,9 +85,9 @@ def plot_resting_state(A): ...@@ -85,9 +85,9 @@ def plot_resting_state(A):
gs1 = gridspec.GridSpec(1, 3) gs1 = gridspec.GridSpec(1, 3)
gs1.update(left=0.06, right=0.72, top=0.95, wspace=0.4, bottom=0.35) gs1.update(left=0.06, right=0.72, top=0.95, wspace=0.4, bottom=0.35)
axes['A'] = pl.subplot(gs1[:-1, :1]) axes['A'] = pl.subplot(gs1[:1, :1])
axes['B'] = pl.subplot(gs1[:-1, 1:2]) axes['B'] = pl.subplot(gs1[:1, 1:2])
axes['C'] = pl.subplot(gs1[:-1, 2:]) axes['C'] = pl.subplot(gs1[:1, 2:])
gs2 = gridspec.GridSpec(3, 1) gs2 = gridspec.GridSpec(3, 1)
gs2.update(left=0.78, right=0.95, top=0.95, bottom=0.35) gs2.update(left=0.78, right=0.95, top=0.95, bottom=0.35)
......
...@@ -85,9 +85,9 @@ def plot_resting_state(A): ...@@ -85,9 +85,9 @@ def plot_resting_state(A):
gs1 = gridspec.GridSpec(1, 3) gs1 = gridspec.GridSpec(1, 3)
gs1.update(left=0.06, right=0.72, top=0.95, wspace=0.4, bottom=0.35) gs1.update(left=0.06, right=0.72, top=0.95, wspace=0.4, bottom=0.35)
axes['A'] = pl.subplot(gs1[:-1, :1]) axes['A'] = pl.subplot(gs1[:1, :1])
axes['B'] = pl.subplot(gs1[:-1, 1:2]) axes['B'] = pl.subplot(gs1[:1, 1:2])
axes['C'] = pl.subplot(gs1[:-1, 2:]) axes['C'] = pl.subplot(gs1[:1, 2:])
gs2 = gridspec.GridSpec(3, 1) gs2 = gridspec.GridSpec(3, 1)
gs2.update(left=0.78, right=0.95, top=0.95, bottom=0.35) gs2.update(left=0.78, right=0.95, top=0.95, bottom=0.35)
......
%% Cell type:markdown id:b1331599 tags: %% Cell type:markdown id:b1331599 tags:
# Down-scaled multi-area model # Down-scaled multi-area model
%% Cell type:markdown id:b952d0ea tags: %% Cell type:markdown id:b952d0ea tags:
#### Notebook structure <a class="anchor" id="toc"></a> #### Notebook structure <a class="anchor" id="toc"></a>
* [S0. Configuration](#section_0) * [S0. Configuration](#section_0)
* [S1. Paramters specification](#section_1) * [S1. Paramters specification](#section_1)
* [1.1. Parameters to tune](#section_1_1) * [1.1. Parameters to tune](#section_1_1)
* [1.2. Default parameters](#section_1_2) * [1.2. Default parameters](#section_1_2)
* [S2. Multi-area model instantiation and simulation](#section_2) * [S2. Multi-area model instantiation and simulation](#section_2)
* [2.1. Insantiate a multi-area model](#section_2_1) * [2.1. Insantiate a multi-area model](#section_2_1)
* [2.2. Predict firing rates from theory](#section_2_2) * [2.2. Predict firing rates from theory](#section_2_2)
* [2.3. Extract interarea connectivity](#section_2_3) * [2.3. Extract interarea connectivity](#section_2_3)
* [2.4. Run the simulation](#section_2_4) * [2.4. Run the simulation](#section_2_4)
* [S3. Simulation results validation and connection extraction](#section_3) * [S3. Simulation results validation and connection extraction](#section_3)
* [S4. Data loading and processing](#section_4) * [S4. Data loading and processing](#section_4)
* [S5. Simulation results visualization](#section_5) * [S5. Simulation results visualization](#section_5)
* [5.1. Instantaneous and mean firing rate across all populations](#section_5_1) * [5.1. Instantaneous and mean firing rate across all populations](#section_5_1)
* [5.2. Raster plot of spiking activity for single area](#section_5_2) * [5.2. Raster plot of spiking activity for single area](#section_5_2)
* [5.3. Population-averaged firing rate](#section_5_3) * [5.3. Population-averaged firing rate](#section_5_3)
* [5.4 Time-averaged population rates](#section_5_4) * [5.4 Time-averaged population rates](#section_5_4)
* [5.5. Average pairwise correlation coefficients of spiking activity](#section_5_5) * [5.5. Average pairwise correlation coefficients of spiking activity](#section_5_5)
* [5.6. Irregularity of spiking activity](#section_5_6) * [5.6. Irregularity of spiking activity](#section_5_6)
* [5.7. Time series of population- and area-averaged firing rates](#section_5_7) * [5.7. Time series of population- and area-averaged firing rates](#section_5_7)
%% Cell type:markdown id:bd3d4b0e tags: %% Cell type:markdown id:bd3d4b0e tags:
<br> <br>
%% Cell type:markdown id:d782e527 tags: %% Cell type:markdown id:d782e527 tags:
## S0. Configuration <a class="anchor" id="section_0"></a> ## S0. Configuration <a class="anchor" id="section_0"></a>
%% Cell type:code id:06e7f341-c226-4782-a913-4868027d7d06 tags: %% Cell type:code id:06e7f341-c226-4782-a913-4868027d7d06 tags:
``` python ``` python
# Create config file # Create config file
with open('config.py', 'w') as fp: with open('config.py', 'w') as fp:
fp.write( fp.write(
'''import os '''import os
base_path = os.path.abspath(".") base_path = os.path.abspath(".")
data_path = os.path.abspath("simulations") data_path = os.path.abspath("simulations")
jobscript_template = "python {base_path}/run_simulation.py {label}" jobscript_template = "python {base_path}/run_simulation.py {label}"
submit_cmd = "bash -c" submit_cmd = "bash -c"
''') ''')
``` ```
%% Cell type:code id:96517739 tags: %% Cell type:code id:96517739 tags:
``` python ``` python
# Import dependencies # Import dependencies
%matplotlib inline %matplotlib inline
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
import numpy as np import numpy as np
import os import os
import nest import nest
from IPython.display import display, HTML from IPython.display import display, HTML
import json import json
# Import the MultiAreaModel class # Import the MultiAreaModel class
from multiarea_model import MultiAreaModel from multiarea_model import MultiAreaModel
from multiarea_model import Analysis from multiarea_model import Analysis
from config import base_path, data_path from config import base_path, data_path
import sys import sys
sys.path.append('./figures') sys.path.append('./figures')
from MAM2EBRAINS import plot_instan_mean_firing_rate from MAM2EBRAINS import plot_instan_mean_firing_rate
``` ```
%% Output %% Output
-- N E S T -- -- N E S T --
Copyright (C) 2004 The NEST Initiative Copyright (C) 2004 The NEST Initiative
Version: 3.4 Version: 3.4
Built: May 17 2023 20:48:31 Built: May 17 2023 20:48:31
This program is provided AS IS and comes with This program is provided AS IS and comes with
NO WARRANTY. See the file LICENSE for details. NO WARRANTY. See the file LICENSE for details.
Problems or suggestions? Problems or suggestions?
Visit https://www.nest-simulator.org Visit https://www.nest-simulator.org
Type 'nest.help()' to find out more about NEST. Type 'nest.help()' to find out more about NEST.
%% Cell type:code id:7e07b0d0 tags: %% Cell type:code id:7e07b0d0 tags:
``` python ``` python
!pip install nested_dict dicthash !pip install nested_dict dicthash
``` ```
%% Output %% Output
Requirement already satisfied: nested_dict in /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/site-packages (1.61) Requirement already satisfied: nested_dict in /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/site-packages (1.61)
Requirement already satisfied: dicthash in /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/site-packages (0.0.2) Requirement already satisfied: dicthash in /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/site-packages (0.0.2)
Requirement already satisfied: future in /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/site-packages (from dicthash) (0.18.2) Requirement already satisfied: future in /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/site-packages (from dicthash) (0.18.2)
%% Cell type:code id:1d440c07-9b69-4e52-8573-26b13493bc5a tags: %% Cell type:code id:1d440c07-9b69-4e52-8573-26b13493bc5a tags:
``` python ``` python
# Jupyter notebook display format setting # Jupyter notebook display format setting
style = """ style = """
<style> <style>
table {float:left} table {float:left}
</style> </style>
""" """
display(HTML(style)) display(HTML(style))
``` ```
%% Output %% Output
%% Cell type:markdown id:27160ba8 tags: %% Cell type:markdown id:27160ba8 tags:
Go back to [Notebook structure](#toc) Go back to [Notebook structure](#toc)
%% Cell type:markdown id:565be233 tags: %% Cell type:markdown id:565be233 tags:
<br> <br>
%% Cell type:markdown id:df83f5ea-1c4b-44d3-9926-01786aa46e14 tags: %% Cell type:markdown id:df83f5ea-1c4b-44d3-9926-01786aa46e14 tags:
## S1. Paramters specification <a class="anchor" id="section_1"></a> ## S1. Paramters specification <a class="anchor" id="section_1"></a>
%% Cell type:markdown id:30655817 tags: %% Cell type:markdown id:30655817 tags:
### 1.1. Parameters to tune <a class="anchor" id="section_1_1"></a> ### 1.1. Parameters to tune <a class="anchor" id="section_1_1"></a>
%% Cell type:markdown id:4f67c1ba tags: %% Cell type:markdown id:4f67c1ba tags:
|Parameter |Default value |Value range/options |Value assigned |Description | |Parameter |Default value |Value range/options |Value assigned |Description |
|:----------------------------:|:-----------------------:|:--------------------------------------------------------------------:|:------------------:|:-----------:| |:----------------------------:|:-----------------------:|:--------------------------------------------------------------------:|:------------------:|:-----------:|
|scale_down_to |1. |(0, 1.] |0.005 |$^1$ | |scale_down_to |1. |(0, 1.] |0.005 |$^1$ |
|cc_weights_factor |1. |(0, 1.] |1. |$^2$ | |cc_weights_factor |1. |(0, 1.] |1. |$^2$ |
|areas_simulated |complete_area_list |All sublists of complete_area_list |complete_area_list |$^3$ | |areas_simulated |complete_area_list |All sublists of complete_area_list |complete_area_list |$^3$ |
|replace_non_simulated_areas |None |None, 'hom_poisson_stat', 'het_poisson_stat', 'het_current_nonstat' |'het_poisson_stat' |$^4$ | |replace_non_simulated_areas |None |None, 'hom_poisson_stat', 'het_poisson_stat', 'het_current_nonstat' |'het_poisson_stat' |$^4$ |
%% Cell type:markdown id:a2161477 tags: %% Cell type:markdown id:a2161477 tags:
1. `scale_down_to` <br> 1. `scale_down_to` <br>
`scale_down_to` is the down-scaling factor which defines the the ratio of the full scale multi-area model being down-scaled to a model with fewer neurons and indegrees so as to be simulated on machines with lower computational ability and the simulation results can be obtained within relative shorter period of time. <br> Its deafualt value if `1.` meaning full scale simulation. <br> In the pre-set downscale version, it's set as `0.005`, where the numer of neurons and indegrees are both scaled down to 0.5% of its full scale amount, where the model can usually be simulated on a local machine. <br> **Warning**: This will not yield reasonable dynamical results from the network and is only meant to demonstrate the simulation workflow <br> `scale_down_to` is the down-scaling factor which defines the the ratio of the full scale multi-area model being down-scaled to a model with fewer neurons and indegrees so as to be simulated on machines with lower computational ability and the simulation results can be obtained within relative shorter period of time. <br> Its deafualt value if `1.` meaning full scale simulation. <br> In the pre-set downscale version, it's set as `0.005`, where the numer of neurons and indegrees are both scaled down to 0.5% of its full scale amount, where the model can usually be simulated on a local machine. <br> **Warning**: This will not yield reasonable dynamical results from the network and is only meant to demonstrate the simulation workflow <br>
2. `cc_weights_factor` <br> 2. `cc_weights_factor` <br>
This scaling factor controls the cortico-cortical synaptic strength. <br> By default it's set as `1.0`, where the inter-area synaptic strength is the same as the intra-areal. <br> **Important**: This factor changes the network activity from ground state to metastable state. <br> This scaling factor controls the cortico-cortical synaptic strength. <br> By default it's set as `1.0`, where the inter-area synaptic strength is the same as the intra-areal. <br> **Important**: This factor changes the network activity from ground state to metastable state. <br>
3. `areas_simulated` <br> 3. `areas_simulated` <br>
This parameter specifies the cortical areas included in the simulation process. Its default value is `complete_area_list` meaning all the areas in the complete_area_list will be actually simulated. <br> This parameter specifies the cortical areas included in the simulation process. Its default value is `complete_area_list` meaning all the areas in the complete_area_list will be actually simulated. <br>
complete_area_list = ['V1', 'V2', 'VP', 'V3', 'V3A', 'MT', 'V4t', 'V4', 'VOT', 'MSTd', 'PIP', 'PO', 'DP', 'MIP', 'MDP', 'VIP', 'LIP', 'PITv', 'PITd', 'MSTl', 'CITv', 'CITd', 'FEF', 'TF', 'AITv', 'FST', '7a', 'STPp', 'STPa', '46', 'AITd', 'TH'] <br> complete_area_list = ['V1', 'V2', 'VP', 'V3', 'V3A', 'MT', 'V4t', 'V4', 'VOT', 'MSTd', 'PIP', 'PO', 'DP', 'MIP', 'MDP', 'VIP', 'LIP', 'PITv', 'PITd', 'MSTl', 'CITv', 'CITd', 'FEF', 'TF', 'AITv', 'FST', '7a', 'STPp', 'STPa', '46', 'AITd', 'TH'] <br>
The value assigned to simulation_areas can be any sublist of the compete_area_list specifying areas a user want to include in his/her simulation. <br> The value assigned to simulation_areas can be any sublist of the compete_area_list specifying areas a user want to include in his/her simulation. <br>
4. `replace_non_simulated_areas` <br> 4. `replace_non_simulated_areas` <br>
The paramter `replace_non_simulated_areas` defines how non-simulated areas will be replaced. <br> It's set as `None` by default when the parameter areas_simulated is set as full_area_list where all areas will be simulated so that no areas need to be replaced. <br> Other options are: `'hom_poisson_stat'`, `'het_poisson_stat'`, and `'het_current_nonstat'`. `'hom_poisson_stat'` is a manually set parameter which can be tuned. When it's set as 'het_poisson_stat' or 'het_current_nonstat', the data to replace the cortico-cortical input is loaded from 'replace_cc_input_source' which is the firing rates of our full scale simulation results. The differenc between 'het_poisson_stat' and 'het_current_nonstat' is that 'het_poisson_stat' is the mean of the time-series firing rate so that it's static, yet 'het_current_nonstat' is time-varying specific current, which is varying by time. The paramter `replace_non_simulated_areas` defines how non-simulated areas will be replaced. <br> It's set as `None` by default when the parameter areas_simulated is set as full_area_list where all areas will be simulated so that no areas need to be replaced. <br> Other options are: `'hom_poisson_stat'`, `'het_poisson_stat'`, and `'het_current_nonstat'`. `'hom_poisson_stat'` is a manually set parameter which can be tuned. When it's set as 'het_poisson_stat' or 'het_current_nonstat', the data to replace the cortico-cortical input is loaded from 'replace_cc_input_source' which is the firing rates of our full scale simulation results. The differenc between 'het_poisson_stat' and 'het_current_nonstat' is that 'het_poisson_stat' is the mean of the time-series firing rate so that it's static, yet 'het_current_nonstat' is time-varying specific current, which is varying by time.
%% Cell type:code id:60265d52 tags: %% Cell type:code id:60265d52 tags:
``` python ``` python
# Downscaling factor # Downscaling factor
scale_down_to = 0.005 # Change it to 1. for running the fullscale network scale_down_to = 0.005 # Change it to 1. for running the fullscale network
# Scaling factor for cortico-cortical connections (chi) # Scaling factor for cortico-cortical connections (chi)
cc_weights_factor = 1. cc_weights_factor = 1.
# Cortical areas included in the simulation # Cortical areas included in the simulation
areas_simulated = ['V1', 'V2', 'VP', 'V3', 'V3A', 'MT', 'V4t', 'V4', 'VOT', 'MSTd', 'PIP', 'PO', 'DP', 'MIP', 'MDP', 'VIP', 'LIP', 'PITv', 'PITd', 'MSTl', 'CITv', 'CITd', 'FEF', 'TF', 'AITv', 'FST', '7a', 'STPp', 'STPa', '46', 'AITd', 'TH'] areas_simulated = ['V1', 'V2', 'VP', 'V3', 'V3A', 'MT', 'V4t', 'V4', 'VOT', 'MSTd', 'PIP', 'PO', 'DP', 'MIP', 'MDP', 'VIP', 'LIP', 'PITv', 'PITd', 'MSTl', 'CITv', 'CITd', 'FEF', 'TF', 'AITv', 'FST', '7a', 'STPp', 'STPa', '46', 'AITd', 'TH']
# Firing rates used to replace the non-simulated areas # Firing rates used to replace the non-simulated areas
replace_non_simulated_areas = 'het_poisson_stat' replace_non_simulated_areas = 'het_poisson_stat'
``` ```
%% Cell type:markdown id:de11b07f tags: %% Cell type:markdown id:de11b07f tags:
### 1.2. Default parameters <a class="anchor" id="section_1_2"></a> ### 1.2. Default parameters <a class="anchor" id="section_1_2"></a>
We try our best not to confuse users with too many parameters. However, if you want to change more parameters and explore the model, you can do so by passing a dictionary to the `default_params` argument of the `MultiAreaModel` class. We try our best not to confuse users with too many parameters. However, if you want to change more parameters and explore the model, you can do so by passing a dictionary to the `default_params` argument of the `MultiAreaModel` class.
%% Cell type:code id:6e4bed8d tags: %% Cell type:code id:6e4bed8d tags:
``` python ``` python
# Connection parameters # Connection parameters
conn_params = { conn_params = {
'replace_non_simulated_areas': 'het_poisson_stat', # Whether to replace non-simulated areas by Poisson sources with the same global rate, by default: None 'replace_non_simulated_areas': 'het_poisson_stat', # Whether to replace non-simulated areas by Poisson sources with the same global rate, by default: None
'g': -11., # It sets the relative inhibitory synaptic strength, by default: -16. 'g': -11., # It sets the relative inhibitory synaptic strength, by default: -16.
'K_stable': 'K_stable.npy', # Whether to apply the stabilization method of Schuecker, Schmidt et al. (2017), by default: None 'K_stable': 'K_stable.npy', # Whether to apply the stabilization method of Schuecker, Schmidt et al. (2017), by default: None
'fac_nu_ext_TH': 1.2, # Increase the external input to 2/3E and 5E in area TH 'fac_nu_ext_TH': 1.2, # Increase the external input to 2/3E and 5E in area TH
'fac_nu_ext_5E': 1.125, # Increase the external Poisson indegree onto 5E 'fac_nu_ext_5E': 1.125, # Increase the external Poisson indegree onto 5E
'fac_nu_ext_6E': 1.41666667, # Increase the external Poisson indegree onto 6E 'fac_nu_ext_6E': 1.41666667, # Increase the external Poisson indegree onto 6E
'av_indegree_V1': 3950. # Adjust the average indegree in V1 based on monkey data 'av_indegree_V1': 3950. # Adjust the average indegree in V1 based on monkey data
} }
# Input parameters # Input parameters
input_params = { input_params = {
'rate_ext': 10. # Rate of the Poissonian spike generator (in spikes/s) 'rate_ext': 10. # Rate of the Poissonian spike generator (in spikes/s)
} }
# Neuron parameters # Neuron parameters
neuron_params = { neuron_params = {
'V0_mean': -150., # Mean for the distribution of initial membrane potentials, by default: -100. 'V0_mean': -150., # Mean for the distribution of initial membrane potentials, by default: -100.
'V0_sd': 50.} # Standard deviation for the distribution of initial membrane potentials, by default: 50. 'V0_sd': 50.} # Standard deviation for the distribution of initial membrane potentials, by default: 50.
# Network parameters # Network parameters
network_params = { network_params = {
'N_scaling': scale_down_to, # Scaling of population sizes, by default: 1. 'N_scaling': scale_down_to, # Scaling of population sizes, by default: 1.
'K_scaling': scale_down_to, # Scaling of indegrees, by default: 1. 'K_scaling': scale_down_to, # Scaling of indegrees, by default: 1.
'fullscale_rates': 'tests/fullscale_rates.json', # Absolute path to the file holding full-scale rates for scaling synaptic weights, by default: None 'fullscale_rates': 'tests/fullscale_rates.json', # Absolute path to the file holding full-scale rates for scaling synaptic weights, by default: None
'input_params': input_params, # Input parameters 'input_params': input_params, # Input parameters
'connection_params': conn_params, # Connection parameters 'connection_params': conn_params, # Connection parameters
'neuron_params': neuron_params # Neuron parameters 'neuron_params': neuron_params # Neuron parameters
} }
# Simulation parameters # Simulation parameters
sim_params = { sim_params = {
'areas_simulated': areas_simulated, 'areas_simulated': areas_simulated,
't_sim': 2000., # Simulated time (in ms), by default: 10.0 't_sim': 2000., # Simulated time (in ms), by default: 10.0
'num_processes': 1, # The number of MPI processes, by default: 1 'num_processes': 1, # The number of MPI processes, by default: 1
'local_num_threads': 1, # The number of threads per MPI process, by default: 1 'local_num_threads': 1, # The number of threads per MPI process, by default: 1
'recording_dict': {'record_vm': False}, 'recording_dict': {'record_vm': False},
'rng_seed': 1 # global random seed 'rng_seed': 1 # global random seed
} }
# Theory paramters (theory_params) # Theory paramters (theory_params)
theory_params = { theory_params = {
'dt': 0.1 # The time step of the mean-field theory integration, by default: 0.01 'dt': 0.1 # The time step of the mean-field theory integration, by default: 0.01
} }
``` ```
%% Cell type:markdown id:1472e9c5 tags: %% Cell type:markdown id:1472e9c5 tags:
Go back to [Notebook structure](#toc) Go back to [Notebook structure](#toc)
%% Cell type:markdown id:c532a861-824f-4713-a311-590aef8b6134 tags: %% Cell type:markdown id:c532a861-824f-4713-a311-590aef8b6134 tags:
<br> <br>
%% Cell type:markdown id:de4a6703 tags: %% Cell type:markdown id:de4a6703 tags:
## S2. Multi-area model instantiation and simulation <a class="anchor" id="section_2"></a> ## S2. Multi-area model instantiation and simulation <a class="anchor" id="section_2"></a>
%% Cell type:markdown id:1fd58841 tags: %% Cell type:markdown id:1fd58841 tags:
### 2.1. Insantiate a multi-area model <a class="anchor" id="section_2_1"></a> ### 2.1. Insantiate a multi-area model <a class="anchor" id="section_2_1"></a>
%% Cell type:code id:ab25f9f8 tags: %% Cell type:code id:ab25f9f8 tags:
``` python ``` python
M = MultiAreaModel(network_params, M = MultiAreaModel(network_params,
simulation=True, simulation=True,
sim_spec=sim_params, sim_spec=sim_params,
theory=True, theory=True,
theory_spec=theory_params) theory_spec=theory_params)
``` ```
%% Output %% Output
Initializing network from dictionary. Initializing network from dictionary.
RAND_DATA_LABEL 588 RAND_DATA_LABEL 588
/srv/main-spack-instance-2302/spack/opt/spack/linux-ubuntu20.04-x86_64/gcc-10.3.0/py-numpy-1.21.6-6fewtq7oarp3vtwlxrrcofz5sxwt55s7/lib/python3.8/site-packages/numpy/core/fromnumeric.py:3440: RuntimeWarning:Mean of empty slice. /srv/main-spack-instance-2302/spack/opt/spack/linux-ubuntu20.04-x86_64/gcc-10.3.0/py-numpy-1.21.6-6fewtq7oarp3vtwlxrrcofz5sxwt55s7/lib/python3.8/site-packages/numpy/core/fromnumeric.py:3440: RuntimeWarning:Mean of empty slice.
/srv/main-spack-instance-2302/spack/opt/spack/linux-ubuntu20.04-x86_64/gcc-10.3.0/py-numpy-1.21.6-6fewtq7oarp3vtwlxrrcofz5sxwt55s7/lib/python3.8/site-packages/numpy/core/_methods.py:189: RuntimeWarning:invalid value encountered in double_scalars /srv/main-spack-instance-2302/spack/opt/spack/linux-ubuntu20.04-x86_64/gcc-10.3.0/py-numpy-1.21.6-6fewtq7oarp3vtwlxrrcofz5sxwt55s7/lib/python3.8/site-packages/numpy/core/_methods.py:189: RuntimeWarning:invalid value encountered in double_scalars
Error in library("aod") : there is no package called ‘aod’ Error in library("aod") : there is no package called ‘aod’
Execution halted Execution halted
No R installation or IndexError, taking hard-coded SLN fit parameters. No R installation or IndexError, taking hard-coded SLN fit parameters.
======================================== ========================================
Customized parameters Customized parameters
-------------------- --------------------
{'K_scaling': 0.005, {'K_scaling': 0.005,
'N_scaling': 0.005, 'N_scaling': 0.005,
'connection_params': {'K_stable': 'K_stable.npy', 'connection_params': {'K_stable': 'K_stable.npy',
'av_indegree_V1': 3950.0, 'av_indegree_V1': 3950.0,
'fac_nu_ext_5E': 1.125, 'fac_nu_ext_5E': 1.125,
'fac_nu_ext_6E': 1.41666667, 'fac_nu_ext_6E': 1.41666667,
'fac_nu_ext_TH': 1.2, 'fac_nu_ext_TH': 1.2,
'g': -11.0, 'g': -11.0,
'replace_non_simulated_areas': 'het_poisson_stat'}, 'replace_non_simulated_areas': 'het_poisson_stat'},
'fullscale_rates': 'tests/fullscale_rates.json', 'fullscale_rates': 'tests/fullscale_rates.json',
'input_params': {'rate_ext': 10.0}, 'input_params': {'rate_ext': 10.0},
'neuron_params': {'V0_mean': -150.0, 'V0_sd': 50.0}} 'neuron_params': {'V0_mean': -150.0, 'V0_sd': 50.0}}
======================================== ========================================
/srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/dicthash/dicthash.py:47: UserWarning:Float too small for safe conversion tointeger. Rounding down to zero. /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/dicthash/dicthash.py:47: UserWarning:Float too small for safe conversion tointeger. Rounding down to zero.
Simulation label: 27d81076e6d6e9e591684be053078477 Simulation label: 27d81076e6d6e9e591684be053078477
Copied files. Copied files.
Initialized simulation class. Initialized simulation class.
%% Cell type:markdown id:91649c30 tags: %% Cell type:markdown id:91649c30 tags:
### 2.2. Predict firing rates from theory <a class="anchor" id="section_2_2"></a> ### 2.2. Predict firing rates from theory <a class="anchor" id="section_2_2"></a>
%% Cell type:code id:6a7ddf0e tags: %% Cell type:code id:6a7ddf0e tags:
``` python ``` python
p, r = M.theory.integrate_siegert() p, r = M.theory.integrate_siegert()
print("Mean-field theory predicts an average " print("Mean-field theory predicts an average "
"firing rate of {0:.3f} spikes/s across all populations.".format(np.mean(r[:, -1]))) "firing rate of {0:.3f} spikes/s across all populations.".format(np.mean(r[:, -1])))
``` ```
%% Output %% Output
Iteration: 0 Iteration: 0
Mean-field theory predicts an average firing rate of 29.588 spikes/s across all populations. Mean-field theory predicts an average firing rate of 29.588 spikes/s across all populations.
%% Cell type:markdown id:2062ddf3 tags: %% Cell type:markdown id:2062ddf3 tags:
### 2.3. Extract interarea connectivity <a class="anchor" id="section_2_3"></a> ### 2.3. Extract interarea connectivity <a class="anchor" id="section_2_3"></a>
%% Cell type:markdown id:8a7c09e0 tags: %% Cell type:markdown id:8a7c09e0 tags:
The connectivity and neuron numbers are stored in the attributes of the model class. Neuron numbers are stored in `M.N` as a dictionary (and in `M.N_vec` as an array), indegrees in `M.K` as a dictionary (and in `M.K_matrix` as an array). Number of synapses can also be access via `M.synapses` (and in `M.syn_matrix` as an array). <br> The connectivity and neuron numbers are stored in the attributes of the model class. Neuron numbers are stored in `M.N` as a dictionary (and in `M.N_vec` as an array), indegrees in `M.K` as a dictionary (and in `M.K_matrix` as an array). Number of synapses can also be access via `M.synapses` (and in `M.syn_matrix` as an array). <br>
%% Cell type:markdown id:b7396606 tags: %% Cell type:markdown id:b7396606 tags:
#### 2.3.1 Node indegrees #### 2.3.1 Node indegrees
%% Cell type:code id:6316ac24 tags: %% Cell type:code id:6316ac24 tags:
``` python ``` python
# Dictionary of nodes indegrees organized as: # Dictionary of nodes indegrees organized as:
# {<source_area>: {<source_pop>: {<target_area>: {<target_pop>: indegree_values}}}} # {<source_area>: {<source_pop>: {<target_area>: {<target_pop>: indegree_values}}}}
# M.K # M.K
``` ```
%% Cell type:markdown id:253a2aba tags: %% Cell type:markdown id:253a2aba tags:
#### 2.3.2 Synapses #### 2.3.2 Synapses
%% Cell type:code id:445a722a tags: %% Cell type:code id:445a722a tags:
``` python ``` python
# Dictionary of synapses that target neurons receive, it is organized as: # Dictionary of synapses that target neurons receive, it is organized as:
# {<source_area>: {<source_pop>: {<target_area>: {<target_pop>: number_of_synapses}}}} # {<source_area>: {<source_pop>: {<target_area>: {<target_pop>: number_of_synapses}}}}
# M.synapses # M.synapses
``` ```
%% Cell type:markdown id:e67f37e9-ec8d-4bb1-bd21-45e966f47ab6 tags: %% Cell type:markdown id:e67f37e9-ec8d-4bb1-bd21-45e966f47ab6 tags:
Go back to [Notebook structure](#toc) Go back to [Notebook structure](#toc)
%% Cell type:markdown id:0c1cad59-81d0-4e24-ac33-13c4ca8c6dec tags: %% Cell type:markdown id:0c1cad59-81d0-4e24-ac33-13c4ca8c6dec tags:
### 2.4. Run the simulation <a class="anchor" id="section_2_4"></a> ### 2.4. Run the simulation <a class="anchor" id="section_2_4"></a>
%% Cell type:code id:15778e9c tags: %% Cell type:code id:15778e9c tags:
``` python ``` python
# run the simulation, depending on the model parameter and downscale ratio, the running time varies largely. # run the simulation, depending on the model parameter and downscale ratio, the running time varies largely.
M.simulation.simulate() M.simulation.simulate()
``` ```
%% Output %% Output
Prepared simulation in 0.00 seconds. Prepared simulation in 0.00 seconds.
Rank 0: created area V1 with 0 local nodes Rank 0: created area V1 with 0 local nodes
Memory after V1 : 1911.94 MB Memory after V1 : 1911.94 MB
Rank 0: created area V2 with 0 local nodes Rank 0: created area V2 with 0 local nodes
Memory after V2 : 1938.54 MB Memory after V2 : 1938.54 MB
Rank 0: created area VP with 0 local nodes Rank 0: created area VP with 0 local nodes
Memory after VP : 1967.84 MB Memory after VP : 1967.84 MB
Rank 0: created area V3 with 0 local nodes Rank 0: created area V3 with 0 local nodes
Memory after V3 : 1996.10 MB Memory after V3 : 1996.10 MB
Rank 0: created area V3A with 0 local nodes Rank 0: created area V3A with 0 local nodes
Memory after V3A : 2016.04 MB Memory after V3A : 2016.04 MB
Rank 0: created area MT with 0 local nodes Rank 0: created area MT with 0 local nodes
Memory after MT : 2041.55 MB Memory after MT : 2041.55 MB
Rank 0: created area V4t with 0 local nodes Rank 0: created area V4t with 0 local nodes
Memory after V4t : 2066.49 MB Memory after V4t : 2066.49 MB
Rank 0: created area V4 with 0 local nodes Rank 0: created area V4 with 0 local nodes
Memory after V4 : 2093.56 MB Memory after V4 : 2093.56 MB
Rank 0: created area VOT with 0 local nodes Rank 0: created area VOT with 0 local nodes
Memory after VOT : 2118.79 MB Memory after VOT : 2118.79 MB
Rank 0: created area MSTd with 0 local nodes Rank 0: created area MSTd with 0 local nodes
Memory after MSTd : 2140.26 MB Memory after MSTd : 2140.26 MB
Rank 0: created area PIP with 0 local nodes Rank 0: created area PIP with 0 local nodes
Memory after PIP : 2161.65 MB Memory after PIP : 2161.65 MB
Rank 0: created area PO with 0 local nodes Rank 0: created area PO with 0 local nodes
Memory after PO : 2183.16 MB Memory after PO : 2183.16 MB
Rank 0: created area DP with 0 local nodes Rank 0: created area DP with 0 local nodes
Memory after DP : 2203.39 MB Memory after DP : 2203.39 MB
Rank 0: created area MIP with 0 local nodes Rank 0: created area MIP with 0 local nodes
Memory after MIP : 2224.92 MB Memory after MIP : 2224.92 MB
Rank 0: created area MDP with 0 local nodes Rank 0: created area MDP with 0 local nodes
Memory after MDP : 2246.44 MB Memory after MDP : 2246.44 MB
Rank 0: created area VIP with 0 local nodes Rank 0: created area VIP with 0 local nodes
Memory after VIP : 2268.37 MB Memory after VIP : 2268.37 MB
Rank 0: created area LIP with 0 local nodes Rank 0: created area LIP with 0 local nodes
Memory after LIP : 2292.32 MB Memory after LIP : 2292.32 MB
Rank 0: created area PITv with 0 local nodes Rank 0: created area PITv with 0 local nodes
Memory after PITv : 2317.62 MB Memory after PITv : 2317.62 MB
Rank 0: created area PITd with 0 local nodes Rank 0: created area PITd with 0 local nodes
Memory after PITd : 2342.83 MB Memory after PITd : 2342.83 MB
Rank 0: created area MSTl with 0 local nodes Rank 0: created area MSTl with 0 local nodes
Memory after MSTl : 2364.29 MB Memory after MSTl : 2364.29 MB
Rank 0: created area CITv with 0 local nodes Rank 0: created area CITv with 0 local nodes
Memory after CITv : 2383.51 MB Memory after CITv : 2383.51 MB
Rank 0: created area CITd with 0 local nodes Rank 0: created area CITd with 0 local nodes
Memory after CITd : 2402.80 MB Memory after CITd : 2402.80 MB
Rank 0: created area FEF with 0 local nodes Rank 0: created area FEF with 0 local nodes
Memory after FEF : 2424.31 MB Memory after FEF : 2424.31 MB
Rank 0: created area TF with 0 local nodes Rank 0: created area TF with 0 local nodes
Memory after TF : 2439.88 MB Memory after TF : 2439.88 MB
Rank 0: created area AITv with 0 local nodes Rank 0: created area AITv with 0 local nodes
Memory after AITv : 2462.55 MB Memory after AITv : 2462.55 MB
Rank 0: created area FST with 0 local nodes Rank 0: created area FST with 0 local nodes
Memory after FST : 2479.28 MB Memory after FST : 2479.28 MB
Rank 0: created area 7a with 0 local nodes Rank 0: created area 7a with 0 local nodes
Memory after 7a : 2500.42 MB Memory after 7a : 2500.42 MB
Rank 0: created area STPp with 0 local nodes Rank 0: created area STPp with 0 local nodes
Memory after STPp : 2519.25 MB Memory after STPp : 2519.25 MB
Rank 0: created area STPa with 0 local nodes Rank 0: created area STPa with 0 local nodes
Memory after STPa : 2538.40 MB Memory after STPa : 2538.40 MB
Rank 0: created area 46 with 0 local nodes Rank 0: created area 46 with 0 local nodes
Memory after 46 : 2553.76 MB Memory after 46 : 2553.76 MB
Rank 0: created area AITd with 0 local nodes Rank 0: created area AITd with 0 local nodes
Memory after AITd : 2576.32 MB Memory after AITd : 2576.32 MB
Rank 0: created area TH with 0 local nodes Rank 0: created area TH with 0 local nodes
Memory after TH : 2589.07 MB Memory after TH : 2589.07 MB
Created areas and internal connections in 2.19 seconds. Created areas and internal connections in 2.19 seconds.
Created cortico-cortical connections in 22.42 seconds. Created cortico-cortical connections in 22.42 seconds.
Simulated network in 74.96 seconds. Simulated network in 74.96 seconds.
%% Cell type:markdown id:fd6e3232 tags: %% Cell type:markdown id:fd6e3232 tags:
Go back to [Notebook structure](#toc) Go back to [Notebook structure](#toc)
%% Cell type:markdown id:4003c5a5-4a6f-49c5-be17-09f1bc68c411 tags: %% Cell type:markdown id:4003c5a5-4a6f-49c5-be17-09f1bc68c411 tags:
<br> <br>
%% Cell type:markdown id:28e071f8 tags: %% Cell type:markdown id:28e071f8 tags:
## S3. Simulation results validation and connection extraction <a class="anchor" id="section_3"></a> ## S3. Simulation results validation and connection extraction <a class="anchor" id="section_3"></a>
%% Cell type:markdown id:89c7b7cf tags: %% Cell type:markdown id:89c7b7cf tags:
### 3.1 Test if the correct number of synapses has been created ### 3.1 Test if the correct number of synapses has been created
%% Cell type:code id:dc3b1820 tags: %% Cell type:code id:dc3b1820 tags:
``` python ``` python
# # Uncomment the lines in this code cell below to test if the number of synapses created by NEST matches the expected values # # Uncomment the lines in this code cell below to test if the number of synapses created by NEST matches the expected values
# print("Testing synapse numbers") # print("Testing synapse numbers")
# for target_area_name in M.area_list: # for target_area_name in M.area_list:
# target_area = M.simulation.areas[M.simulation.areas.index(target_area_name)] # target_area = M.simulation.areas[M.simulation.areas.index(target_area_name)]
# for source_area_name in M.area_list: # for source_area_name in M.area_list:
# source_area = M.simulation.areas[M.simulation.areas.index(source_area_name)] # source_area = M.simulation.areas[M.simulation.areas.index(source_area_name)]
# for target_pop in M.structure[target_area.name]: # for target_pop in M.structure[target_area.name]:
# target_nodes = target_area.gids[target_pop] # target_nodes = target_area.gids[target_pop]
# for source_pop in M.structure[source_area.name]: # for source_pop in M.structure[source_area.name]:
# source_nodes = source_area.gids[source_pop] # source_nodes = source_area.gids[source_pop]
# created_syn = nest.GetConnections(source=source_nodes, # created_syn = nest.GetConnections(source=source_nodes,
# target=target_nodes) # target=target_nodes)
# syn = M.synapses[target_area.name][target_pop][source_area.name][source_pop] # syn = M.synapses[target_area.name][target_pop][source_area.name][source_pop]
# assert(len(created_syn) == int(syn)) # assert(len(created_syn) == int(syn))
``` ```
%% Cell type:markdown id:57401110 tags: %% Cell type:markdown id:57401110 tags:
### 3.2 Extract connections information ### 3.2 Extract connections information
**Warning**: Memory explosion <br> **Warning**: Memory explosion <br>
To obtain the connections information, you can extract the lists of connected sources and targets. Moreover, you can access additional synaptic details, such as synaptic weights and delays. To obtain the connections information, you can extract the lists of connected sources and targets. Moreover, you can access additional synaptic details, such as synaptic weights and delays.
%% Cell type:code id:e7eb052e tags: %% Cell type:code id:e7eb052e tags:
``` python ``` python
# conns = nest.GetConnections() # conns = nest.GetConnections()
# conns_sparse_matrix = conns.get(['source', 'target', 'weight']) # conns_sparse_matrix = conns.get(['source', 'target', 'weight'])
# srcs = conns_sparse_matrix['source'] # srcs = conns_sparse_matrix['source']
# tgts = conns_sparse_matrix['target'] # tgts = conns_sparse_matrix['target']
# weights = conns_sparse_matrix['weight'] # weights = conns_sparse_matrix['weight']
``` ```
%% Cell type:markdown id:ef4b2e4b tags: %% Cell type:markdown id:ef4b2e4b tags:
You can determine the area and subpopulation to which the neuron ID ranges belong by referring to the file `network_gids.txt`, which is automatically generated during network creation. You can determine the area and subpopulation to which the neuron ID ranges belong by referring to the file `network_gids.txt`, which is automatically generated during network creation.
%% Cell type:code id:902f2800 tags: %% Cell type:code id:902f2800 tags:
``` python ``` python
# # Open the file using a with statement # # Open the file using a with statement
# with open(os.path.join(M.simulation.data_dir,"recordings/network_gids.txt"), "r") as file: # with open(os.path.join(M.simulation.data_dir,"recordings/network_gids.txt"), "r") as file:
# # Read the contents of the file # # Read the contents of the file
# gids = file.read() # gids = file.read()
# # Print the contents # # Print the contents
# print(gids) # print(gids)
``` ```
%% Cell type:markdown id:b1320ab1 tags: %% Cell type:markdown id:b1320ab1 tags:
Go back to [Notebook structure](#toc) Go back to [Notebook structure](#toc)
%% Cell type:markdown id:529b1ade tags: %% Cell type:markdown id:529b1ade tags:
<br> <br>
%% Cell type:markdown id:57ff902c-d6ce-4f96-9e4f-8e3e7166ab66 tags: %% Cell type:markdown id:57ff902c-d6ce-4f96-9e4f-8e3e7166ab66 tags:
## S4. Data loading and processing <a class="anchor" id="section_4"></a> ## S4. Data loading and processing <a class="anchor" id="section_4"></a>
%% Cell type:code id:f5b58845-4d1a-430f-83f4-402fdf918aef tags: %% Cell type:code id:f5b58845-4d1a-430f-83f4-402fdf918aef tags:
``` python ``` python
label_spikes = M.simulation.label label_spikes = M.simulation.label
label = M.simulation.label label = M.simulation.label
``` ```
%% Cell type:code id:6607a73d-1c74-4848-9603-081ad0e7cae8 tags: %% Cell type:code id:6607a73d-1c74-4848-9603-081ad0e7cae8 tags:
``` python ``` python
""" """
Analysis class. Analysis class.
An instance of the analysis class for the given network and simulation. An instance of the analysis class for the given network and simulation.
Can be created as a member class of a multiarea_model instance or standalone. Can be created as a member class of a multiarea_model instance or standalone.
Parameters Parameters
---------- ----------
network : MultiAreaModel network : MultiAreaModel
An instance of the multiarea_model class that specifies An instance of the multiarea_model class that specifies
the network to be analyzed. the network to be analyzed.
simulation : Simulation simulation : Simulation
An instance of the simulation class that specifies An instance of the simulation class that specifies
the simulation to be analyzed. the simulation to be analyzed.
data_list : list of strings {'spikes', vm'}, optional data_list : list of strings {'spikes', vm'}, optional
Specifies which type of data is to load. Defaults to ['spikes']. Specifies which type of data is to load. Defaults to ['spikes'].
load_areas : list of strings with area names, optional load_areas : list of strings with area names, optional
Specifies the areas for which data is to be loaded. Specifies the areas for which data is to be loaded.
Default value is None and leads to loading of data for all Default value is None and leads to loading of data for all
simulated areas. simulated areas.
""" """
# Instantiate an analysis class and load spike data # Instantiate an analysis class and load spike data
A = Analysis(network=M, A = Analysis(network=M,
simulation=M.simulation, simulation=M.simulation,
data_list=['spikes'], data_list=['spikes'],
load_areas=None) load_areas=None)
``` ```
%% Output %% Output
loading spikes loading spikes
%% Cell type:code id:c471e9c8-b1e4-43e4-a6a1-443b8b8963be tags: %% Cell type:code id:c471e9c8-b1e4-43e4-a6a1-443b8b8963be tags:
``` python ``` python
# load spike data and calculate instantaneous and mean firing rates # load spike data and calculate instantaneous and mean firing rates
data = np.loadtxt(M.simulation.data_dir + '/recordings/' + M.simulation.label + "-spikes-1-0.dat", skiprows=3) data = np.loadtxt(M.simulation.data_dir + '/recordings/' + M.simulation.label + "-spikes-1-0.dat", skiprows=3)
tsteps, spikecount = np.unique(data[:,1], return_counts=True) tsteps, spikecount = np.unique(data[:,1], return_counts=True)
rate = spikecount / M.simulation.params['dt'] * 1e3 / np.sum(M.N_vec) rate = spikecount / M.simulation.params['dt'] * 1e3 / np.sum(M.N_vec)
``` ```
%% Cell type:code id:1870cf34-ee62-4614-bc25-c36bc9a7377c tags: %% Cell type:code id:1870cf34-ee62-4614-bc25-c36bc9a7377c tags:
``` python ``` python
""" """
Calculate time-averaged population rates and store them in member pop_rates. Calculate time-averaged population rates and store them in member pop_rates.
If the rates had previously been stored with the same If the rates had previously been stored with the same
parameters, they are loaded from file. parameters, they are loaded from file.
Parameters Parameters
---------- ----------
t_min : float, optional t_min : float, optional
Minimal time in ms of the simulation to take into account Minimal time in ms of the simulation to take into account
for the calculation. Defaults to 500 ms. for the calculation. Defaults to 500 ms.
t_max : float, optional t_max : float, optional
Maximal time in ms of the simulation to take into account Maximal time in ms of the simulation to take into account
for the calculation. Defaults to the simulation time. for the calculation. Defaults to the simulation time.
compute_stat : bool, optional compute_stat : bool, optional
If set to true, the mean and variance of the population rate If set to true, the mean and variance of the population rate
is calculated. Defaults to False. is calculated. Defaults to False.
Caution: Setting to True slows down the computation. Caution: Setting to True slows down the computation.
areas : list, optional areas : list, optional
Which areas to include in the calculcation. Which areas to include in the calculcation.
Defaults to all loaded areas. Defaults to all loaded areas.
pops : list or {'complete'}, optional pops : list or {'complete'}, optional
Which populations to include in the calculation. Which populations to include in the calculation.
If set to 'complete', all populations the respective areas If set to 'complete', all populations the respective areas
are included. Defaults to 'complete'. are included. Defaults to 'complete'.
""" """
A.create_pop_rates() A.create_pop_rates()
print("Computing population rates done") print("Computing population rates done")
``` ```
%% Output %% Output
Loading data from file Loading data from file
Computing population rates done Computing population rates done
%% Cell type:code id:50b7df89 tags: %% Cell type:code id:50b7df89 tags:
``` python ``` python
""" """
Calculate synchrony as the coefficient of variation of the population rate Calculate synchrony as the coefficient of variation of the population rate
and store in member synchrony. Uses helper function synchrony. and store in member synchrony. Uses helper function synchrony.
If the synchrony has previously been stored with the If the synchrony has previously been stored with the
same parameters, they are loaded from file. same parameters, they are loaded from file.
Parameters Parameters
---------- ----------
t_min : float, optional t_min : float, optional
Minimal time in ms of the simulation to take into account Minimal time in ms of the simulation to take into account
for the calculation. Defaults to 500 ms. for the calculation. Defaults to 500 ms.
t_max : float, optional t_max : float, optional
Maximal time in ms of the simulation to take into account Maximal time in ms of the simulation to take into account
for the calculation. Defaults to the simulation time. for the calculation. Defaults to the simulation time.
areas : list, optional areas : list, optional
Which areas to include in the calculcation. Which areas to include in the calculcation.
Defaults to all loaded areas. Defaults to all loaded areas.
pops : list or {'complete'}, optional pops : list or {'complete'}, optional
Which populations to include in the calculation. Which populations to include in the calculation.
If set to 'complete', all populations the respective areas If set to 'complete', all populations the respective areas
are included. Defaults to 'complete'. are included. Defaults to 'complete'.
resolution : float, optional resolution : float, optional
Resolution of the population rate. Defaults to 1 ms. Resolution of the population rate. Defaults to 1 ms.
""" """
A.create_synchrony() A.create_synchrony()
print("Computing synchrony done") print("Computing synchrony done")
``` ```
%% Output %% Output
Loading data from file Loading data from file
Computing synchrony done Computing synchrony done
%% Cell type:code id:d43b493c tags: %% Cell type:code id:d43b493c tags:
``` python ``` python
""" """
Calculate poulation-averaged LvR (see Shinomoto et al. 2009) and Calculate poulation-averaged LvR (see Shinomoto et al. 2009) and
store as member pop_LvR. Uses helper function LvR. store as member pop_LvR. Uses helper function LvR.
Parameters Parameters
---------- ----------
t_min : float, optional t_min : float, optional
Minimal time in ms of the simulation to take into account Minimal time in ms of the simulation to take into account
for the calculation. Defaults to 500 ms. for the calculation. Defaults to 500 ms.
t_max : float, optional t_max : float, optional
Maximal time in ms of the simulation to take into account Maximal time in ms of the simulation to take into account
for the calculation. Defaults to the simulation time. for the calculation. Defaults to the simulation time.
areas : list, optional areas : list, optional
Which areas to include in the calculcation. Which areas to include in the calculcation.
Defaults to all loaded areas. Defaults to all loaded areas.
pops : list or {'complete'}, optional pops : list or {'complete'}, optional
Which populations to include in the calculation. Which populations to include in the calculation.
If set to 'complete', all populations the respective areas If set to 'complete', all populations the respective areas
are included. Defaults to 'complete'. are included. Defaults to 'complete'.
""" """
A.create_pop_LvR() A.create_pop_LvR()
print("Computing population LvR done") print("Computing population LvR done")
``` ```
%% Output %% Output
Loading data from file Loading data from file
Computing population LvR done Computing population LvR done
%% Cell type:code id:401ece2d-47c8-4775-80ae-92a8e432520c tags: %% Cell type:code id:401ece2d-47c8-4775-80ae-92a8e432520c tags:
``` python ``` python
""" """
Calculate time series of population- and area-averaged firing rates. Calculate time series of population- and area-averaged firing rates.
Uses ah.pop_rate_time_series. Uses ah.pop_rate_time_series.
If the rates have previously been stored with the If the rates have previously been stored with the
same parameters, they are loaded from file. same parameters, they are loaded from file.
Parameters Parameters
---------- ----------
t_min : float, optional t_min : float, optional
Minimal time in ms of the simulation to take into account Minimal time in ms of the simulation to take into account
for the calculation. Defaults to 500 ms. for the calculation. Defaults to 500 ms.
t_max : float, optional t_max : float, optional
Maximal time in ms of the simulation to take into account Maximal time in ms of the simulation to take into account
for the calculation. Defaults to the simulation time. for the calculation. Defaults to the simulation time.
areas : list, optional areas : list, optional
Which areas to include in the calculcation. Which areas to include in the calculcation.
Defaults to all loaded areas. Defaults to all loaded areas.
pops : list or {'complete'}, optional pops : list or {'complete'}, optional
Which populations to include in the calculation. Which populations to include in the calculation.
If set to 'complete', all populations the respective areas If set to 'complete', all populations the respective areas
are included. Defaults to 'complete'. are included. Defaults to 'complete'.
kernel : {'gauss_time_window', 'alpha_time_window', 'rect_time_window'}, optional kernel : {'gauss_time_window', 'alpha_time_window', 'rect_time_window'}, optional
Specifies the kernel to be convolved with the spike histogram. Specifies the kernel to be convolved with the spike histogram.
Defaults to 'binned', which corresponds to no convolution. Defaults to 'binned', which corresponds to no convolution.
resolution: float, optional resolution: float, optional
Width of the convolution kernel. Specifically it correponds to: Width of the convolution kernel. Specifically it correponds to:
- 'binned' : bin width of the histogram - 'binned' : bin width of the histogram
- 'gauss_time_window' : sigma - 'gauss_time_window' : sigma
- 'alpha_time_window' : time constant of the alpha function - 'alpha_time_window' : time constant of the alpha function
- 'rect_time_window' : width of the moving rectangular function - 'rect_time_window' : width of the moving rectangular function
""" """
A.create_rate_time_series() A.create_rate_time_series()
print("Computing rate time series done") print("Computing rate time series done")
``` ```
%% Output %% Output
Loading data from file Loading data from file
Loading data from file Loading data from file
Computing rate time series done Computing rate time series done
%% Cell type:code id:fa3ea20e-e456-4608-a711-e2c320bcaf91 tags: %% Cell type:code id:fa3ea20e-e456-4608-a711-e2c320bcaf91 tags:
``` python ``` python
A.save() A.save()
``` ```
%% Output %% Output
pop_LvR pop_LvR
pop_rates pop_rates
synchrony synchrony
%% Cell type:markdown id:2da9728d-4481-4a15-b810-d125e39cbe4e tags: %% Cell type:markdown id:2da9728d-4481-4a15-b810-d125e39cbe4e tags:
Go back to [Notebook structure](#toc) Go back to [Notebook structure](#toc)
%% Cell type:markdown id:4d43d223-a62e-448a-a7ea-8379b8be8e86 tags: %% Cell type:markdown id:4d43d223-a62e-448a-a7ea-8379b8be8e86 tags:
<br> <br>
%% Cell type:markdown id:bb71c922 tags: %% Cell type:markdown id:bb71c922 tags:
## S5. Simulation results visualziation <a class="anchor" id="section_5"></a> ## S5. Simulation results visualziation <a class="anchor" id="section_5"></a>
%% Cell type:code id:3659d4de-cb23-4980-bddb-0245cd4d1bb7 tags: %% Cell type:code id:3659d4de-cb23-4980-bddb-0245cd4d1bb7 tags:
``` python ``` python
# # visualization settings # # visualization settings
# icolor = 'r' # icolor = 'r'
# ecolor = 'b' # ecolor = 'b'
# population_labels = ['2/3E', '2/3I', '4E', '4I', '5E', '5I', '6E', '6I'] # population_labels = ['2/3E', '2/3I', '4E', '4I', '5E', '5I', '6E', '6I']
``` ```
%% Cell type:markdown id:38ddd973 tags: %% Cell type:markdown id:38ddd973 tags:
### 5.1. Instantaneous and mean firing rate across all populations <a class="anchor" id="section_5_1"></a> ### 5.1. Instantaneous and mean firing rate across all populations <a class="anchor" id="section_5_1"></a>
%% Cell type:code id:bea30fc8 tags: %% Cell type:code id:bea30fc8 tags:
``` python ``` python
plot_instan_mean_firing_rate(tsteps, rate, sim_params) # plot_instan_mean_firing_rate(tsteps, rate, sim_params)
``` ```
%% Output %% Output
--------------------------------------------------------------------------- ---------------------------------------------------------------------------
CalledProcessError Traceback (most recent call last) CalledProcessError Traceback (most recent call last)
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/texmanager.py:233, in TexManager._run_checked_subprocess(self, command, tex, cwd) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/texmanager.py:233, in TexManager._run_checked_subprocess(self, command, tex, cwd)
232 try: 232 try:
--> 233 report = subprocess.check_output( --> 233 report = subprocess.check_output(
234 command, cwd=cwd if cwd is not None else self.texcache, 234 command, cwd=cwd if cwd is not None else self.texcache,
235 stderr=subprocess.STDOUT) 235 stderr=subprocess.STDOUT)
236 except FileNotFoundError as exc: 236 except FileNotFoundError as exc:
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/subprocess.py:415, in check_output(timeout, *popenargs, **kwargs) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/subprocess.py:415, in check_output(timeout, *popenargs, **kwargs)
413 kwargs['input'] = empty 413 kwargs['input'] = empty
--> 415 return run(*popenargs, stdout=PIPE, timeout=timeout, check=True, --> 415 return run(*popenargs, stdout=PIPE, timeout=timeout, check=True,
416 **kwargs).stdout 416 **kwargs).stdout
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/subprocess.py:516, in run(input, capture_output, timeout, check, *popenargs, **kwargs) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/._view/6axslmv6jvf4v2nte3uwlayg4vhsjoha/lib/python3.8/subprocess.py:516, in run(input, capture_output, timeout, check, *popenargs, **kwargs)
515 if check and retcode: 515 if check and retcode:
--> 516 raise CalledProcessError(retcode, process.args, --> 516 raise CalledProcessError(retcode, process.args,
517 output=stdout, stderr=stderr) 517 output=stdout, stderr=stderr)
518 return CompletedProcess(process.args, retcode, stdout, stderr) 518 return CompletedProcess(process.args, retcode, stdout, stderr)
CalledProcessError: Command '['latex', '-interaction=nonstopmode', '--halt-on-error', '../540885da22f8d5147f0088a11e43701a.tex']' returned non-zero exit status 1. CalledProcessError: Command '['latex', '-interaction=nonstopmode', '--halt-on-error', '../540885da22f8d5147f0088a11e43701a.tex']' returned non-zero exit status 1.
The above exception was the direct cause of the following exception: The above exception was the direct cause of the following exception:
RuntimeError Traceback (most recent call last) RuntimeError Traceback (most recent call last)
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/IPython/core/formatters.py:339, in BaseFormatter.__call__(self, obj) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/IPython/core/formatters.py:339, in BaseFormatter.__call__(self, obj)
337 pass 337 pass
338 else: 338 else:
--> 339 return printer(obj) --> 339 return printer(obj)
340 # Finally look for special method names 340 # Finally look for special method names
341 method = get_real_method(obj, self.print_method) 341 method = get_real_method(obj, self.print_method)
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/IPython/core/pylabtools.py:151, in print_figure(fig, fmt, bbox_inches, base64, **kwargs) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/IPython/core/pylabtools.py:151, in print_figure(fig, fmt, bbox_inches, base64, **kwargs)
148 from matplotlib.backend_bases import FigureCanvasBase 148 from matplotlib.backend_bases import FigureCanvasBase
149 FigureCanvasBase(fig) 149 FigureCanvasBase(fig)
--> 151 fig.canvas.print_figure(bytes_io, **kw) --> 151 fig.canvas.print_figure(bytes_io, **kw)
152 data = bytes_io.getvalue() 152 data = bytes_io.getvalue()
153 if fmt == 'svg': 153 if fmt == 'svg':
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/backend_bases.py:2295, in FigureCanvasBase.print_figure(self, filename, dpi, facecolor, edgecolor, orientation, format, bbox_inches, pad_inches, bbox_extra_artists, backend, **kwargs) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/backend_bases.py:2295, in FigureCanvasBase.print_figure(self, filename, dpi, facecolor, edgecolor, orientation, format, bbox_inches, pad_inches, bbox_extra_artists, backend, **kwargs)
2289 renderer = _get_renderer( 2289 renderer = _get_renderer(
2290 self.figure, 2290 self.figure,
2291 functools.partial( 2291 functools.partial(
2292 print_method, orientation=orientation) 2292 print_method, orientation=orientation)
2293 ) 2293 )
2294 with getattr(renderer, "_draw_disabled", nullcontext)(): 2294 with getattr(renderer, "_draw_disabled", nullcontext)():
-> 2295 self.figure.draw(renderer) -> 2295 self.figure.draw(renderer)
2297 if bbox_inches: 2297 if bbox_inches:
2298 if bbox_inches == "tight": 2298 if bbox_inches == "tight":
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/artist.py:73, in _finalize_rasterization.<locals>.draw_wrapper(artist, renderer, *args, **kwargs) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/artist.py:73, in _finalize_rasterization.<locals>.draw_wrapper(artist, renderer, *args, **kwargs)
71 @wraps(draw) 71 @wraps(draw)
72 def draw_wrapper(artist, renderer, *args, **kwargs): 72 def draw_wrapper(artist, renderer, *args, **kwargs):
---> 73 result = draw(artist, renderer, *args, **kwargs) ---> 73 result = draw(artist, renderer, *args, **kwargs)
74 if renderer._rasterizing: 74 if renderer._rasterizing:
75 renderer.stop_rasterizing() 75 renderer.stop_rasterizing()
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/artist.py:50, in allow_rasterization.<locals>.draw_wrapper(artist, renderer) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/artist.py:50, in allow_rasterization.<locals>.draw_wrapper(artist, renderer)
47 if artist.get_agg_filter() is not None: 47 if artist.get_agg_filter() is not None:
48 renderer.start_filter() 48 renderer.start_filter()
---> 50 return draw(artist, renderer) ---> 50 return draw(artist, renderer)
51 finally: 51 finally:
52 if artist.get_agg_filter() is not None: 52 if artist.get_agg_filter() is not None:
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/figure.py:2837, in Figure.draw(self, renderer) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/figure.py:2837, in Figure.draw(self, renderer)
2834 # ValueError can occur when resizing a window. 2834 # ValueError can occur when resizing a window.
2836 self.patch.draw(renderer) 2836 self.patch.draw(renderer)
-> 2837 mimage._draw_list_compositing_images( -> 2837 mimage._draw_list_compositing_images(
2838 renderer, self, artists, self.suppressComposite) 2838 renderer, self, artists, self.suppressComposite)
2840 for sfig in self.subfigs: 2840 for sfig in self.subfigs:
2841 sfig.draw(renderer) 2841 sfig.draw(renderer)
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/image.py:132, in _draw_list_compositing_images(renderer, parent, artists, suppress_composite) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/image.py:132, in _draw_list_compositing_images(renderer, parent, artists, suppress_composite)
130 if not_composite or not has_images: 130 if not_composite or not has_images:
131 for a in artists: 131 for a in artists:
--> 132 a.draw(renderer) --> 132 a.draw(renderer)
133 else: 133 else:
134 # Composite any adjacent images together 134 # Composite any adjacent images together
135 image_group = [] 135 image_group = []
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/artist.py:50, in allow_rasterization.<locals>.draw_wrapper(artist, renderer) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/artist.py:50, in allow_rasterization.<locals>.draw_wrapper(artist, renderer)
47 if artist.get_agg_filter() is not None: 47 if artist.get_agg_filter() is not None:
48 renderer.start_filter() 48 renderer.start_filter()
---> 50 return draw(artist, renderer) ---> 50 return draw(artist, renderer)
51 finally: 51 finally:
52 if artist.get_agg_filter() is not None: 52 if artist.get_agg_filter() is not None:
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/axes/_base.py:3055, in _AxesBase.draw(self, renderer) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/axes/_base.py:3055, in _AxesBase.draw(self, renderer)
3052 for spine in self.spines.values(): 3052 for spine in self.spines.values():
3053 artists.remove(spine) 3053 artists.remove(spine)
-> 3055 self._update_title_position(renderer) -> 3055 self._update_title_position(renderer)
3057 if not self.axison: 3057 if not self.axison:
3058 for _axis in self._get_axis_list(): 3058 for _axis in self._get_axis_list():
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/axes/_base.py:3005, in _AxesBase._update_title_position(self, renderer) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/axes/_base.py:3005, in _AxesBase._update_title_position(self, renderer)
3003 _log.debug('top of Axes not in the figure, so title not moved') 3003 _log.debug('top of Axes not in the figure, so title not moved')
3004 return 3004 return
-> 3005 if title.get_window_extent(renderer).ymin < top: -> 3005 if title.get_window_extent(renderer).ymin < top:
3006 _, y = self.transAxes.inverted().transform((0, top)) 3006 _, y = self.transAxes.inverted().transform((0, top))
3007 title.set_position((x, y)) 3007 title.set_position((x, y))
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/text.py:910, in Text.get_window_extent(self, renderer, dpi) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/text.py:910, in Text.get_window_extent(self, renderer, dpi)
907 raise RuntimeError('Cannot get window extent w/o renderer') 907 raise RuntimeError('Cannot get window extent w/o renderer')
909 with cbook._setattr_cm(self.figure, dpi=dpi): 909 with cbook._setattr_cm(self.figure, dpi=dpi):
--> 910 bbox, info, descent = self._get_layout(self._renderer) --> 910 bbox, info, descent = self._get_layout(self._renderer)
911 x, y = self.get_unitless_position() 911 x, y = self.get_unitless_position()
912 x, y = self.get_transform().transform((x, y)) 912 x, y = self.get_transform().transform((x, y))
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/text.py:309, in Text._get_layout(self, renderer) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/text.py:309, in Text._get_layout(self, renderer)
306 ys = [] 306 ys = []
308 # Full vertical extent of font, including ascenders and descenders: 308 # Full vertical extent of font, including ascenders and descenders:
--> 309 _, lp_h, lp_d = renderer.get_text_width_height_descent( --> 309 _, lp_h, lp_d = renderer.get_text_width_height_descent(
310 "lp", self._fontproperties, 310 "lp", self._fontproperties,
311 ismath="TeX" if self.get_usetex() else False) 311 ismath="TeX" if self.get_usetex() else False)
312 min_dy = (lp_h - lp_d) * self._linespacing 312 min_dy = (lp_h - lp_d) * self._linespacing
314 for i, line in enumerate(lines): 314 for i, line in enumerate(lines):
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/backends/backend_agg.py:259, in RendererAgg.get_text_width_height_descent(self, s, prop, ismath) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/backends/backend_agg.py:259, in RendererAgg.get_text_width_height_descent(self, s, prop, ismath)
257 texmanager = self.get_texmanager() 257 texmanager = self.get_texmanager()
258 fontsize = prop.get_size_in_points() 258 fontsize = prop.get_size_in_points()
--> 259 w, h, d = texmanager.get_text_width_height_descent( --> 259 w, h, d = texmanager.get_text_width_height_descent(
260 s, fontsize, renderer=self) 260 s, fontsize, renderer=self)
261 return w, h, d 261 return w, h, d
263 if ismath: 263 if ismath:
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/texmanager.py:335, in TexManager.get_text_width_height_descent(self, tex, fontsize, renderer) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/texmanager.py:335, in TexManager.get_text_width_height_descent(self, tex, fontsize, renderer)
333 if tex.strip() == '': 333 if tex.strip() == '':
334 return 0, 0, 0 334 return 0, 0, 0
--> 335 dvifile = self.make_dvi(tex, fontsize) --> 335 dvifile = self.make_dvi(tex, fontsize)
336 dpi_fraction = renderer.points_to_pixels(1.) if renderer else 1 336 dpi_fraction = renderer.points_to_pixels(1.) if renderer else 1
337 with dviread.Dvi(dvifile, 72 * dpi_fraction) as dvi: 337 with dviread.Dvi(dvifile, 72 * dpi_fraction) as dvi:
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/texmanager.py:271, in TexManager.make_dvi(self, tex, fontsize) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/texmanager.py:271, in TexManager.make_dvi(self, tex, fontsize)
262 # Generate the dvi in a temporary directory to avoid race 262 # Generate the dvi in a temporary directory to avoid race
263 # conditions e.g. if multiple processes try to process the same tex 263 # conditions e.g. if multiple processes try to process the same tex
264 # string at the same time. Having tmpdir be a subdirectory of the 264 # string at the same time. Having tmpdir be a subdirectory of the
(...) (...)
268 # the absolute path may contain characters (e.g. ~) that TeX does 268 # the absolute path may contain characters (e.g. ~) that TeX does
269 # not support.) 269 # not support.)
270 with TemporaryDirectory(dir=Path(dvifile).parent) as tmpdir: 270 with TemporaryDirectory(dir=Path(dvifile).parent) as tmpdir:
--> 271 self._run_checked_subprocess( --> 271 self._run_checked_subprocess(
272 ["latex", "-interaction=nonstopmode", "--halt-on-error", 272 ["latex", "-interaction=nonstopmode", "--halt-on-error",
273 f"../{texfile.name}"], tex, cwd=tmpdir) 273 f"../{texfile.name}"], tex, cwd=tmpdir)
274 (Path(tmpdir) / Path(dvifile).name).replace(dvifile) 274 (Path(tmpdir) / Path(dvifile).name).replace(dvifile)
275 return dvifile 275 return dvifile
File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/texmanager.py:241, in TexManager._run_checked_subprocess(self, command, tex, cwd) File /srv/main-spack-instance-2302/spack/var/spack/environments/ebrains-23-02/.spack-env/view/lib/python3.8/site-packages/matplotlib/texmanager.py:241, in TexManager._run_checked_subprocess(self, command, tex, cwd)
237 raise RuntimeError( 237 raise RuntimeError(
238 'Failed to process string with tex because {} could not be ' 238 'Failed to process string with tex because {} could not be '
239 'found'.format(command[0])) from exc 239 'found'.format(command[0])) from exc
240 except subprocess.CalledProcessError as exc: 240 except subprocess.CalledProcessError as exc:
--> 241 raise RuntimeError( --> 241 raise RuntimeError(
242 '{prog} was not able to process the following string:\n' 242 '{prog} was not able to process the following string:\n'
243 '{tex!r}\n\n' 243 '{tex!r}\n\n'
244 'Here is the full report generated by {prog}:\n' 244 'Here is the full report generated by {prog}:\n'
245 '{exc}\n\n'.format( 245 '{exc}\n\n'.format(
246 prog=command[0], 246 prog=command[0],
247 tex=tex.encode('unicode_escape'), 247 tex=tex.encode('unicode_escape'),
248 exc=exc.output.decode('utf-8'))) from exc 248 exc=exc.output.decode('utf-8'))) from exc
249 _log.debug(report) 249 _log.debug(report)
250 return report 250 return report
RuntimeError: latex was not able to process the following string: RuntimeError: latex was not able to process the following string:
b'lp' b'lp'
Here is the full report generated by latex: Here is the full report generated by latex:
This is pdfTeX, Version 3.14159265-2.6-1.40.20 (TeX Live 2019/Debian) (preloaded format=latex) This is pdfTeX, Version 3.14159265-2.6-1.40.20 (TeX Live 2019/Debian) (preloaded format=latex)
restricted \write18 enabled. restricted \write18 enabled.
entering extended mode entering extended mode
(../540885da22f8d5147f0088a11e43701a.tex (../540885da22f8d5147f0088a11e43701a.tex
LaTeX2e <2020-02-02> patch level 2 LaTeX2e <2020-02-02> patch level 2
L3 programming layer <2020-02-14> L3 programming layer <2020-02-14>
(/usr/share/texlive/texmf-dist/tex/latex/base/article.cls (/usr/share/texlive/texmf-dist/tex/latex/base/article.cls
Document Class: article 2019/12/20 v1.4l Standard LaTeX document class Document Class: article 2019/12/20 v1.4l Standard LaTeX document class
(/usr/share/texlive/texmf-dist/tex/latex/base/size10.clo)) (/usr/share/texlive/texmf-dist/tex/latex/base/size10.clo))
(/usr/share/texlive/texmf-dist/tex/latex/type1cm/type1cm.sty) (/usr/share/texlive/texmf-dist/tex/latex/type1cm/type1cm.sty)
! LaTeX Error: File `type1ec.sty' not found. ! LaTeX Error: File `type1ec.sty' not found.
Type X to quit or <RETURN> to proceed, Type X to quit or <RETURN> to proceed,
or enter new name. (Default extension: sty) or enter new name. (Default extension: sty)
Enter file name: Enter file name:
! Emergency stop. ! Emergency stop.
<read *> <read *>
l.6 \usepackage l.6 \usepackage
[utf8]{inputenc}^^M [utf8]{inputenc}^^M
No pages of output. No pages of output.
Transcript written on 540885da22f8d5147f0088a11e43701a.log. Transcript written on 540885da22f8d5147f0088a11e43701a.log.
%% Cell type:markdown id:e91c436e-db94-4cd7-a531-29c032efeeae tags: %% Cell type:markdown id:e91c436e-db94-4cd7-a531-29c032efeeae tags:
### 5.2 Resting state plots <a class="anchor" id="section_5_2"></a> ### 5.2 Resting state plots <a class="anchor" id="section_5_2"></a>
%% Cell type:markdown id:aeae56a4 tags: %% Cell type:markdown id:aeae56a4 tags:
**Fig 5. Resting state of the model with χ =1.9.** (A-C) Raster plot of spiking activity of 3% of the neurons in area V1 (A), V2 (B), and FEF (C). Blue: excitatory neurons, red: inhibitory neurons. (D-F) Spiking statistics across all 32 areas for the respective populations shown as area-averaged box plots. Crosses: medians, boxes: interquartile range (IQR), whiskers extend to the most extremeobservat ions within 1.5×IQR beyond the IQR. (D) Population-averaged firing rates. (E) Average pairwise correlation coefficients of spiking activity. (F) Irregularity measured by revised local variation LvR averaged across neurons. (G) Area-averaged firing rates, shown as raw binned spike histograms with 1ms bin width (gray) and convolved histograms, with aGaussian kernel (black) of optimal width. **Fig 5. Resting state of the model with χ =1.9.** (A-C) Raster plot of spiking activity of 3% of the neurons in area V1 (A), V2 (B), and FEF (C). Blue: excitatory neurons, red: inhibitory neurons. (D-F) Spiking statistics across all 32 areas for the respective populations shown as area-averaged box plots. Crosses: medians, boxes: interquartile range (IQR), whiskers extend to the most extremeobservat ions within 1.5×IQR beyond the IQR. (D) Population-averaged firing rates. (E) Average pairwise correlation coefficients of spiking activity. (F) Irregularity measured by revised local variation LvR averaged across neurons. (G) Area-averaged firing rates, shown as raw binned spike histograms with 1ms bin width (gray) and convolved histograms, with aGaussian kernel (black) of optimal width.
%% Cell type:code id:ae19bcc3 tags: %% Cell type:code id:ae19bcc3 tags:
``` python ``` python
from MAM2EBRAINS import plot_resting_state from MAM2EBRAINS import plot_resting_state
plot_resting_state(A) plot_resting_state(A)
``` ```
%% Output
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
Input In [29], in <module>
1 from MAM2EBRAINS import plot_resting_state
----> 2 plot_resting_state(A)
TypeError: plot_resting_state() takes 0 positional arguments but 1 was given
%% Cell type:markdown id:3ef52a7c tags: %% Cell type:markdown id:3ef52a7c tags:
### 5.2 Raster plot of spiking activity for single area <a class="anchor" id="section_5_2"></a> ### 5.2 Raster plot of spiking activity for single area <a class="anchor" id="section_5_2"></a>
Fig. 3 (A-C) Raster plot of spiking activity of 3% of the neurons in area V1 (A), V2 (B), and FEF (C). Blue: excitatory neurons, red: inhibitory neurons. (D-F) Spiking statistics across all 32 areas for the respective populations shown as area-averaged box plots. Crosses: medians, boxes: interquartile range (IQR), whiskers extend to the most extremeobservat ions within 1.5×IQR beyond the IQR. Fig. 3 (A-C) Raster plot of spiking activity of 3% of the neurons in area V1 (A), V2 (B), and FEF (C). Blue: excitatory neurons, red: inhibitory neurons. (D-F) Spiking statistics across all 32 areas for the respective populations shown as area-averaged box plots. Crosses: medians, boxes: interquartile range (IQR), whiskers extend to the most extremeobservat ions within 1.5×IQR beyond the IQR.
%% Cell type:code id:1da18fee tags: %% Cell type:code id:1da18fee tags:
``` python ``` python
M2E.plot_raster_plot(A) M2E.plot_raster_plot(A)
``` ```
%% Cell type:code id:7d5a5e32-bd12-4e91-a65d-91d279edc450 tags: %% Cell type:code id:7d5a5e32-bd12-4e91-a65d-91d279edc450 tags:
``` python ``` python
# load spike data # load spike data
spike_data = A.spike_data spike_data = A.spike_data
``` ```
%% Cell type:code id:7947fea3-ba4c-4b1d-94fc-16614e4e4a11 tags: %% Cell type:code id:7947fea3-ba4c-4b1d-94fc-16614e4e4a11 tags:
``` python ``` python
# plotting raster plot of spiking activity for single area # plotting raster plot of spiking activity for single area
from matplotlib import gridspec from matplotlib import gridspec
# axes = {} # axes = {}
# gs1 = gridspec.GridSpec(1, 3) # gs1 = gridspec.GridSpec(1, 3)
# # gs1.update(left=0.06, right=0.72, top=0.95, wspace=0.4, bottom=0.35) # # gs1.update(left=0.06, right=0.72, top=0.95, wspace=0.4, bottom=0.35)
# axes['A'] = plt.subplot(gs1[:1, :1], figsize=(16,9), gridspec_kw={'height_ratios': [1, 2]}) # axes['A'] = plt.subplot(gs1[:1, :1], figsize=(16,9), gridspec_kw={'height_ratios': [1, 2]})
# axes['B'] = plt.subplot(gs1[:1, 1:2]) # axes['B'] = plt.subplot(gs1[:1, 1:2])
# axes['C'] = plt.subplot(gs1[:1, 2:]) # axes['C'] = plt.subplot(gs1[:1, 2:])
f = plt.figure(figsize=(10,3)) f = plt.figure(figsize=(10,3))
sub = [131, 132, 133] sub = [131, 132, 133]
areas = ['V1', 'V2', 'FEF'] areas = ['V1', 'V2', 'FEF']
labels = ['A', 'B', 'C'] labels = ['A', 'B', 'C']
t_min = 0. t_min = 0.
t_max = 500. t_max = 500.
# t_min = 3000. # t_min = 3000.
# t_max = 3500. # t_max = 3500.
# icolor = myred # icolor = myred
# ecolor = myblue # ecolor = myblue
# frac_neurons = 0.03 # frac_neurons = 0.03
frac_neurons = 0.3 frac_neurons = 0.3
for i, area in enumerate(areas): for i, area in enumerate(areas):
# ax = axes[labels[i]] # ax = axes[labels[i]]
# ax = plt.subplot() # ax = plt.subplot()
ax = f.add_subplot(sub[i]) ax = f.add_subplot(sub[i])
if area in spike_data: if area in spike_data:
n_pops = len(spike_data[area]) n_pops = len(spike_data[area])
# Determine number of neurons that will be plotted for this area (for # Determine number of neurons that will be plotted for this area (for
# vertical offset) # vertical offset)
offset = 0 offset = 0
n_to_plot = {} n_to_plot = {}
for pop in M.structure[area]: for pop in M.structure[area]:
n_to_plot[pop] = int(M.N[area][pop] * frac_neurons) n_to_plot[pop] = int(M.N[area][pop] * frac_neurons)
offset = offset + n_to_plot[pop] offset = offset + n_to_plot[pop]
y_max = offset + 1 y_max = offset + 1
prev_pop = '' prev_pop = ''
yticks = [] yticks = []
yticklocs = [] yticklocs = []
for jj, pop in enumerate(M.structure[area]): for jj, pop in enumerate(M.structure[area]):
if pop[0:-1] != prev_pop: if pop[0:-1] != prev_pop:
prev_pop = pop[0:-1] prev_pop = pop[0:-1]
yticks.append('L' + population_labels[jj][0:-1]) yticks.append('L' + population_labels[jj][0:-1])
yticklocs.append(offset - 0.5 * n_to_plot[pop]) yticklocs.append(offset - 0.5 * n_to_plot[pop])
ind = np.where(np.logical_and( ind = np.where(np.logical_and(
spike_data[area][pop][:, 1] <= t_max, spike_data[area][pop][:, 1] >= t_min)) spike_data[area][pop][:, 1] <= t_max, spike_data[area][pop][:, 1] >= t_min))
pop_data = spike_data[area][pop][ind] pop_data = spike_data[area][pop][ind]
pop_neurons = np.unique(pop_data[:, 0]) pop_neurons = np.unique(pop_data[:, 0])
neurons_to_ = np.arange(np.min(spike_data[area][pop][:, 0]), np.min( neurons_to_ = np.arange(np.min(spike_data[area][pop][:, 0]), np.min(
spike_data[area][pop][:, 0]) + n_to_plot[pop], 1) spike_data[area][pop][:, 0]) + n_to_plot[pop], 1)
if pop.find('E') > (-1): if pop.find('E') > (-1):
pcolor = ecolor pcolor = ecolor
else: else:
pcolor = icolor pcolor = icolor
for kk in range(n_to_plot[pop]): for kk in range(n_to_plot[pop]):
spike_times = pop_data[pop_data[:, 0] == neurons_to_[kk], 1] spike_times = pop_data[pop_data[:, 0] == neurons_to_[kk], 1]
_ = ax.plot(spike_times, np.zeros(len(spike_times)) + _ = ax.plot(spike_times, np.zeros(len(spike_times)) +
offset - kk, '.', color=pcolor, markersize=1) offset - kk, '.', color=pcolor, markersize=1)
offset = offset - n_to_plot[pop] offset = offset - n_to_plot[pop]
y_min = offset y_min = offset
ax.set_title(areas[i]) ax.set_title(areas[i])
ax.set_xlim([t_min, t_max]) ax.set_xlim([t_min, t_max])
ax.set_ylim([y_min, y_max]) ax.set_ylim([y_min, y_max])
ax.set_yticklabels(yticks) ax.set_yticklabels(yticks)
ax.set_yticks(yticklocs) ax.set_yticks(yticklocs)
ax.set_xlabel('Time (s)', labelpad=-0.1) ax.set_xlabel('Time (s)', labelpad=-0.1)
ax.set_xticks([t_min, t_min + 250., t_max]) ax.set_xticks([t_min, t_min + 250., t_max])
ax.set_xticklabels([r'$3.$', r'$3.25$', r'$3.5$']) ax.set_xticklabels([r'$3.$', r'$3.25$', r'$3.5$'])
``` ```
%% Cell type:markdown id:019d805e tags: %% Cell type:markdown id:019d805e tags:
### 5.3 Population-averaged firing rates <a class="anchor" id="section_5_3"></a> ### 5.3 Population-averaged firing rates <a class="anchor" id="section_5_3"></a>
Fig 3. (D) Population-averaged firing rates Fig 3. (D) Population-averaged firing rates
%% Cell type:code id:b069bc59-44ae-450a-b0a5-b073951e3604 tags: %% Cell type:code id:b069bc59-44ae-450a-b0a5-b073951e3604 tags:
``` python ``` python
# load data # load data
# stationary firing rates # stationary firing rates
fn = os.path.join(data_path, label, 'Analysis', 'pop_rates.json') fn = os.path.join(data_path, label, 'Analysis', 'pop_rates.json')
with open(fn, 'r') as f: with open(fn, 'r') as f:
pop_rates = json.load(f) pop_rates = json.load(f)
``` ```
%% Cell type:code id:96baac8f-3a0e-4d69-92ac-1d7cb4aff0d8 tags: %% Cell type:code id:96baac8f-3a0e-4d69-92ac-1d7cb4aff0d8 tags:
``` python ``` python
def set_boxplot_props(d): def set_boxplot_props(d):
for i in range(len(d['boxes'])): for i in range(len(d['boxes'])):
if i % 2 == 0: if i % 2 == 0:
d['boxes'][i].set_facecolor(icolor) d['boxes'][i].set_facecolor(icolor)
d['boxes'][i].set_color(icolor) d['boxes'][i].set_color(icolor)
else: else:
d['boxes'][i].set_facecolor(ecolor) d['boxes'][i].set_facecolor(ecolor)
d['boxes'][i].set_color(ecolor) d['boxes'][i].set_color(ecolor)
plt.setp(d['whiskers'], color='k') plt.setp(d['whiskers'], color='k')
plt.setp(d['fliers'], color='k', markerfacecolor='k', marker='+') plt.setp(d['fliers'], color='k', markerfacecolor='k', marker='+')
plt.setp(d['medians'], color='none') plt.setp(d['medians'], color='none')
plt.setp(d['caps'], color='k') plt.setp(d['caps'], color='k')
plt.setp(d['means'], marker='x', color='k', plt.setp(d['means'], marker='x', color='k',
markerfacecolor='k', markeredgecolor='k', markersize=3.) markerfacecolor='k', markeredgecolor='k', markersize=3.)
# print("plotting Population rates") # print("plotting Population rates")
rates = np.zeros((len(M.area_list), 8)) rates = np.zeros((len(M.area_list), 8))
for i, area in enumerate(M.area_list): for i, area in enumerate(M.area_list):
for j, pop in enumerate(M.structure[area][::-1]): for j, pop in enumerate(M.structure[area][::-1]):
rate = pop_rates[area][pop][0] rate = pop_rates[area][pop][0]
if rate == 0.0: if rate == 0.0:
rate = 1e-5 rate = 1e-5
if area == 'TH' and j > 3: # To account for missing layer 4 in TH if area == 'TH' and j > 3: # To account for missing layer 4 in TH
rates[i][j + 2] = rate rates[i][j + 2] = rate
else: else:
rates[i][j] = rate rates[i][j] = rate
rates = np.transpose(rates) rates = np.transpose(rates)
masked_rates = np.ma.masked_where(rates < 1e-4, rates) masked_rates = np.ma.masked_where(rates < 1e-4, rates)
# ax = axes['D'] # ax = axes['D']
ax = plt.subplot() ax = plt.subplot()
d = plt.boxplot(np.transpose(rates), vert=False, d = plt.boxplot(np.transpose(rates), vert=False,
patch_artist=True, whis=1.5, showmeans=True) patch_artist=True, whis=1.5, showmeans=True)
set_boxplot_props(d) set_boxplot_props(d)
ax.plot(np.mean(rates, axis=1), np.arange( ax.plot(np.mean(rates, axis=1), np.arange(
1., len(M.structure['V1']) + 1., 1.), 'x', color='k', markersize=3) 1., len(M.structure['V1']) + 1., 1.), 'x', color='k', markersize=3)
ax.set_yticklabels(population_labels[::-1], size=8) ax.set_yticklabels(population_labels[::-1], size=8)
ax.set_yticks(np.arange(1., len(M.structure['V1']) + 1., 1.)) ax.set_yticks(np.arange(1., len(M.structure['V1']) + 1., 1.))
ax.set_ylim((0., len(M.structure['V1']) + .5)) ax.set_ylim((0., len(M.structure['V1']) + .5))
x_max = 100. x_max = 100.
ax.set_title("Population-averaged firing rates") ax.set_title("Population-averaged firing rates")
ax.set_xlim((-1., x_max)) ax.set_xlim((-1., x_max))
ax.set_xlabel(r'Rate (spikes/s)', labelpad=-0.1) ax.set_xlabel(r'Rate (spikes/s)', labelpad=-0.1)
ax.set_xticks([0., 50.]) ax.set_xticks([0., 50.])
``` ```
%% Cell type:markdown id:473d0882-8e45-4330-bfa2-2c7e1af0dac4 tags: %% Cell type:markdown id:473d0882-8e45-4330-bfa2-2c7e1af0dac4 tags:
### 5.4 Time-averaged population rates <a class="anchor" id="section_5_4"></a> ### 5.4 Time-averaged population rates <a class="anchor" id="section_5_4"></a>
%% Cell type:code id:721d1f03-df25-468d-8075-a807025a9c58 tags: %% Cell type:code id:721d1f03-df25-468d-8075-a807025a9c58 tags:
``` python ``` python
""" """
Plot overview over time-averaged population rates encoded in colors Plot overview over time-averaged population rates encoded in colors
with areas along x-axis and populations along y-axis. with areas along x-axis and populations along y-axis.
Parameters Parameters
---------- ----------
area_list : list, optional area_list : list, optional
Specifies with areas are plotted in which order. Specifies with areas are plotted in which order.
Default to None, leading to plotting of all areas ordered by architectural type. Default to None, leading to plotting of all areas ordered by architectural type.
output : {'pdf', 'png', 'eps'}, optional output : {'pdf', 'png', 'eps'}, optional
If given, the function stores the plot to a file of the given format. If given, the function stores the plot to a file of the given format.
""" """
A.show_rates() A.show_rates()
``` ```
%% Cell type:markdown id:06a595de tags: %% Cell type:markdown id:06a595de tags:
### 5.5 Average pairwise correlation coefficients of spiking activity <a class="anchor" id="section_5_5"></a> ### 5.5 Average pairwise correlation coefficients of spiking activity <a class="anchor" id="section_5_5"></a>
Fig 5. (E) Average pairwise correlation coefficients of spiking activity Fig 5. (E) Average pairwise correlation coefficients of spiking activity
%% Cell type:code id:a8e77836-4c37-4b78-b7c4-5e11bc67b4fa tags: %% Cell type:code id:a8e77836-4c37-4b78-b7c4-5e11bc67b4fa tags:
``` python ``` python
# load data # load data
# correlation coefficients # correlation coefficients
# fn = os.path.join(data_path, label, 'Analysis', 'corrcoeff.json') # fn = os.path.join(data_path, label, 'Analysis', 'corrcoeff.json')
fn = os.path.join(data_path, label, 'Analysis', 'synchrony.json') fn = os.path.join(data_path, label, 'Analysis', 'synchrony.json')
# synchrony.json # synchrony.json
with open(fn, 'r') as f: with open(fn, 'r') as f:
corrcoeff = json.load(f) corrcoeff = json.load(f)
``` ```
%% Cell type:code id:218367da-82ef-47b6-bf15-083ef3d43013 tags: %% Cell type:code id:218367da-82ef-47b6-bf15-083ef3d43013 tags:
``` python ``` python
# print("plotting Synchrony") # print("plotting Synchrony")
syn = np.zeros((len(M.area_list), 8)) syn = np.zeros((len(M.area_list), 8))
for i, area in enumerate(M.area_list): for i, area in enumerate(M.area_list):
for j, pop in enumerate(M.structure[area][::-1]): for j, pop in enumerate(M.structure[area][::-1]):
value = corrcoeff[area][pop] value = corrcoeff[area][pop]
if value == 0.0: if value == 0.0:
value = 1e-5 value = 1e-5
if area == 'TH' and j > 3: # To account for missing layer 4 in TH if area == 'TH' and j > 3: # To account for missing layer 4 in TH
syn[i][j + 2] = value syn[i][j + 2] = value
else: else:
syn[i][j] = value syn[i][j] = value
syn = np.transpose(syn) syn = np.transpose(syn)
masked_syn = np.ma.masked_where(syn < 1e-4, syn) masked_syn = np.ma.masked_where(syn < 1e-4, syn)
# ax = axes['E'] # ax = axes['E']
ax = plt.subplot() ax = plt.subplot()
d = ax.boxplot(np.transpose(syn), vert=False, d = ax.boxplot(np.transpose(syn), vert=False,
patch_artist=True, whis=1.5, showmeans=True) patch_artist=True, whis=1.5, showmeans=True)
set_boxplot_props(d) set_boxplot_props(d)
ax.plot(np.mean(syn, axis=1), np.arange( ax.plot(np.mean(syn, axis=1), np.arange(
1., len(M.structure['V1']) + 1., 1.), 'x', color='k', markersize=3) 1., len(M.structure['V1']) + 1., 1.), 'x', color='k', markersize=3)
ax.set_yticklabels(population_labels[::-1], size=8) ax.set_yticklabels(population_labels[::-1], size=8)
ax.set_yticks(np.arange(1., len(M.structure['V1']) + 1., 1.)) ax.set_yticks(np.arange(1., len(M.structure['V1']) + 1., 1.))
ax.set_ylim((0., len(M.structure['V1']) + .5)) ax.set_ylim((0., len(M.structure['V1']) + .5))
# ax.set_xticks(np.arange(0.0, 0.601, 0.2)) # ax.set_xticks(np.arange(0.0, 0.601, 0.2))
ax.set_xticks(np.arange(0.0, 10.0, 2.0)) ax.set_xticks(np.arange(0.0, 10.0, 2.0))
ax.set_xlabel('Correlation coefficient', labelpad=-0.1) ax.set_xlabel('Correlation coefficient', labelpad=-0.1)
``` ```
%% Cell type:markdown id:a3847e67 tags: %% Cell type:markdown id:a3847e67 tags:
### 5.6 Irregularity of spiking activity <a class="anchor" id="section_5_6"></a> ### 5.6 Irregularity of spiking activity <a class="anchor" id="section_5_6"></a>
Fig 5. (F) Irregularity measured by revised local variation LvR averaged across neurons Fig 5. (F) Irregularity measured by revised local variation LvR averaged across neurons
%% Cell type:code id:65377033-f3c0-4f90-be13-70594cfda292 tags: %% Cell type:code id:65377033-f3c0-4f90-be13-70594cfda292 tags:
``` python ``` python
# load data # load data
# local variance revised (LvR) # local variance revised (LvR)
fn = os.path.join(data_path, label, 'Analysis', 'pop_LvR.json') fn = os.path.join(data_path, label, 'Analysis', 'pop_LvR.json')
with open(fn, 'r') as f: with open(fn, 'r') as f:
pop_LvR = json.load(f) pop_LvR = json.load(f)
``` ```
%% Cell type:code id:d7480a9b tags: %% Cell type:code id:d7480a9b tags:
``` python ``` python
# print("plotting Irregularity") # print("plotting Irregularity")
LvR = np.zeros((len(M.area_list), 8)) LvR = np.zeros((len(M.area_list), 8))
for i, area in enumerate(M.area_list): for i, area in enumerate(M.area_list):
for j, pop in enumerate(M.structure[area][::-1]): for j, pop in enumerate(M.structure[area][::-1]):
value = pop_LvR[area][pop] value = pop_LvR[area][pop]
if value == 0.0: if value == 0.0:
value = 1e-5 value = 1e-5
if area == 'TH' and j > 3: # To account for missing layer 4 in TH if area == 'TH' and j > 3: # To account for missing layer 4 in TH
LvR[i][j + 2] = value LvR[i][j + 2] = value
else: else:
LvR[i][j] = value LvR[i][j] = value
LvR = np.transpose(LvR) LvR = np.transpose(LvR)
masked_LvR = np.ma.masked_where(LvR < 1e-4, LvR) masked_LvR = np.ma.masked_where(LvR < 1e-4, LvR)
# ax = axes['F'] # ax = axes['F']
ax = plt.subplot() ax = plt.subplot()
d = ax.boxplot(np.transpose(LvR), vert=False, d = ax.boxplot(np.transpose(LvR), vert=False,
patch_artist=True, whis=1.5, showmeans=True) patch_artist=True, whis=1.5, showmeans=True)
set_boxplot_props(d) set_boxplot_props(d)
ax.plot(np.mean(LvR, axis=1), np.arange( ax.plot(np.mean(LvR, axis=1), np.arange(
1., len(M.structure['V1']) + 1., 1.), 'x', color='k', markersize=3) 1., len(M.structure['V1']) + 1., 1.), 'x', color='k', markersize=3)
ax.set_yticklabels(population_labels[::-1], size=8) ax.set_yticklabels(population_labels[::-1], size=8)
ax.set_yticks(np.arange(1., len(M.structure['V1']) + 1., 1.)) ax.set_yticks(np.arange(1., len(M.structure['V1']) + 1., 1.))
ax.set_ylim((0., len(M.structure['V1']) + .5)) ax.set_ylim((0., len(M.structure['V1']) + .5))
x_max = 1.9 x_max = 1.9
ax.set_xlim((0., x_max)) ax.set_xlim((0., x_max))
ax.set_xlabel('Irregularity', labelpad=-0.1) ax.set_xlabel('Irregularity', labelpad=-0.1)
ax.set_xticks([0., 1., 2.]) ax.set_xticks([0., 1., 2.])
# axes['G'].spines['right'].set_color('none') # axes['G'].spines['right'].set_color('none')
# axes['G'].spines['left'].set_color('none') # axes['G'].spines['left'].set_color('none')
# axes['G'].spines['top'].set_color('none') # axes['G'].spines['top'].set_color('none')
# axes['G'].spines['bottom'].set_color('none') # axes['G'].spines['bottom'].set_color('none')
# axes['G'].yaxis.set_ticks_position("none") # axes['G'].yaxis.set_ticks_position("none")
# axes['G'].xaxis.set_ticks_position("none") # axes['G'].xaxis.set_ticks_position("none")
# axes['G'].set_xticks([]) # axes['G'].set_xticks([])
# axes['G'].set_yticks([]) # axes['G'].set_yticks([])
``` ```
%% Cell type:markdown id:90ae8f4c tags: %% Cell type:markdown id:90ae8f4c tags:
### 5.7 Time series of area-averaged firing rates <a class="anchor" id="section_5_7"></a> ### 5.7 Time series of area-averaged firing rates <a class="anchor" id="section_5_7"></a>
Area-averaged firing rates, shown as raw binned spike histograms with 1ms bin width (gray) and convolved histograms, with a Gaussian kernel (black) of optimal width. Area-averaged firing rates, shown as raw binned spike histograms with 1ms bin width (gray) and convolved histograms, with a Gaussian kernel (black) of optimal width.
%% Cell type:code id:0308d50a-1906-4860-9194-7f8664bd1f9d tags: %% Cell type:code id:0308d50a-1906-4860-9194-7f8664bd1f9d tags:
``` python ``` python
# load data # load data
# time series of firing rates # time series of firing rates
rate_time_series = {} rate_time_series = {}
for area in areas: for area in areas:
fn = os.path.join(data_path, label, fn = os.path.join(data_path, label,
'Analysis', 'Analysis',
'rate_time_series_full', 'rate_time_series_full',
'rate_time_series_full_{}.npy'.format(area)) 'rate_time_series_full_{}.npy'.format(area))
rate_time_series[area] = np.load(fn) rate_time_series[area] = np.load(fn)
# # time series of firing rates convolved with a kernel # # time series of firing rates convolved with a kernel
# rate_time_series_auto_kernel = {} # rate_time_series_auto_kernel = {}
# for area in areas: # for area in areas:
# fn = os.path.join(data_path, label, # fn = os.path.join(data_path, label,
# 'Analysis', # 'Analysis',
# 'rate_time_series_auto_kernel', # 'rate_time_series_auto_kernel',
# 'rate_time_series_auto_kernel_{}.npy'.format(area)) # 'rate_time_series_auto_kernel_{}.npy'.format(area))
# rate_time_series_auto_kernel[area] = np.load(fn) # rate_time_series_auto_kernel[area] = np.load(fn)
``` ```
%% Cell type:code id:4460d823-543a-482b-8ef1-a049e5837af4 tags: %% Cell type:code id:4460d823-543a-482b-8ef1-a049e5837af4 tags:
``` python ``` python
print("Plotting rate time series") print("Plotting rate time series")
pos = axes['G'].get_position() pos = axes['G'].get_position()
ax = [] ax = []
h = pos.y1 - pos.y0 h = pos.y1 - pos.y0
w = pos.x1 - pos.x0 w = pos.x1 - pos.x0
ax.append(pl.axes([pos.x0, pos.y0, w, 0.28 * h])) ax.append(pl.axes([pos.x0, pos.y0, w, 0.28 * h]))
ax.append(pl.axes([pos.x0, pos.y0 + 0.33 * h, w, 0.28 * h])) ax.append(pl.axes([pos.x0, pos.y0 + 0.33 * h, w, 0.28 * h]))
ax.append(pl.axes([pos.x0, pos.y0 + 0.67 * h, w, 0.28 * h])) ax.append(pl.axes([pos.x0, pos.y0 + 0.67 * h, w, 0.28 * h]))
colors = ['0.5', '0.3', '0.0'] colors = ['0.5', '0.3', '0.0']
t_min = 500. t_min = 500.
t_max = 10500. t_max = 10500.
time = np.arange(500., t_max) time = np.arange(500., t_max)
for i, area in enumerate(areas[::-1]): for i, area in enumerate(areas[::-1]):
ax[i].spines['right'].set_color('none') ax[i].spines['right'].set_color('none')
ax[i].spines['top'].set_color('none') ax[i].spines['top'].set_color('none')
ax[i].yaxis.set_ticks_position("left") ax[i].yaxis.set_ticks_position("left")
ax[i].xaxis.set_ticks_position("none") ax[i].xaxis.set_ticks_position("none")
binned_spikes = rate_time_series[area][np.where( binned_spikes = rate_time_series[area][np.where(
np.logical_and(time >= t_min, time < t_max))] np.logical_and(time >= t_min, time < t_max))]
ax[i].plot(time, binned_spikes, color=colors[0], label=area) ax[i].plot(time, binned_spikes, color=colors[0], label=area)
rate = rate_time_series_auto_kernel[area] rate = rate_time_series_auto_kernel[area]
ax[i].plot(time, rate, color=colors[2], label=area) ax[i].plot(time, rate, color=colors[2], label=area)
ax[i].set_xlim((500., t_max)) ax[i].set_xlim((500., t_max))
ax[i].text(0.8, 0.7, area, transform=ax[i].transAxes) ax[i].text(0.8, 0.7, area, transform=ax[i].transAxes)
if i > 0: if i > 0:
ax[i].spines['bottom'].set_color('none') ax[i].spines['bottom'].set_color('none')
ax[i].set_xticks([]) ax[i].set_xticks([])
ax[i].set_yticks([0., 30.]) ax[i].set_yticks([0., 30.])
else: else:
ax[i].set_xticks([1000., 5000., 10000.]) ax[i].set_xticks([1000., 5000., 10000.])
ax[i].set_xticklabels([r'$1.$', r'$5.$', r'$10.$']) ax[i].set_xticklabels([r'$1.$', r'$5.$', r'$10.$'])
ax[i].set_yticks([0., 5.]) ax[i].set_yticks([0., 5.])
if i == 1: if i == 1:
ax[i].set_ylabel(r'Rate (spikes/s)') ax[i].set_ylabel(r'Rate (spikes/s)')
ax[0].set_xlabel('Time (s)', labelpad=-0.05) ax[0].set_xlabel('Time (s)', labelpad=-0.05)
``` ```
%% Cell type:markdown id:ef74ca3e-98dc-49c9-a4a0-2c640e29b1d9 tags: %% Cell type:markdown id:ef74ca3e-98dc-49c9-a4a0-2c640e29b1d9 tags:
Go back to [Notebook structure](#toc) Go back to [Notebook structure](#toc)
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment