|
| 1 | +How to perform inter-session alignment |
| 2 | +====================================== |
| 3 | + |
| 4 | +In this how-to, we will assess and correct changes in probe position across multiple experimental sessions |
| 5 | +using 'inter-session alignment'. |
| 6 | + |
| 7 | +This is often valuable chronic-recording experiments, where the goal is to track units across sessions |
| 8 | + |
| 9 | +A full tutorial, including details on the many settings for this procedure, can be found `here <TODO: ADD LINK (INTERNAL)>`_. |
| 10 | + |
| 11 | +Running inter-session alignment |
| 12 | +------------------------------- |
| 13 | + |
| 14 | +In SpikeInterface, it is recommended to perform inter-session alignment |
| 15 | +following within-session motion correction (if used) and before whitening. |
| 16 | + |
| 17 | +Preprocessed recordings should be stored in a list before being passed |
| 18 | +to the `session_alignment` functions: |
| 19 | + |
| 20 | +.. code:: python |
| 21 | +
|
| 22 | + recordings_list = [prepro_session_1, prepro_session_2, ...] |
| 23 | +
|
| 24 | +Here, we will simulate such an experiment by generating a pair of sessions in |
| 25 | +which the probe is displaced 200 micrometers (um) along its y-axis (depth). |
| 26 | + |
| 27 | +The first step involves running all required imports: |
| 28 | + |
| 29 | +.. code:: python |
| 30 | +
|
| 31 | + from spikeinterface.generation.session_displacement_generator import generate_session_displacement_recordings |
| 32 | + import spikeinterface.full as si |
| 33 | + from spikeinterface.preprocessing.inter_session_alignment import ( # TODO: should add all of the below to spikeinterface.full? (CHECK) |
| 34 | + session_alignment, |
| 35 | + ) |
| 36 | + from spikeinterface.widgets import plot_session_alignment, plot_activity_histogram_2d |
| 37 | + import matplotlib.pyplot as plt |
| 38 | +
|
| 39 | +
|
| 40 | +and then generating the test recordings: |
| 41 | + |
| 42 | +.. code:: python |
| 43 | +
|
| 44 | + recordings_list, _ = generate_session_displacement_recordings( # TODO: add to spikeinterface.full ? |
| 45 | + num_units=5, |
| 46 | + recording_durations=[10, 10], |
| 47 | + recording_shifts=((0, 0), (0, 200)), # (x offset, y offset) pairs |
| 48 | + ) |
| 49 | +
|
| 50 | +We won't explicitly preprocess these recordings in this how-to, but you can imagine |
| 51 | +preprocessing steps have already been run (e.g. filtering, common reference etc.). |
| 52 | + |
| 53 | +To run inter-session alignment, we need to detect peaks and compute the peak locations, |
| 54 | +as the location of firing neurons are used as anchors to align the sessions. |
| 55 | + |
| 56 | +If you are **running inter-session alignment following motion correction**, the peaks will |
| 57 | +already be detected and localised. In this case, please jump to the |
| 58 | +:ref:`alignment guide <_with_motion_correction>`. |
| 59 | + |
| 60 | +In this section we will imagine motion correction was not run, so we need to compute the peaks: |
| 61 | + |
| 62 | +.. code:: python |
| 63 | +
|
| 64 | + peaks_list, peak_locations_list = session_alignment.compute_peaks_locations_for_session_alignment( |
| 65 | + recordings_list, |
| 66 | + detect_kwargs={"method": "locally_exclusive"}, |
| 67 | + localize_peaks_kwargs={"method": "grid_convolution"}, |
| 68 | + ) |
| 69 | +
|
| 70 | +The peak locations (before correction) can be visualised with the plotting function: |
| 71 | + |
| 72 | +.. code:: python |
| 73 | +
|
| 74 | + plot_session_alignment( |
| 75 | + recordings_list, |
| 76 | + peaks_list, |
| 77 | + peak_locations_list, |
| 78 | + drift_raster_map_kwargs={"clim":(-250, 0)} # fix the color limit across plots for easy comparison |
| 79 | + ) |
| 80 | + plt.show() |
| 81 | +
|
| 82 | +Now, we are ready to perform inter-session alignment. There are many options associated |
| 83 | +with this method—the simplest way to edit these is to fetch the default options |
| 84 | +using the getter function as below: |
| 85 | + |
| 86 | +.. code:: python |
| 87 | +
|
| 88 | + estimate_histogram_kwargs = session_alignment.get_estimate_histogram_kwargs() |
| 89 | + estimate_histogram_kwargs["histogram_type"] = "activity_2d" # TODO: RENAME |
| 90 | +
|
| 91 | + corrected_recordings_list, extra_info = session_alignment.align_sessions( |
| 92 | + recordings_list, |
| 93 | + peaks_list, |
| 94 | + peak_locations_list, |
| 95 | + estimate_histogram_kwargs=estimate_histogram_kwargs |
| 96 | + ) |
| 97 | +
|
| 98 | +To assess the performance of inter-session alignment, `plot_session_alignment()` |
| 99 | +will plot both the original and corrected recordings: |
| 100 | + |
| 101 | +.. code:: python |
| 102 | +
|
| 103 | + plot_session_alignment( # TODO: is this signature confusing? |
| 104 | + recordings_list, |
| 105 | + peaks_list, |
| 106 | + peak_locations_list, |
| 107 | + extra_info["session_histogram_list"], |
| 108 | + **extra_info["corrected"], |
| 109 | + spatial_bin_centers=extra_info["spatial_bin_centers"], |
| 110 | + drift_raster_map_kwargs={"clim":(-250, 0)} |
| 111 | + ) |
| 112 | + plt.show() |
| 113 | +
|
| 114 | +As we have used 2d histograms for alignment, we can also plot these with ``plot_activity_histogram_2d()``. |
| 115 | + |
| 116 | +.. _with_motion_correction: |
| 117 | + |
| 118 | +Inter-session alignment after motion correction |
| 119 | +----------------------------------------------- |
| 120 | + |
| 121 | +If motion correction has already been performed, it is possible to reuse the |
| 122 | +previously computed peaks and peak locations, avoiding the need for re-computation. |
| 123 | +We will use the special function `align_sessions_after_motion_correction()` for this case. |
| 124 | + |
| 125 | +Critically, the last preprocessing step prior to inter-session alignment should be motion correction, |
| 126 | +so the correction for inter-session displacement will be **added directly to the motion correction**. |
| 127 | +This is beneficial as it avoids interpolating the data (i.e. shifting the traces) more than once. |
| 128 | + |
| 129 | +.. admonition:: Warning |
| 130 | + :class: warning |
| 131 | + |
| 132 | + To ensure that inter-session alignment adds the displacement directly to the motion-corrected recording |
| 133 | + to avoid repeated interpolation, motion correction must be the final operation applied to the recording |
| 134 | + prior to inter-session alignment. |
| 135 | + |
| 136 | + You can verify this by confirming the recording is an ``InterpolateMotionRecording`` with: |
| 137 | + |
| 138 | + .. code:: python |
| 139 | +
|
| 140 | + type(recording) # quick check, should print `InterpolateMotionRecording` |
| 141 | +
|
| 142 | + from spikeinterace.sortingcomponents.motion.motion_utils import InterpolateMotionRecording |
| 143 | +
|
| 144 | + assert isinstance(recording, InterpolateMotionRecording) # error if not true |
| 145 | +
|
| 146 | + ``align_sessions_after_motion_correction()`` will raise an error if the passed recordings |
| 147 | + are not all `InterpolateMotionRecordings`. |
| 148 | + |
| 149 | +Again, let's create some test data. We can create a recording with motion errors, |
| 150 | +then split it in two to simulate two separate sessions: |
| 151 | + |
| 152 | +.. code:: python |
| 153 | +
|
| 154 | + # Generate the recording with motion artefact |
| 155 | + motion_recording, _ = si.generate_drifting_recording(duration=100) |
| 156 | + total_duration = motion_recording.get_duration() |
| 157 | + split_time = total_duration / 2 |
| 158 | +
|
| 159 | + # Split in two to simulate two sessions |
| 160 | + recording_part1 = motion_recording.time_slice(start_time=0, end_time=split_time) |
| 161 | + recording_part2 = motion_recording.time_slice(start_time=split_time, end_time=total_duration) |
| 162 | +
|
| 163 | +
|
| 164 | +Next, motion correction is performed, storing the results in a list: |
| 165 | + |
| 166 | +.. code:: python |
| 167 | +
|
| 168 | + # perform motion correction on each session, storing the outputs in lists |
| 169 | + recordings_list = [] |
| 170 | + motion_info_list = [] |
| 171 | + for recording in [recording_part1, recording_part2]: |
| 172 | +
|
| 173 | + rec, motion_info = si.correct_motion(recording, output_motion_info=True, preset="rigid_fast") |
| 174 | +
|
| 175 | + recordings_list.append(rec) |
| 176 | + motion_info_list.append(motion_info) |
| 177 | +
|
| 178 | +Now, we are ready to use ``align_sessions_after_motion_correction()`` |
| 179 | +to align the motion-corrected sessions. |
| 180 | + |
| 181 | +This function should always be used for aligning motion-corrected sessions, |
| 182 | +as it ensures the alignment parameters are properly matched. |
| 183 | + |
| 184 | +We can pass any arguments directly to ``align_sessions`` using the ``align_sessions_kwargs`` argument: |
| 185 | + |
| 186 | +.. code:: python |
| 187 | +
|
| 188 | + estimate_histogram_kwargs = session_alignment.get_estimate_histogram_kwargs() |
| 189 | + estimate_histogram_kwargs["histogram_type"] = "activity_2d" # TODO: RENAME |
| 190 | +
|
| 191 | + align_sessions_kwargs = {"estimate_histogram_kwargs": estimate_histogram_kwargs} |
| 192 | +
|
| 193 | + corrected_recordings_list, extra_info = session_alignment.align_sessions_after_motion_correction( |
| 194 | + recordings_list, motion_info_list, align_sessions_kwargs |
| 195 | + ) |
| 196 | +
|
| 197 | +As above, the inter-session alignment can be assessed using ``plot_session_alignment()``. |
0 commit comments