|
| 1 | +""" |
| 2 | +Example for usecases.rst |
| 3 | +""" |
| 4 | + |
| 5 | +from itertools import cycle |
| 6 | +import numpy as np |
| 7 | +from quantities import ms, mV, kHz |
| 8 | +import matplotlib.pyplot as plt |
| 9 | +from neo import Block, Segment, ChannelView, Group, SpikeTrain, AnalogSignal |
| 10 | + |
| 11 | +store_signals = False |
| 12 | + |
| 13 | +block = Block(name="probe data", tetrode_ids=["Tetrode #1", "Tetrode #2"]) |
| 14 | +block.segments = [ |
| 15 | + Segment(name="trial #1", index=0), |
| 16 | + Segment(name="trial #2", index=1), |
| 17 | + Segment(name="trial #3", index=2), |
| 18 | +] |
| 19 | + |
| 20 | +n_units = {"Tetrode #1": 2, "Tetrode #2": 5} |
| 21 | + |
| 22 | +# Create a group for each neuron, annotate each group with the tetrode from which it was recorded |
| 23 | +groups = [] |
| 24 | +counter = 0 |
| 25 | +for tetrode_id, n in n_units.items(): |
| 26 | + groups.extend([Group(name=f"neuron #{counter + i + 1}", tetrode_id=tetrode_id) for i in range(n)]) |
| 27 | + counter += n |
| 28 | +block.groups.extend(groups) |
| 29 | + |
| 30 | +iter_group = cycle(groups) |
| 31 | + |
| 32 | +# Create dummy data, one segment at a time |
| 33 | +for segment in block.segments: |
| 34 | + |
| 35 | + # create two 4-channel AnalogSignals with dummy data |
| 36 | + signals = { |
| 37 | + "Tetrode #1": AnalogSignal(np.random.rand(1000, 4) * mV, sampling_rate=10 * kHz, tetrode_id="Tetrode #1"), |
| 38 | + "Tetrode #2": AnalogSignal(np.random.rand(1000, 4) * mV, sampling_rate=10 * kHz, tetrode_id="Tetrode #2"), |
| 39 | + } |
| 40 | + if store_signals: |
| 41 | + segment.analogsignals.extend(signals.values()) |
| 42 | + |
| 43 | + # create spike trains with dummy data |
| 44 | + # we will pretend the spikes have been extracted from the dummy signal |
| 45 | + for tetrode_id in ("Tetrode #1", "Tetrode #2"): |
| 46 | + for i in range(n_units[tetrode_id]): |
| 47 | + spiketrain = SpikeTrain(np.random.uniform(0, 100, size=30) * ms, t_stop=100 * ms) |
| 48 | + # assign each spiketrain to the appropriate segment |
| 49 | + segment.spiketrains.append(spiketrain) |
| 50 | + # assign each spiketrain to a given neuron |
| 51 | + current_group = next(iter_group) |
| 52 | + current_group.add(spiketrain) |
| 53 | + if store_signals: |
| 54 | + # add to the group a reference to the signal from which the spikes were obtained |
| 55 | + # this does not give a 1:1 correspondance between spike trains and signals, |
| 56 | + # for that we could use additional groups (and have groups of groups) |
| 57 | + current_group.add(signals[tetrode_id]) |
| 58 | + |
| 59 | + |
| 60 | +# Now plot the data |
| 61 | + |
| 62 | +# .. by trial |
| 63 | +plt.figure() |
| 64 | +for seg in block.segments: |
| 65 | + print(f"Analyzing segment {seg.index}") |
| 66 | + stlist = [st - st.t_start for st in seg.spiketrains] |
| 67 | + plt.subplot(len(block.segments), 1, seg.index + 1) |
| 68 | + count, bins = np.histogram(stlist) |
| 69 | + plt.bar(bins[:-1], count, width=bins[1] - bins[0]) |
| 70 | + plt.title(f"PSTH in segment {seg.index}") |
| 71 | +plt.show() |
| 72 | + |
| 73 | +# ..by neuron |
| 74 | + |
| 75 | +plt.figure() |
| 76 | +for i, group in enumerate(block.groups): |
| 77 | + stlist = [st - st.t_start for st in group.spiketrains] |
| 78 | + plt.subplot(len(block.groups), 1, i + 1) |
| 79 | + count, bins = np.histogram(stlist) |
| 80 | + plt.bar(bins[:-1], count, width=bins[1] - bins[0]) |
| 81 | + plt.title(f"PSTH of unit {group.name}") |
| 82 | +plt.show() |
| 83 | + |
| 84 | +# ..by tetrode |
| 85 | + |
| 86 | +plt.figure() |
| 87 | +for i, tetrode_id in enumerate(block.annotations["tetrode_ids"]): |
| 88 | + stlist = [] |
| 89 | + for unit in block.filter(objects=Group, tetrode_id=tetrode_id): |
| 90 | + stlist.extend([st - st.t_start for st in unit.spiketrains]) |
| 91 | + plt.subplot(2, 1, i + 1) |
| 92 | + count, bins = np.histogram(stlist) |
| 93 | + plt.bar(bins[:-1], count, width=bins[1] - bins[0]) |
| 94 | + plt.title(f"PSTH blend of tetrode {tetrode_id}") |
| 95 | +plt.show() |
0 commit comments