NestIO compatible with Nest 3.2

[morales-gregorio]

Hi all,

@rgutzen @jasperalbers and I were looking into the NestIO today and we came up with this edits that make the current NestIO compatible with both Nest 2.x and 3.x versions.

We incorporated all the changes suggested by @essink into this PR.

Essentially we force the IO to ignore all the headers (which raises a warning that it did so) as suggested by @rgutzen, which I also wrote into PR #739

Additionally, we found that the current NestIO decides whether to read neo.AnalogSignals or neo.SpikeTrain based on the file name alone. This is problematic since Nest 3.2 (maybe others too) saves spikes to the .dat format by default, which the IO expected to be analogsignals. I fixed this issue by adding a new keyword target_object which then determines what kind of object will be read. This does not allow to read analogsignals and spikes into the same segment (two different calls are needed), but this seems a rather good practice in my opinion, since spiketrains and analogsignales are rarely combined into the same segment.

I tested the code to load a .dat file generated with Nest 3.2 and it works, I can add a proper test somewhere if you want.

Best,
Aitor

[morales-gregorio]

Looking into the capability to deal with lists it only seems to complicate things and it is not capable of dealing with multiple files of the same type, which is something I would expect. I would suggest getting rid of the capability to deal with lists altogether and stick to single files as most other IOs do. This then allows to remove code further and if someone really wants to have the analogsignals and spikes together it is just a matter of running this two times.

[Moritz-Alexander-Kern]

Hi all,
thanks for taking care of this.
We created a NEST-elephant tutorial and we ran into this issue with the nestio and NEST 3.3.

Since this was tested with NEST 3.2, maybe add this info to the docstring(?):

python-neo/neo/io/nestio.py

Line 4 in 50abfbc

Tested with NEST2.10.0

[JuliaSprenger]

Hi @morales-gregorio Do you have a small test file for this format that could be added to gin?

[morales-gregorio]

Just rebased the branch and some chaos has followed. My local history looks as one would expect, with this branch now splitting from the current master. I am not sure why Github interprets that all the previous commits are part of this PR :/

[morales-gregorio]

@JuliaSprenger I created a small NEST (3.3) simulation where I record both spikes and analogsignals. The results and simulation script are here: test_NEST_simulation_and_results.zip.

The current intended use is therefore something like this:

import neo
import quantities as pq
import matplotlib.pyplot as plt
# Load spikes
io = neo.NestIO(['sr1-142-0.dat', 'sr2-154-0.dat'], target_object='SpikeTrain')
bl = io.read_block([], t_start=0*pq.s, t_stop=20*pq.s)

# Plot spikes
plt.figure(figsize=(4, 3))
for i, st in enumerate(bl.segments[0].spiketrains):
    if st.annotations['source_file'] == 'sr1-142-0.dat':
        color = 'b'
    elif st.annotations['source_file'] == 'sr2-154-0.dat':
        color = 'r'
    plt.scatter(st.times, [i]*len(st.times), marker='|', c=color)
plt.xlabel('Time (ms)')
plt.ylabel('Neurons')
plt.savefig('nest_spikes.png', dpi=300)

# Load anasigs
io = neo.NestIO(['vm1-141-0.dat', 'vm2-153-0.dat'], target_object='AnalogSignal')
bl = io.read_block([], t_start=0*pq.s, t_stop=20*pq.s)

# Plot anasigs
plt.figure(figsize=(4, 3))
for i, ana in enumerate(bl.segments[0].analogsignals):
    if ana.annotations['source_file'] == 'vm1-141-0.dat':
        color = 'b'
    elif ana.annotations['source_file'] == 'vm2-153-0.dat':
        color = 'r'

    plt.plot(ana.times, 8*i + ana.magnitude, c=color, alpha=0.5)

plt.xlabel('Time (ms)')
plt.ylabel('Neurons')
plt.savefig('nest_anasigs.png', dpi=300)

Which creates these two plots in this case:

[heplesser]

Hello everyone!

After discussion with @mdenker and @apdavison yesterday, I provide here a little script that creates a network consisting of

• aeif_cond_alpha neurons (to show multimeter recording)
• iaf_psc_alpha_ps neurons (to show precise time recording)
• weight_recorder (to show weight recording)

Running the script will generate a number of spike recorder, voltmeter, multimeter and weight recorder output files with different settings. Note in particular the effect of "time_in_steps": If false (default), time is recorded in ms under heading "time_ms", if true time will be recorded in columns "time_step" and "time_offset", where the actual spike time is given by time_step * resolution - time_offset. Unfortunately, the resolution is not stored in the file, so I suggest to not support these files for now.

Do not rely on ordering of columns or rows. In parallel simulations, there will be one file per virtual process, enumerated by the last part of the file name.

NEST 3 files can be safely identified by the

# NEST version: 3.7.0-rc1
# RecordingBackendASCII version: 2

header. NEST 2 files do not have these headers. I would suggest to use the header for identification, also to be prepared for future versions of the ASCII backend.

nest.ResetKernel()
nest.overwrite_files = True

g = nest.Create('poisson_generator', params={'rate': 100000})
n = nest.Create('aeif_cond_alpha', 10)
p = nest.Create('iaf_psc_alpha_ps', 5)

rnt = nest.Create('spike_recorder', params={'record_to': 'ascii',
                                            'label': 'aeif_spikes_times'})
rns = nest.Create('spike_recorder', params={'record_to': 'ascii',
                                           'time_in_steps': True,
                                            'label': 'aeif_spikes_steps'})

rpt = nest.Create('spike_recorder', params={'record_to': 'ascii',
                                            'label': 'precise_spikes_times'})
rps = nest.Create('spike_recorder', params={'record_to': 'ascii',
                                           'time_in_steps': True,
                                            'label': 'precise_spikes_steps'})

vt = nest.Create('voltmeter', params={'record_to': 'ascii',
                                      'label': 'voltmeter_times'})
vs = nest.Create('voltmeter', params={'time_in_steps': True,
                                      'label': 'voltmeter_steps'})

mm1t = nest.Create('multimeter', params={'record_from': ('g_ex', 'g_in', 'V_m', 'w'),
                                         'record_to': 'ascii',
                                            'label': 'multimeter_1ms'})
mmrt = nest.Create('multimeter', params={'record_from': ('g_ex', 'g_in', 'V_m', 'w'),
                                         'record_to': 'ascii',
                                         'interval': nest.resolution,
                                            'label': 'multimeter_res'})

w = nest.Create('weight_recorder', params={'record_to': 'ascii',
                                            'label': 'weight_rec'})

nest.SetDefaults('stdp_synapse', {'weight_recorder': w})

nest.Connect(g, n)
nest.Connect(n, n, syn_spec={'synapse_model': 'stdp_synapse', 'weight': 100.})
nest.Connect(n, p, syn_spec={'weight': 1000.})
nest.Connect(n, rnt + rns)
nest.Connect(p, rpt + rps)
nest.Connect(vt + vs + mm1t + mmrt, n)

nest.Simulate(100)

[morales-gregorio]

Hi @heplesser ! Thanks for pushing this forward and coming up with good test examples.

I think we should not introduce a NEST dependency to the neo tests, I would instead suggest to create several (small) files for multiple NEST versions. The files should have descriptive names like "nest3.7_voltmeter_ascii.dat". What do you think?

[heplesser]

Hi @heplesser ! Thanks for pushing this forward and coming up with good test examples.

I think we should not introduce a NEST dependency to the neo tests, I would instead suggest to create several (small) files for multiple NEST versions. The files should have descriptive names like "nest3.7_voltmeter_ascii.dat". What do you think?

@morales-gregorio I agree, otherwise we could easily end up with circular dependencies and all-to-heavy test setups. Frozen example files for testing should work nicely. It should suffice to base the versioning on the ASCII Recording Backend version, which changes more slowly than the NEST version.

[morales-gregorio]

Hi @heplesser ! Thanks for pushing this forward and coming up with good test examples.
I think we should not introduce a NEST dependency to the neo tests, I would instead suggest to create several (small) files for multiple NEST versions. The files should have descriptive names like "nest3.7_voltmeter_ascii.dat". What do you think?

@morales-gregorio I agree, otherwise we could easily end up with circular dependencies and all-to-heavy test setups. Frozen example files for testing should work nicely. It should suffice to base the versioning on the ASCII Recording Backend version, which changes more slowly than the NEST version.

@heplesser could you create those files? 😇

[Moritz-Alexander-Kern]

Hey @morales-gregorio ,

the current test files for testing NEST-IO are located on NeuralEnsemble /
ephy_testing_data in nest here: https://gin.g-node.org/NeuralEnsemble/ephy_testing_data/src/master/nest

The contents are:

0gid-1time-1256-0.gdf
0gid-1time-2gex-1262-0.dat
0gid-1time-2gex-3Vm-1261-0.dat
0gid-1time-2Vm-1259-0.dat
0gid-1time-2Vm-3gex-4gin-1260-0.dat
0gid-1time-2Vm-3Iex-4Iin-1264-0.dat
0gid-1time_in_steps-1258-0.gdf
0gid-1time_in_steps-2Vm-1263-0.dat
0time-1255-0.gdf
0time_in_steps-1257-0.gdf
brunel-delta-nest_mod.py  # <--- script that generated the files
N1-0gid-1time-2Vm-1265-0.dat
N1-0time-1Vm-1266-0.dat
N1-0Vm-1267-0.dat

I have copied and formatted the above script into generate_data.py:

"""
Script to generate test data for the nestIO.

This script creates a network simulation using the NEST simulator (https://www.nest-simulator.org/).
It generates output files for spike recorders, voltmeters, multimeters and weight recorders with
different settings. Note the effect of "time_in_steps": If false (default), time is recorded in ms
under heading "time_ms"; if true, time will be recorded in columns "time_step" and "time_offset",
where the actual spike time is given by time_step * resolution - time_offset. The resolution is not
stored in the file. In parallel simulations, there will be one file per virtual process, enumerated
by the last part of the file name.

Usage:
    python generate_data.py


"""

import nest


def main():
    """
    Executes the main logic of the script.

    This function sets up a simulation using the NEST simulator (https://www.nest-simulator.org/)
    for simulating neural networks. It creates a network consisting of aeif_cond_alpha neurons
    (to show multimeter recording), iaf_psc_alpha_ps neurons (to show precise time recording),
    and weight_recorder (to show weight recording). The simulation is run for a specified duration.

    Note: Do not rely on ordering of columns or rows. In parallel simulations, there will be
    one file per virtual process, enumerated by the last part of the file name.

    Parameters:
        None

    Returns:
        None

    """
    nest.ResetKernel()
    nest.overwrite_files = True

    g = nest.Create("poisson_generator", params={"rate": 100000})
    n = nest.Create("aeif_cond_alpha", 10)
    p = nest.Create("iaf_psc_alpha_ps", 5)

    rnt = nest.Create(
        "spike_recorder",
        params={"record_to": "ascii", "label": "aeif_spikes_times"},
    )
    rns = nest.Create(
        "spike_recorder",
        params={
            "record_to": "ascii",
            "time_in_steps": True,
            "label": "aeif_spikes_steps",
        },
    )

    rpt = nest.Create(
        "spike_recorder",
        params={"record_to": "ascii", "label": "precise_spikes_times"},
    )
    rps = nest.Create(
        "spike_recorder",
        params={
            "record_to": "ascii",
            "time_in_steps": True,
            "label": "precise_spikes_steps",
        },
    )

    vt = nest.Create("voltmeter", params={"record_to": "ascii", "label": "voltmeter_times"})
    vs = nest.Create("voltmeter", params={"time_in_steps": True, "label": "voltmeter_steps"})

    mm1t = nest.Create(
        "multimeter",
        params={
            "record_from": ("g_ex", "g_in", "V_m", "w"),
            "record_to": "ascii",
            "label": "multimeter_1ms",
        },
    )
    mmrt = nest.Create(
        "multimeter",
        params={
            "record_from": ("g_ex", "g_in", "V_m", "w"),
            "record_to": "ascii",
            "interval": nest.resolution,
            "label": "multimeter_res",
        },
    )

    w = nest.Create("weight_recorder", params={"record_to": "ascii", "label": "weight_rec"})

    nest.SetDefaults("stdp_synapse", {"weight_recorder": w})

    nest.Connect(g, n)
    nest.Connect(n, n, syn_spec={"synapse_model": "stdp_synapse", "weight": 100.0})
    nest.Connect(n, p, syn_spec={"weight": 1000.0})
    nest.Connect(n, rnt + rns)
    nest.Connect(p, rpt + rps)
    nest.Connect(vt + vs + mm1t + mmrt, n)

    nest.Simulate(100)


if __name__ == "__main__":
    main()

which generates: (.zip attached)

aeif_spikes_steps-18-0.dat
aeif_spikes_times-17-0.dat
generate_data.py
multimeter_1ms-23-0.dat
multimeter_res-24-0.dat
precise_spikes_steps-20-0.dat
precise_spikes_times-19-0.dat
voltmeter_times-21-0.dat
weight_rec-25-0.dat

nestIO_test_data.zip

Refer to this readme here: https://gin.g-node.org/NeuralEnsemble/ephy_testing_data to add the files.

[apdavison]

There are 8 failing tests. I've looked into the first of these, and it seems to be that the case of id_column=None is not being correctly handled.

[mdenker]

Short update: A lot of info has been compiled on this PR, this should be sufficient to generate a decent NestIO for 3.x NEST releases. I am currently working on a refactor of this IO on https://github.com/INM-6/python-neo/tree/fork/NestIO_3.2 to implement header handling, integrate the new example data provided by Moritz and HEP, fix failing tests and adding new ones. I will create a PR against this PR once done.

[mdenker]

@heplesser: If one was to implement support for reading step-wise recording (times given as steps and offset) right away: Is there a possible sensible default value for the recording resolution? Or is it better to require this value explicitly from the user?

[samuelgarcia]

@mdenker : tell us when your branch will be ready enough to be merged.

[alejoe91]

@mdenker any update on this? Otherwise let's postpone it to 0.15.0 and move on with the release

[mdenker]

I'm afraid it's gonna be 0.15.0....