update required package versions and adjust codebase

README.md

@@ -89,7 +89,7 @@ The script produces a figure saved in three different graphics file formats.
 Contains the file `reachgraspio.py`, which contains the loading routine specific to the Reach-to-Grasp experiments in this repository. This loading routine merges the recorded data with metadata information from the odML files into a common Neo object. It is recommended that this loading routine is used in combination with the odML and Neo libraries (see below) to work on the data.
 
 ### Python dependencies
-The provided code is based on the Python libraries Neo[1], odML[2], and Elephant[3], that are required by the example scripts and the reachgraspio loading routine. In short, Neo provides the data model, generic Blackrock loading routines, and APIs used to load the data; odML provides an API to handle the metadata files; and Elephant is a library for the analysis of neuronal data based on the Neo data model that is used by the example script for filtering raw signals to obtain offline filtered LFPs. To run the provided code install the Python dependencies listed in `code/requirements.txt`. Note however, that the versions specified in the `code/requirements.txt` are fixed to match the `reachgraspio.py` loading routine (see Updates below) and more recent versions might exist.
+The provided code is based on the Python libraries Neo[1], odML[2], and Elephant[3], that are required by the example scripts and the reachgraspio loading routine. In short, Neo provides the data model, generic Blackrock loading routines, and APIs used to load the data; odML provides an API to handle the metadata files; and Elephant is a library for the analysis of neuronal data based on the Neo data model that is used by the example script for filtering raw signals to obtain offline filtered LFPs. To run the provided code install the Python dependencies listed in `code/requirements.txt`. Note however, that the versions specified in the `code/requirements.txt` are fixed to match the provided version of the`reachgraspio.py` loading routine (see Updates below) and more recent versions might exist.
 * [1] https://github.com/NeuralEnsemble/python-neo
 * [2] https://github.com/G-Node/python-odml
 * [3] https://github.com/NeuralEnsemble/elephant

code/data_overview_1.py

@@ -110,7 +110,7 @@ def get_arraygrid(signals, chosen_el):
         connector_aligned_id = signals.array_annotations['connector_aligned_ids'][sig_idx]
         x, y = int((connector_aligned_id -1)// 10), int((connector_aligned_id - 1) % 10)
 
-        if signals.array_annotations['channel_ids'][sig_idx] == chosen_el:
+        if np.asarray(signals.array_annotations['channel_ids'][sig_idx], dtype=int) == chosen_el:
             array_grid[x, y] = -0.7
         elif rejections[sig_idx]:
             array_grid[x, y] = -0.35
@@ -529,7 +529,8 @@ assert len(el_raw_sig) == 1
 el_raw_sig = el_raw_sig[0]
 
 # plotting raw signal trace of chosen electrode
-chosen_el_idx = np.where(el_raw_sig.array_annotations['channel_ids'] == chosen_el[monkey])[0][0]
+chids = np.asarray(el_raw_sig.array_annotations['channel_ids'], dtype=int)
+chosen_el_idx = np.where(chids == chosen_el[monkey])[0][0]
 ax4.plot(el_raw_sig.times.rescale(plotting_time_unit),
          el_raw_sig[:, chosen_el_idx].squeeze().rescale(raw_signal_unit),
          color='k')
@@ -591,8 +592,8 @@ ax4.text(timebar_xmin + 0.25 * pq.s, timebar_ypos + timebar_labeloffset,
 # =============================================================================
 # get behavioral signals
 ainp_signals = [nsig for nsig in trial_segment.analogsignals if not nsig.annotations['neural_signal']][0]
-
-force_channel_idx = np.where(ainp_signals.array_annotations['channel_ids'] == 141)[0][0]
+chids = np.asarray(ainp_signals.array_annotations['channel_ids'], dtype=int)
+force_channel_idx = np.where(chids == 141)[0][0]
 ainp_trialz_signals = [a for a in trialz_seg.analogsignals if not a.annotations['neural_signal']]
 assert len(ainp_trialz_signals)
 ainp_trialz = ainp_trialz_signals[0][:, force_channel_idx]
@@ -631,14 +632,14 @@ fploff_trz = eventz[stopidx].rescale(plotting_time_unit)
 # plotting grip force and object displacement
 ai_legend = []
 ai_legend_txt = []
-for chidx, chid in enumerate(ainp_signals.array_annotations['channel_ids']):
+for chidx, chid in enumerate(np.asarray(ainp_signals.array_annotations['channel_ids'], dtype=int)):
     ainp = ainp_signals[:, chidx]
-    if ainp.array_annotations['channel_ids'][0] in trialx_chids:
+    if int(ainp.array_annotations['channel_ids'][0]) in trialx_chids:
         ainp_times = ainp.times.rescale(plotting_time_unit)
         mask = (ainp_times > sr) & (ainp_times < stop)
         ainp_ampli = stats.zscore(ainp.magnitude[mask])
 
-        if ainp.array_annotations['channel_ids'][0] != 143:
+        if int(ainp.array_annotations['channel_ids'][0]) != 143:
             color = 'gray'
             ai_legend_txt.append('grip force')
         else:
@@ -648,7 +649,7 @@ for chidx, chid in enumerate(ainp_signals.array_annotations['channel_ids']):
             ax5a.plot(ainp_times[mask], ainp_ampli, color=color)[0])
 
     # get force load of this trial for next plot
-    elif ainp.array_annotations['channel_ids'][0] == 141:
+    elif int(ainp.array_annotations['channel_ids'][0]) == 141:
         ainp_times = ainp.times.rescale(plotting_time_unit)
         mask = (ainp_times > fplon) & (ainp_times < fploff)
         force_av_01 = np.mean(ainp.rescale(behav_signal_unit).magnitude[mask])

code/data_overview_2.py

@@ -141,7 +141,8 @@ seg.analogsignals = [raw_signal]
 seg = load_segment(seg, load_wavefroms=True)
 
 # only keep the chosen electrode signal in the AnalogSignal object
-mask = np.isin(seg.analogsignals[0].array_annotations['channel_ids'], chosen_els[monkey])
+mask = np.isin(np.asarray(seg.analogsignals[0].array_annotations['channel_ids'], dtype=int),
+               chosen_els[monkey])
 
 seg.analogsignals[0] = seg.analogsignals[0][:, mask]
 
@@ -247,7 +248,8 @@ selected_trial = cut_segments[trial_index]
 for el_idx, electrode_id in enumerate(chosen_els[monkey]):
 
     # PLOT ANALOGSIGNALS in upper plot
-    chosen_el_idx = np.where(cut_segments[0].analogsignals[0].array_annotations['channel_ids'] == electrode_id)[0][0]
+    chids = np.asarray(cut_segments[0].analogsignals[0].array_annotations['channel_ids'], dtype=int)
+    chosen_el_idx = np.where(chids == electrode_id)[0][0]
     anasig = selected_trial.analogsignals[0][:, chosen_el_idx]
     ax1.plot(anasig.times.rescale(time_unit),
              np.asarray(anasig.rescale(lfp_unit))
@@ -300,7 +302,8 @@ ax3.set_xlabel('time [%s]' % time_unit, fontdict=fontdict_axis)
 
 
 # plot data for each chosen trial
-chosen_el_idx = np.where(cut_segments[0].analogsignals[0].array_annotations['channel_ids'] == chosen_el[monkey])[0][0]
+chids = np.asarray(cut_segments[0].analogsignals[0].array_annotations['channel_ids'], dtype=int)
+chosen_el_idx = np.where(chids == chosen_el[monkey])[0][0]
 for trial_idx, trial_id in enumerate(trial_indexes):
     trial_spikes = cut_segments[trial_id].filter(channel_id=chosen_el[monkey], objects='SpikeTrain')
     trial_type = cut_segments[trial_id].annotations['trialtype']

code/example.py

@@ -105,9 +105,12 @@ for anasig in data_segment.analogsignals:
         continue
 
     # identify nsx source of signals in this AnalogSignal object
-    nsx = np.unique(anasig.array_annotations['nsx'])
-    assert len(nsx) == 1, 'Different nsx sources in AnalogSignal'
-    nsx = nsx[0]
+    if 'nsx' in anasig.annotations:
+        nsx = anasig.annotations['nsx']
+    else:
+        nsx = np.unique(anasig.array_annotations['nsx'])
+        assert len(nsx) == 1, 'Different nsx sources in AnalogSignal'
+        nsx = nsx[0]
 
     if nsx == 2:
         # AnalogSignal is LFP from ns2
@@ -121,14 +124,19 @@ for anasig in data_segment.analogsignals:
 # Create LFP signal by filtering raw signal if not present already
 if filtered_anasig is None:
     # Use the Elephant library to filter the signal, filter only target channel
-    target_channel_index = np.where(target_channel_id == raw_anasig.array_annotations['channel_ids'])[0]
+    channel_ids = np.asarray(raw_anasig.array_annotations['channel_ids'], dtype=int)
+    target_channel_index = np.where(target_channel_id == channel_ids)[0]
     raw_signal = raw_anasig[:, target_channel_index]
     f_anasig = butter(
         raw_signal,
         highpass_freq=None,
         lowpass_freq=250 * pq.Hz,
         order=4)
-    f_anasig.name = f'LFP (offline filtered ns{anasig.array_annotations["nsx"][0]})'
+    if 'nsx' in anasig.annotations:
+        nsx = anasig.annotations['nsx']
+    else:
+        nsx = anasig.array_annotations["nsx"][0]
+    f_anasig.name = f'LFP (offline filtered ns{nsx})'
     # Attach all offline filtered LFPs to the segment of data
     data_segment.analogsignals.append(f_anasig)
 
@@ -221,8 +229,12 @@ nsx_colors = {2: 'k', 5: 'r', 6: 'b'}
 for i, anasig in enumerate(trial_segment.analogsignals):
     # only visualize neural data
     if anasig.annotations['neural_signal']:
-        nsx = anasig.array_annotations['nsx'][0]
-        target_channel_index = np.where(anasig.array_annotations['channel_ids'] == target_channel_id)[0]
+        if 'nsx' in anasig.annotations:
+            nsx = anasig.annotations['nsx']
+        else:
+            nsx = anasig.array_annotations['nsx'][0]
+        channel_ids = np.asarray(anasig.array_annotations['channel_ids'], dtype=int)
+        target_channel_index = np.where(channel_ids == target_channel_id)[0]
         target_signal = anasig[:, target_channel_index]
         plt.plot(
             target_signal.times.rescale(time_unit),

code/reachgraspio/reachgraspio.py

@@ -884,7 +884,7 @@ class ReachGraspIO(BlackrockIO):
         """
 
         if self.odmldoc:
-            chids = asig.array_annotations['channel_ids']
+            chids = np.asarray(asig.array_annotations['channel_ids'], dtype=int)
             neural_chids = [chid in self.avail_electrode_ids for chid in chids]
 
             if not any(neural_chids):
@@ -894,7 +894,11 @@ class ReachGraspIO(BlackrockIO):
 
                 # Annotate filter settings from odML
                 nchan = asig.shape[-1]
-                filter = 'Filter_ns%i' % asig.array_annotations['nsx'][0]
+                if 'nsx' in asig.annotations:
+                    nsx = asig.annotations['nsx']
+                else:
+                    nsx = asig.array_annotations['nsx'][0]
+                filter = 'Filter_ns%i' % nsx
                 sec = self.odmldoc['Cerebus']['NeuralSignalProcessor']['NeuralSignals'][filter]
                 props = sec.properties
                 hi_pass_freq = np.full((nchan), pq.Quantity(props['HighPassFreq'].values[0],
@@ -905,7 +909,7 @@ class ReachGraspIO(BlackrockIO):
                 lo_pass_order = np.zeros_like(lo_pass_freq)
                 filter_type = np.empty((nchan), np.str)
                 for chidx in range(nchan):
-                    filter = 'Filter_ns%i' % asig.array_annotations['nsx'][chidx]
+                    filter = 'Filter_ns%i' % nsx
                     sec = self.odmldoc['Cerebus']['NeuralSignalProcessor']['NeuralSignals'][filter]
                     hi_pass_freq[chidx] = pq.Quantity(
                         sec.properties['HighPassFreq'].values[0],
@@ -970,8 +974,8 @@ class ReachGraspIO(BlackrockIO):
             for seg in block.segments:
                 for obj in seg.analogsignals:
                     if 'neural_signal' in obj.annotations and obj.annotations[
-                        'neural_signal'] and 'channel_ids' in obj.array_annotations:
-                        chids = obj.array_annotations['channel_ids']
+                       'neural_signal'] and 'channel_ids' in obj.array_annotations:
+                        chids = np.asarray(obj.array_annotations['channel_ids'], dtype=int)
                         ca_ids, *coordinates = self.__convert_chids_and_coordinates(chids)
                         obj.array_annotations.update(dict(connector_aligned_ids=ca_ids,
                                                           coordinates_x=coordinates[0],
@@ -1032,21 +1036,23 @@ class ReachGraspIO(BlackrockIO):
             arraygrid_row_num=sec.properties['GridRows'].values,
             arraygrid_col_num=sec.properties['GridColumns'].values)
 
-        secs = self.odmldoc['UtahArray']['Array'].sections
         bl.annotate(avail_electrode_ids=self.avail_electrode_ids)
 
-        # TODO: add list of behavioral channels
-        # bl.annotate(avail_behavsig_indexes=[])
-
     def __correct_filter_shifts(self, asig):
         if self.odmldoc and asig.annotations['neural_signal']:
             # assert all signals are originating from same nsx file
-            if len(np.unique(asig.array_annotations['nsx'])) > 1:
+            if 'nsx' in asig.array_annotations and \
+                    len(np.unique(asig.array_annotations['nsx'])) > 1:
                 raise ValueError('Multiple nsx file origins (%s) in single AnalogSignal'
                                  ''.format(asig.array_annotations['nsx']))
 
+            if 'nsx' in asig.annotations:
+                nsx = asig.annotations['nsx']
+            else:
+                nsx = asig.array_annotations['nsx'][0]
+
             # Get and correct for shifts
-            filter_name = 'Filter_ns%i' % asig.array_annotations['nsx'][0]  # use nsx of 1st signal
+            filter_name = 'Filter_ns%i' % nsx  # use nsx of 1st signal
             sec = self.odmldoc['Cerebus']['NeuralSignalProcessor']['NeuralSignals'][filter_name]
             shift = pq.Quantity(
                 sec.properties['EstimatedShift'].values[0],

code/requirements.txt

@@ -1,12 +1,12 @@
-numpy>=1.8.2
-quantities>=0.10.1
-scipy>=0.14.0
-six
+numpy>=1.8.2, <=1.22.3
+quantities>=0.10.1, <=0.13.0
+scipy>=0.14.0, <=1.8.0
+six; python_version >= '3.3'
 enum34; python_version < '3.4'
 
 # Plotting module required for example scripts
 matplotlib>=3.0.3
 
-neo==0.9.0
-elephant==0.9.0
-odml==1.5.1
+neo>=0.9.0, <0.11
+elephant>=0.9.0, <0.12
+odml>=1.5.1, <1.6