more examples touchups

Author: zm711

examples/plot_igorio.py

@@ -16,6 +16,8 @@
 # Downloaded from Human Brain Project Collaboratory
 # Digital Reconstruction of Neocortical Microcircuitry (nmc-portal)
 # http://microcircuits.epfl.ch/#/animal/8ecde7d1-b2d2-11e4-b949-6003088da632
+
+
 datafile_url = "https://microcircuits.epfl.ch/data/released_data/B95.zip"
 filename_zip = "B95.zip"
 filename = "grouped_ephys/B95/B95_Ch0_IDRest_107.ibw"

examples/plot_imageseq.py

@@ -19,6 +19,7 @@
 # Now we need to generate some data
 # We will just make a nice box and then we can attach this 
 # ImageSequence to a variety of ROIs
+# our ImageSequence will be 50 frames of 100x100 pixel images
 
 l = []
 for frame in range(50):
@@ -28,6 +29,10 @@
         for x in range(100):
             l[frame][y].append(random.randint(0, 50))
 
+#####################################################################
+# we then make our image sequence and pull out our results from the
+# image_seq
+
 image_seq = ImageSequence(l, sampling_rate=500 * pq.Hz, spatial_scale="m", units="V")
 
 result = image_seq.signal_from_region(
@@ -44,5 +49,5 @@
     plt.plot(result[i].times, result[i])
     plt.xlabel("seconde")
     plt.ylabel("valeur")
-
+plt.tight_layout()
 plt.show()

examples/plot_multi_tetrode_example.py

@@ -2,13 +2,21 @@
 Analyzing and Plotting Data with Neo Structures
 ===============================================
 """
+######################################################
+# First we import some packages. Since we are making simulated
+# data we will import quite a few neo features as well as use
+# quantities to provide units
 
 from itertools import cycle
 import numpy as np
 from quantities import ms, mV, kHz
 import matplotlib.pyplot as plt
 from neo import Block, Segment, ChannelView, Group, SpikeTrain, AnalogSignal
 
+##########################################################################
+# For Neo we start with a block of data that will contain segments of data
+# so we will create a block of probe data that has a couple tetrodes
+# Then we will load in 3 segments (for examples trials of a stimulus)
 store_signals = False
 
 block = Block(name="probe data", tetrode_ids=["Tetrode #1", "Tetrode #2"])
@@ -18,8 +26,13 @@
     Segment(name="trial #3", index=2),
 ]
 
+# we will decide how many units each tetrode has found. If only science was this easy
 n_units = {"Tetrode #1": 2, "Tetrode #2": 5}
 
+##################################################################################
+# Neo can also have groups. Groups are structures within a block that can cross segments
+# for example we could group a neuron across trials or across probes. 
+
 # Create a group for each neuron, annotate each group with the tetrode from which it was recorded
 groups = []
 counter = 0
@@ -30,19 +43,28 @@
 
 iter_group = cycle(groups)
 
+##########################################################################################
+# Segments are also containers of data. Segments can hold raw signal data like an AnalogSignal
+# Segments can also hold spiketrain data (in a SpikeTrain). It can also hold event data (which
+# we are not show in this example)
+
+
 # Create dummy data, one segment at a time
 for segment in block.segments:
 
-    # create two 4-channel AnalogSignals with dummy data
+    # create two 4-channel AnalogSignals with simulated data (because we have two tetrodes!)
+    # note that the AnalogSignal with have numpy array-like data with units and sampling rates
+    # Neo keeps track of these units while also giving you the flexibility of treating the raw data
+    # like a numpy array
     signals = {
         "Tetrode #1": AnalogSignal(np.random.rand(1000, 4) * mV, sampling_rate=10 * kHz, tetrode_id="Tetrode #1"),
         "Tetrode #2": AnalogSignal(np.random.rand(1000, 4) * mV, sampling_rate=10 * kHz, tetrode_id="Tetrode #2"),
     }
     if store_signals:
         segment.analogsignals.extend(signals.values())
 
-    # create spike trains with dummy data
-    # we will pretend the spikes have been extracted from the dummy signal
+    # create spike trains with simulated data
+    # we will pretend the spikes have been extracted from the simulated signal
     for tetrode_id in ("Tetrode #1", "Tetrode #2"):
         for i in range(n_units[tetrode_id]):
             spiketrain = SpikeTrain(np.random.uniform(0, 100, size=30) * ms, t_stop=100 * ms)
@@ -60,6 +82,12 @@
 
 ###################################################
 # Now we will plot the data
+# Neo doesn't provide it's own plotting functions, but
+# since its data can be treated like numpy arrays
+# it is easy to use standard packages like matplotlib
+# for all your plotting needs
+# We do a classic in neuroscience and show various ways 
+# to plot a PSTH (Peristimulus histogram)
 
 ###################################################
 # .. by trial
@@ -71,6 +99,7 @@
     count, bins = np.histogram(stlist)
     plt.bar(bins[:-1], count, width=bins[1] - bins[0])
     plt.title(f"PSTH in segment {seg.index}")
+plt.tight_layout()
 plt.show()
 
 ####################################################
@@ -83,6 +112,7 @@
     count, bins = np.histogram(stlist)
     plt.bar(bins[:-1], count, width=bins[1] - bins[0])
     plt.title(f"PSTH of unit {group.name}")
+plt.tight_layout()
 plt.show()
 
 ###########################################################
@@ -97,4 +127,5 @@
     count, bins = np.histogram(stlist)
     plt.bar(bins[:-1], count, width=bins[1] - bins[0])
     plt.title(f"PSTH blend of tetrode {tetrode_id}")
+plt.tight_layout()
 plt.show()

examples/plot_read_files_neo_rawio.py

@@ -11,21 +11,24 @@
 import urllib
 from neo.rawio import PlexonRawIO
 
-url_repo = "https://web.gin.g-node.org/NeuralEnsemble/ephy_testing_data/raw/master/"
-
 ##############################################################
 # Get Plexon files
-# We will be pulling these files down, but if you have local file
+
+# We will be pulling these files down, but if you have a local file
 # then all you need to do is specify the file location on your
-# computer
+# computer. NeuralEnsemble keeps a wide variety of freely accesible, small
+# test files that can be used. So for this example we will take advantage of that
+# fact.
 
+url_repo = "https://web.gin.g-node.org/NeuralEnsemble/ephy_testing_data/raw/master/"
 distantfile = url_repo + "plexon/File_plexon_3.plx"
 localfile = "File_plexon_3.plx"
 urllib.request.urlretrieve(distantfile, localfile)
 
 ###############################################################
 # Create a reader
-# All it takes to create a reader is giving the filename
+
+# All it takes to create a reader is giving the filename (or dirname)
 # Then we need to do the slow step of parsing the header with the
 # `parse_header` function. This collects metadata as well as
 # make all future steps much faster for us
@@ -42,56 +45,70 @@
 # we can use None to mean look at all channels. We also need to 
 # specify the block of data (block_index) as well as the segment
 # (seg_index). Then we give the index start and stop. Since we 
-# often thing in time to go from time to index would just require
+# often think in time: to go from time to index would just require
 # the sample rate (so index = time / sampling_rate)
+
 channel_indexes = None  # could be channel_indexes = [0]
 raw_sigs = reader.get_analogsignal_chunk(
     block_index=0, seg_index=0, i_start=1024, i_stop=2048, channel_indexes=channel_indexes
 )
 
 # raw_sigs are not voltages so to convert to voltages we do the follwing
 float_sigs = reader.rescale_signal_raw_to_float(raw_sigs, dtype="float64")
-# each rawio gives you access to lots of information about your data
+
+# We can see that the shapes are the same, but that the datatypes
+# are different once we've rescaled our data
+print("Raw Data: ", raw_sigs.shape, raw_sigs.dtype)
+print("Scaled Data: ", float_sigs.shape, float_sigs.dtype)
+
+###############################################################
+# Each rawio gives you access to lots of information about your data
+# some of this information comes from functions
+# other information is stored as metadata in the reader.header
+
 sampling_rate = reader.get_signal_sampling_rate()
+# Like above we need to indicate the block and segment
 t_start = reader.get_signal_t_start(block_index=0, seg_index=0)
 units = reader.header["signal_channels"][0]["units"]
+
 # and we can display all of this information
-print(raw_sigs.shape, raw_sigs.dtype)
-print(float_sigs.shape, float_sigs.dtype)
 print(f"{sampling_rate=}, {t_start=}, {units=}")
 
 
 ####################################################################
-# Some rawio's and file formats all provide information about spikes
-# If an rawio can't read this data it will raise an error, but lucky
+# Some rawio's and file formats also provide information about spikes
+# If a rawio can't read this data it will raise an error, but lucky
 # for us PlexonRawIO does have spikes data!!
 
 # Count units and spikes per unit
 nb_unit = reader.spike_channels_count()
-print("nb_unit", nb_unit)
+print(f"nb_unit: {nb_unit}\n") # nb_unit stands for number of units
+print("spike_channel_index     nb_spike")
 for spike_channel_index in range(nb_unit):
     nb_spike = reader.spike_count(block_index=0, seg_index=0, spike_channel_index=spike_channel_index)
-    print(f"{spike_channel_index=}\n{nb_spike}\n")
+    print(f"{spike_channel_index}: {nb_spike}\n")
 
 # Read spike times and rescale (just like analogsignal above)
 spike_timestamps = reader.get_spike_timestamps(
     block_index=0, seg_index=0, spike_channel_index=0, t_start=0.0, t_stop=10.0
 )
-print(f"{spike_timestamps.shape=}\n{spike_timestamps.dtype=}\n{spike_timestamps[:5]=}")
+
+print(f"{spike_timestamps.shape=}\n{spike_timestamps.dtype=}\n{spike_timestamps[:5]=}\n")
 spike_times = reader.rescale_spike_timestamp(spike_timestamps, dtype="float64")
-print(f"{spike_times.shape=}\n{spike_times.dtype=}\n{spike_times[:5]}")
+print(f"{spike_times.shape=}\n{spike_times.dtype=}\n{spike_times[:5]}\n")
 
 #######################################################################
 # Some file formats can also give waveform information. We are lucky
-# again. Our file has waveform data!!
+# again our file has waveform data!! We forms are a 3d dataset of 
+# (nb_spike, nb_channel, nb_sample)
 
 # Read spike waveforms
 raw_waveforms = reader.get_spike_raw_waveforms(
     block_index=0, seg_index=0, spike_channel_index=0, t_start=0.0, t_stop=10.0
 )
-print(f"{raw_waveforms.shape=}\n{raw_waveforms.dtype=}\n{raw_waveforms[0, 0, :4]=}")
+print(f"{raw_waveforms.shape=}\n{raw_waveforms.dtype=}\n{raw_waveforms[0, 0, :4]=}\n")
 float_waveforms = reader.rescale_waveforms_to_float(raw_waveforms, dtype="float32", spike_channel_index=0)
-print(f"{float_waveforms.shape=}\n{float_waveforms.dtype=}{float_waveforms[0, 0, :4]=}")
+print(f"{float_waveforms.shape=}\n{float_waveforms.dtype=}{float_waveforms[0, 0, :4]=}\n")
 
 #########################################################################
 # RawIOs can also read event timestamps. But looks like our luck ran out
@@ -102,18 +119,36 @@
 localfile = "File_plexon_2.plx"
 urllib.request.urlretrieve(distantfile, localfile)
 
-# Count events per channel
+#########################################################################
+# Since this is a new file we need to read initialize our reader as well
+# as parse the header
+
 reader = PlexonRawIO(filename="File_plexon_2.plx")
 reader.parse_header()
+# if we look at this header we see it is different than the header above
+print(reader.header)
+
+###########################################################################
+# Now let's look at some event data. This could be things like stimuli applied
+# during the course of an ephys recording
+
 nb_event_channel = reader.event_channels_count()
-print("nb_event_channel:", nb_event_channel)
+print(f"nb_event_channel: {nb_event_channel}")
+# now iterate through the channels 
 for chan_index in range(nb_event_channel):
     nb_event = reader.event_count(block_index=0, seg_index=0, event_channel_index=chan_index)
-    print("chan_index", chan_index, "nb_event", nb_event)
+    print(f"chan_index: {chan_index} nb_event: {nb_event}\n")
+
+
+###############################################################################
+# Finally we can get our actual event timestamps. Some file formats provide the
+# real timestamps (timestamps in s/ms) others have raw timestamps (in samples)
+# so we can do the same style of functions as above. Get the raw timestamps
+# and convert to real times with a rescale function.
 
 ev_timestamps, ev_durations, ev_labels = reader.get_event_timestamps(
     block_index=0, seg_index=0, event_channel_index=0, t_start=None, t_stop=None
 )
-print(f"{ev_timestamps=}\n{ev_durations=}\n{ev_labels=}")
+print(f"{ev_timestamps=}\n{ev_durations=}\n{ev_labels=}\n")
 ev_times = reader.rescale_event_timestamp(ev_timestamps, dtype="float64")
-print(f"{ev_times=}")
+print(f"{ev_times=}\n")