diff --git a/examples/300_pointneurons/build_network.py b/examples/300_pointneurons/build_network.py new file mode 100644 index 0000000..ee8a562 --- /dev/null +++ b/examples/300_pointneurons/build_network.py @@ -0,0 +1,157 @@ +import os +import numpy as np + +from bmtk.builder import NetworkBuilder +from bmtk.builder.bionet import SWCReader +from bmtk.utils.io.spike_trains import PoissonSpikesGenerator +from bmtk.builder.aux.node_params import positions_columinar, xiter_random + +build_recurrent_edges = True + +print('Building internal network') +# List of non-virtual cell models +cell_models = [ + { + 'model_name': 'Scnn1a', 'ei': 'e', + 'model_template': 'nest:iaf_psc_alpha', + 'dynamics_params': '472363762_point.json' + }, + { + 'model_name': 'Rorb', 'ei': 'e', + 'model_template': 'nest:iaf_psc_alpha', + 'dynamics_params': '473863510_point.json' + }, + { + 'model_name': 'Nr5a1', 'ei': 'e', + 'model_template': 'nest:iaf_psc_alpha', + 'dynamics_params': '473863035_point.json' + }, + { + 'model_name': 'PV1', 'ei': 'i', + 'model_template': 'nest:iaf_psc_alpha', + 'dynamics_params': '472912177_point.json' + }, + { + 'model_name': 'PV2', 'ei': 'i', + 'model_template': 'nest:iaf_psc_alpha', + 'dynamics_params': '473862421_point.json' + } +] + +''' +morphologies = {p['model_name']: SWCReader(os.path.join('../shared_components/morphologies', + '{}.swc'.format(p['morphology']))) + for p in cell_models} +def build_edges(src, trg, sections=['basal', 'apical'], dist_range=[50.0, 150.0]): + """Function used to randomly assign a synaptic location based on the section (soma, basal, apical) and an + arc-length dist_range from the soma. This function should be passed into the network and called during the build + process. + + :param src: source cell (dict) + :param trg: target cell (dict) + :param sections: list of target cell sections to synapse onto + :param dist_range: range (distance from soma center) to place + :return: + """ + # Get morphology and soma center for the target cell + swc_reader = morphologies[trg['model_name']] + target_coords = [trg['x'], trg['y'], trg['z']] + + sec_ids, sec_xs = swc_reader.choose_sections(sections, dist_range) # randomly choose sec_ids + coords = swc_reader.get_coord(sec_ids, sec_xs, soma_center=target_coords) # get coords of sec_ids + dist = swc_reader.get_dist(sec_ids) + swctype = swc_reader.get_type(sec_ids) + return sec_ids, sec_xs, coords[0][0], coords[0][1], coords[0][2], dist[0], swctype[0] +''' + +# Build a network of 300 biophysical cells to simulate +internal = NetworkBuilder("internal") +for i, model_props in enumerate(cell_models): + n_cells = 80 if model_props['ei'] == 'e' else 30 # 80% excitatory, 20% inhib + + # Randomly get positions uniformly distributed in a column + positions = positions_columinar(N=n_cells, center=[0, 10.0, 0], max_radius=50.0, height=200.0) + + internal.add_nodes(N=n_cells, + x=positions[:, 0], y=positions[:, 1], z=positions[:, 2], + rotation_angle_yaxis=xiter_random(N=n_cells, min_x=0.0, max_x=2 * np.pi), # randomly rotate y axis + model_type='point_process', + **model_props) + + + +if build_recurrent_edges: + def n_connections(src, trg, prob=0.5, min_syns=2, max_syns=7): + return 0 if np.random.uniform() > prob else np.random.randint(min_syns, max_syns) + + # exc --> exc connections + cm = internal.add_edges(source={'ei': 'e'}, target={'ei': 'e'}, + connection_rule=lambda s, t: np.random.binomial(1, 0.2), + dynamics_params='ExcToExc.json', + model_template='static_synapse', + # syn_weight=2.5, + delay=2.0) + cm.add_properties('syn_weight', rule=2.5, dtypes=np.float) + + # exc --> inh connections + cm = internal.add_edges(source={'ei': 'e'}, target={'ei': 'i'}, + connection_rule=lambda s, t: np.random.binomial(1, 0.5), + dynamics_params='ExcToInh.json', + model_template='static_synapse', + # syn_weight=7.0, + delay=2.0) + cm.add_properties('syn_weight', rule=7.0, dtypes=np.float) + + # inh --> exc connections + cm = internal.add_edges(source={'ei': 'i'}, target={'ei': 'e'}, + connection_rule=lambda s, t: np.random.binomial(1, 0.5), + #connection_rule=lambda *_: np.random.randint(0, 4), + dynamics_params='InhToExc.json', + model_template='static_synapse', + # syn_weight=-7.5, + delay=2.0) + cm.add_properties('syn_weight', rule=-7.5, dtypes=np.float) + + # inh --> inh connections + cm = internal.add_edges(source={'ei': 'i'}, target={'ei': 'i'}, + connection_rule=lambda s, t: np.random.binomial(1, 0.5), + dynamics_params='InhToInh.json', + model_template='static_synapse', + # syn_weight=-3.0, + delay=2.0) + cm.add_properties('syn_weight', rule=-3.0, dtypes=np.float) + + +internal.build() + +print('Saving internal') +internal.save(output_dir='network') + + +print('Building external connections') +external = NetworkBuilder("external") +external.add_nodes(N=100, model_type='virtual', ei='e') +cm = external.add_edges(target=internal.nodes(ei='e'), source=external.nodes(), + connection_rule=lambda *_: np.random.binomial(1, .7), + dynamics_params='ExcToExc.json', + model_template='static_synapse', + # syn_weight=50.0, + delay=2.0) +cm.add_properties('syn_weight', rule=50.0, dtypes=np.float) + + +cm = external.add_edges(target=internal.nodes(ei='i'), source=external.nodes(), + connection_rule=lambda *_: np.random.binomial(1, .7), + dynamics_params='ExcToInh.json', + model_template='static_synapse', + # syn_weight=75.0, + delay=2.0) +cm.add_properties('syn_weight', rule=65.0, dtypes=np.float) + + +external.build() + +print('Saving external') +external.save(output_dir='network') + + diff --git a/examples/300_pointneurons/circuit_config.json b/examples/300_pointneurons/circuit_config.json new file mode 100755 index 0000000..79bf387 --- /dev/null +++ b/examples/300_pointneurons/circuit_config.json @@ -0,0 +1,37 @@ +{ + "manifest": { + "$NETWORK_DIR": "./network", + "$COMPONENT_DIR": "../shared_components/nest_models" + }, + + "components": { + "point_neuron_models_dir": "$COMPONENT_DIR/cell_models", + "synaptic_models_dir": "$COMPONENT_DIR/synaptic_models" + }, + + "networks": { + "nodes": [ + { + "nodes_file": "$NETWORK_DIR/internal_nodes.h5", + "node_types_file": "$NETWORK_DIR/internal_node_types.csv" + }, + { + "nodes_file": "$NETWORK_DIR/external_nodes.h5", + "node_types_file": "$NETWORK_DIR/external_node_types.csv" + } + ], + + "edges": [ + { + "edges_file": "$NETWORK_DIR/internal_internal_edges.h5", + "edge_types_file": "$NETWORK_DIR/internal_internal_edge_types.csv", + "enabled": true + }, + { + "edges_file": "$NETWORK_DIR/external_internal_edges.h5", + "edge_types_file": "$NETWORK_DIR/external_internal_edge_types.csv" + } + ] + + } +} \ No newline at end of file diff --git a/examples/300_pointneurons/config.json b/examples/300_pointneurons/config.json new file mode 100755 index 0000000..890caa3 --- /dev/null +++ b/examples/300_pointneurons/config.json @@ -0,0 +1,4 @@ +{ + "network": "./circuit_config.json", + "simulation": "./simulation_config.json" +} diff --git a/examples/300_pointneurons/generate_spikes.py b/examples/300_pointneurons/generate_spikes.py new file mode 100644 index 0000000..aaf1739 --- /dev/null +++ b/examples/300_pointneurons/generate_spikes.py @@ -0,0 +1,7 @@ +import os +from bmtk.utils.io.spike_trains import PoissonSpikesGenerator + +if not os.path.exists('inputs'): + os.mkdir('inputs') +psg = PoissonSpikesGenerator(range(100), 14.0, tstop=3000.0, precision=3) +psg.to_hdf5('inputs/external_spike_trains.h5') \ No newline at end of file diff --git a/examples/300_pointneurons/inputs/external_spike_trains.h5 b/examples/300_pointneurons/inputs/external_spike_trains.h5 new file mode 100644 index 0000000..106a377 Binary files /dev/null and b/examples/300_pointneurons/inputs/external_spike_trains.h5 differ diff --git a/examples/300_pointneurons/network/external_internal_edge_types.csv b/examples/300_pointneurons/network/external_internal_edge_types.csv new file mode 100644 index 0000000..5e2bf83 --- /dev/null +++ b/examples/300_pointneurons/network/external_internal_edge_types.csv @@ -0,0 +1,3 @@ +edge_type_id target_query source_query delay dynamics_params model_template +100 ei=='e' * 2.0 ExcToExc.json static_synapse +101 ei=='i' * 2.0 ExcToInh.json static_synapse diff --git a/examples/300_pointneurons/network/external_internal_edges.h5 b/examples/300_pointneurons/network/external_internal_edges.h5 new file mode 100644 index 0000000..344d118 Binary files /dev/null and b/examples/300_pointneurons/network/external_internal_edges.h5 differ diff --git a/examples/300_pointneurons/network/external_node_types.csv b/examples/300_pointneurons/network/external_node_types.csv new file mode 100644 index 0000000..07b4593 --- /dev/null +++ b/examples/300_pointneurons/network/external_node_types.csv @@ -0,0 +1,2 @@ +node_type_id model_type ei +100 virtual e diff --git a/examples/300_pointneurons/network/external_nodes.h5 b/examples/300_pointneurons/network/external_nodes.h5 new file mode 100644 index 0000000..d3375c6 Binary files /dev/null and b/examples/300_pointneurons/network/external_nodes.h5 differ diff --git a/examples/300_pointneurons/network/internal_internal_edge_types.csv b/examples/300_pointneurons/network/internal_internal_edge_types.csv new file mode 100644 index 0000000..02a5c35 --- /dev/null +++ b/examples/300_pointneurons/network/internal_internal_edge_types.csv @@ -0,0 +1,5 @@ +edge_type_id target_query source_query delay dynamics_params model_template +100 ei=='e' ei=='e' 2.0 ExcToExc.json static_synapse +101 ei=='i' ei=='e' 2.0 ExcToInh.json static_synapse +102 ei=='e' ei=='i' 2.0 InhToExc.json static_synapse +103 ei=='i' ei=='i' 2.0 InhToInh.json static_synapse diff --git a/examples/300_pointneurons/network/internal_internal_edges.h5 b/examples/300_pointneurons/network/internal_internal_edges.h5 new file mode 100644 index 0000000..3c6390b Binary files /dev/null and b/examples/300_pointneurons/network/internal_internal_edges.h5 differ diff --git a/examples/300_pointneurons/network/internal_node_types.csv b/examples/300_pointneurons/network/internal_node_types.csv new file mode 100644 index 0000000..2209b58 --- /dev/null +++ b/examples/300_pointneurons/network/internal_node_types.csv @@ -0,0 +1,6 @@ +node_type_id ei model_template model_type dynamics_params model_name +104 i nest:iaf_psc_alpha point_process 473862421_point.json PV2 +100 e nest:iaf_psc_alpha point_process 472363762_point.json Scnn1a +101 e nest:iaf_psc_alpha point_process 473863510_point.json Rorb +102 e nest:iaf_psc_alpha point_process 473863035_point.json Nr5a1 +103 i nest:iaf_psc_alpha point_process 472912177_point.json PV1 diff --git a/examples/300_pointneurons/network/internal_nodes.h5 b/examples/300_pointneurons/network/internal_nodes.h5 new file mode 100644 index 0000000..d039746 Binary files /dev/null and b/examples/300_pointneurons/network/internal_nodes.h5 differ diff --git a/examples/300_pointneurons/node_sets.json b/examples/300_pointneurons/node_sets.json new file mode 100644 index 0000000..fc7f686 --- /dev/null +++ b/examples/300_pointneurons/node_sets.json @@ -0,0 +1,7 @@ +{ + "external": { + "population": "external" + }, + + "recorded_cells": [0, 80, 160, 240, 270] +} \ No newline at end of file diff --git a/examples/300_pointneurons/output/membrane_potential.h5 b/examples/300_pointneurons/output/membrane_potential.h5 new file mode 100644 index 0000000..8d26757 Binary files /dev/null and b/examples/300_pointneurons/output/membrane_potential.h5 differ diff --git a/examples/300_pointneurons/output/spikes.h5 b/examples/300_pointneurons/output/spikes.h5 new file mode 100644 index 0000000..06da4dc Binary files /dev/null and b/examples/300_pointneurons/output/spikes.h5 differ diff --git a/examples/300_pointneurons/plot_potential.py b/examples/300_pointneurons/plot_potential.py new file mode 100644 index 0000000..9675650 --- /dev/null +++ b/examples/300_pointneurons/plot_potential.py @@ -0,0 +1,76 @@ +""" +A script for plotting cell variable (eg voltage, calcium) traces from a SONATA hdf5 file. + +To plot a single ouput: + $ python plot_potential.py path/to/cell_var.h5 + +If you want to perform a side-by-side comparison of multiple output files: + $ python plot_potential.py cell_var_0.h5 cell_var_1.h5 ... cell_var_n.h5 +cell_var_0.h5 is used as the base-line for plotting number of gid/subplots and time-traces. + +""" + +import os +import matplotlib.pyplot as plt +from optparse import OptionParser +from bmtk.utils.cell_vars import CellVarsFile + + +def plot_vars(file_names, cell_var='v', gid_list=[]): + """Plots variable traces for a SONATA h5 file. If multiple spike files are specified will do a side-by-side + comparsion for each gid. + + :param file_names: list of cell_var file names + :param cell_var: cell variable to plot trace + :param gid_list: used to set what gid/subplots to show (if empty list just plot all possible gids) + """ + # convert to list if single spike file passed in + file_names = [file_names] if not isinstance(file_names, (tuple, list)) else file_names + assert(len(file_names) > 0) + + # Use bmtk to parse the cell-var files + cell_var_files = [CellVarsFile(fn) for fn in file_names] + + # get first spike file and properties + cvf_base = cell_var_files[0] + xlim = [cvf_base.time_trace[0], cvf_base.time_trace[-1]] # Use the same x-axis across all subplots + gid_list = cvf_base.gids if not gid_list else gid_list # if gid_list is None just get all gids in first file + n_cells = len(cvf_base.gids) + + fig, ax = plt.subplots(n_cells, 1, figsize=(10, 10)) + for subplot, gid in enumerate(gid_list): + for i, cvf in enumerate(cell_var_files): + # plot all traces + ax[subplot].plot(cvf.time_trace, cvf.data(gid, cell_var), label=file_names[i]) + + ax[subplot].yaxis.set_label_position("right") + ax[subplot].set_ylabel('gid {}'.format(gid), fontsize='xx-small') + ax[subplot].set_xlim(xlim) + if subplot + 1 < n_cells: + # remove x-axis labels on all but the last plot + ax[subplot].set_xticklabels([]) + else: + # Use the last plot to get the legend + handles, labels = ax[subplot].get_legend_handles_labels() + fig.legend(handles, labels, loc='upper right') + + plt.show() + + +if __name__ == '__main__': + parser = OptionParser(usage="Usage: python %prog [options] .h5 [.h5 ...]") + parser.add_option('--variable', type=str, dest='cell_var', default='V_m', help='Cell variable to compare (v, cai, etc)') + parser.add_option('--gids', type='string', dest='gids', default=[], + action='callback', + callback=lambda option, opt, value, parser: setattr(parser.values, option.dest, [int(v) for v in value.split(',')]), + help='comma seperated list of gids to plot') + options, args = parser.parse_args() + + if len(args) == 0: + # If no file is specified see if there is a output/membrane_potential.h5 file to plot membrane voltage + if not os.path.exists('output/membrane_potential.h5'): + raise Exception('Please specifiy hdf5 file to read in arguments. Exiting!') + else: + plot_vars('output/membrane_potential.h5', cell_var=options.cell_var) + else: + plot_vars(file_names=args, cell_var=options.cell_var, gid_list=options.gids) diff --git a/examples/300_pointneurons/plot_spikes.py b/examples/300_pointneurons/plot_spikes.py new file mode 100644 index 0000000..999c3a9 --- /dev/null +++ b/examples/300_pointneurons/plot_spikes.py @@ -0,0 +1,3 @@ +from bmtk.analyzer.visualization.spikes import plot_spikes + +plot_spikes('network/internal_nodes.h5', 'network/internal_node_types.csv', 'output/spikes.h5', group_key='model_name') diff --git a/examples/300_pointneurons/run_pointnet.py b/examples/300_pointneurons/run_pointnet.py new file mode 100755 index 0000000..49a015f --- /dev/null +++ b/examples/300_pointneurons/run_pointnet.py @@ -0,0 +1,18 @@ +import os, sys +from bmtk.simulator import pointnet +from bmtk.analyzer.visualization.spikes import plot_spikes + + +def main(config_file): + configure = pointnet.Config.from_json(config_file) + configure.build_env() + + graph = pointnet.PointNetwork.from_config(configure) + sim = pointnet.PointSimulator.from_config(configure, graph) + sim.run() + + # plot_spikes('network/internal_nodes.h5', 'network/internal_node_types.csv', 'output/spikes.h5', group_key='model_name') + + +if __name__ == '__main__': + main('config.json') diff --git a/examples/300_pointneurons/simulation_config.json b/examples/300_pointneurons/simulation_config.json new file mode 100755 index 0000000..0f7e15f --- /dev/null +++ b/examples/300_pointneurons/simulation_config.json @@ -0,0 +1,49 @@ +{ + "manifest": { + "$BASE_DIR": ".", + "$OUTPUT_DIR": "$BASE_DIR/output", + "$INPUT_DIR": "$BASE_DIR/inputs" + }, + + "run": { + "tstop": 1500.0, + "dt": 0.001, + "nsteps_block": 5000 + }, + + "target_simulator":"NEST", + + "network": "$BASE_DIR/circuit_config.json", + + "conditions": { + "celsius": 34.0, + "v_init": -80 + }, + + "node_sets_file": "node_sets.json", + + "inputs": { + "external_spike_trains": { + "input_type": "spikes", + "module": "h5", + "input_file": "$INPUT_DIR/external_spike_trains.h5", + "node_set": "external" + } + }, + + "output": { + "log_file": "log.txt", + "output_dir": "$OUTPUT_DIR", + "spikes_file": "spikes.h5", + "spikes_sort_order": "time" + }, + + "reports": { + "membrane_potential": { + "cells": "recorded_cells", + "variable_name": "V_m", + "module": "membrane_report", + "sections": "multimeter_report" + } + } +} diff --git a/examples/shared_components/nest_models/cell_models/472363762_point.json b/examples/shared_components/nest_models/cell_models/472363762_point.json new file mode 100644 index 0000000..e6154b1 --- /dev/null +++ b/examples/shared_components/nest_models/cell_models/472363762_point.json @@ -0,0 +1,9 @@ +{ + "I_e": 0.0, + "tau_m": 44.9, + "C_m": 239.0, + "t_ref": 3.0, + "E_L": -78.0, + "V_th": -43.0, + "V_reset": -55.0 +} diff --git a/examples/shared_components/nest_models/cell_models/472912177_point.json b/examples/shared_components/nest_models/cell_models/472912177_point.json new file mode 100644 index 0000000..30b9822 --- /dev/null +++ b/examples/shared_components/nest_models/cell_models/472912177_point.json @@ -0,0 +1,9 @@ +{ + "I_e": 0.0, + "tau_m": 22.2, + "C_m": 180.0, + "t_ref": 3.0, + "E_L": -82.0, + "V_th": -35.0, + "V_reset": -50.0 +} diff --git a/examples/shared_components/nest_models/cell_models/473862421_point.json b/examples/shared_components/nest_models/cell_models/473862421_point.json new file mode 100644 index 0000000..6d7e76a --- /dev/null +++ b/examples/shared_components/nest_models/cell_models/473862421_point.json @@ -0,0 +1,9 @@ +{ + "I_e": 0.0, + "tau_m": 12.5, + "C_m": 78.0, + "t_ref": 3.0, + "E_L": -73.0, + "V_th": -37.0, + "V_reset": -55.0 +} diff --git a/examples/shared_components/nest_models/cell_models/473863035_point.json b/examples/shared_components/nest_models/cell_models/473863035_point.json new file mode 100644 index 0000000..db8e5e4 --- /dev/null +++ b/examples/shared_components/nest_models/cell_models/473863035_point.json @@ -0,0 +1,9 @@ +{ + "I_e": 0.0, + "tau_m": 22.1, + "C_m": 117.0, + "t_ref": 3.0, + "E_L": -78.0, + "V_th": -47.0, + "V_reset": -50.0 +} diff --git a/examples/shared_components/nest_models/cell_models/473863510_point.json b/examples/shared_components/nest_models/cell_models/473863510_point.json new file mode 100644 index 0000000..348c569 --- /dev/null +++ b/examples/shared_components/nest_models/cell_models/473863510_point.json @@ -0,0 +1,9 @@ +{ + "I_e": 0.0, + "tau_m": 11.5, + "C_m": 53.0, + "t_ref": 3.0, + "E_L": -72.0, + "V_th": -25.0, + "V_reset": -50.0 +} diff --git a/examples/shared_components/nest_models/synaptic_models/ExcToExc.json b/examples/shared_components/nest_models/synaptic_models/ExcToExc.json new file mode 100644 index 0000000..2c63c08 --- /dev/null +++ b/examples/shared_components/nest_models/synaptic_models/ExcToExc.json @@ -0,0 +1,2 @@ +{ +} diff --git a/examples/shared_components/nest_models/synaptic_models/ExcToInh.json b/examples/shared_components/nest_models/synaptic_models/ExcToInh.json new file mode 100644 index 0000000..2c63c08 --- /dev/null +++ b/examples/shared_components/nest_models/synaptic_models/ExcToInh.json @@ -0,0 +1,2 @@ +{ +} diff --git a/examples/shared_components/nest_models/synaptic_models/InhToExc.json b/examples/shared_components/nest_models/synaptic_models/InhToExc.json new file mode 100644 index 0000000..bfd870e --- /dev/null +++ b/examples/shared_components/nest_models/synaptic_models/InhToExc.json @@ -0,0 +1,3 @@ +{ +} + diff --git a/examples/shared_components/nest_models/synaptic_models/InhToInh.json b/examples/shared_components/nest_models/synaptic_models/InhToInh.json new file mode 100644 index 0000000..bfd870e --- /dev/null +++ b/examples/shared_components/nest_models/synaptic_models/InhToInh.json @@ -0,0 +1,3 @@ +{ +} +