Changed preprocessing appropriately

bci_essentials/paradigm/p300_paradigm.py

@@ -15,6 +15,7 @@ def __init__(
         epoch_start=0,
         epoch_end=0.6,
         buffer_time=0.01,
+        preprocessing_window=2,
     ):
         """
         Parameters
@@ -34,6 +35,9 @@ def __init__(
         buffer_time : float, *optional*
             Defines the time in seconds after an epoch for which we require EEG data to ensure that all EEG is present in that epoch.
             - Default is `0.01`.
+        preprocessing_window : float, *optional*
+            Defines the time in seconds before and after a marker to include additional EEG data in preprocessing.
+            - Default is `2`.
         """
 
         super().__init__(filters)
@@ -54,6 +58,8 @@ def __init__(
 
         self.buffer_time = buffer_time
 
+        self.preprocessing_window = preprocessing_window
+
     def get_eeg_start_and_end_times(self, markers, timestamps):
         """
         Get the start and end times of the EEG data based on the markers.
@@ -129,8 +135,10 @@ def process_markers(self, markers, marker_timestamps, eeg, eeg_timestamps, fsamp
             # Subtract the marker timestamp from the EEG timestamps so that 0 becomes the marker onset
             marker_eeg_timestamps = eeg_timestamps - marker_timestamp
 
-            # Create the epoch time vector
-            epoch_time = np.arange(self.epoch_start, self.epoch_end, 1 / fsample)
+            # Create the epoch time vector using the preprocessing window
+            epoch_time = np.arange(
+                -self.preprocessing_window, self.preprocessing_window, 1 / fsample
+            )
 
             epoch_X = np.zeros((1, n_channels, len(epoch_time)))
 
@@ -150,23 +158,32 @@ def process_markers(self, markers, marker_timestamps, eeg, eeg_timestamps, fsamp
                 epoch_X[0, :, :], fsample, self.lowcut, self.highcut
             )
 
+            # Trim the preprocessed epoch
+            trimmed_epoch_X = epoch_X[
+                :,
+                :,
+                ((epoch_time >= self.epoch_start) & (epoch_time <= self.epoch_end)),
+            ]
+
             # For each flash index in the marker
             for flash_index in flash_indices:
                 if flash_counts[flash_index] == 0:
-                    object_epochs[flash_index] = epoch_X
+                    object_epochs[flash_index] = trimmed_epoch_X
                     flash_counts[flash_index] += 1
                 else:
                     object_epochs[flash_index] = np.concatenate(
-                        (object_epochs[flash_index], epoch_X), axis=0
+                        (object_epochs[flash_index], trimmed_epoch_X), axis=0
                     )
                     flash_counts[flash_index] += 1
 
         # Average all epochs for each object
-        object_epochs_mean = [np.zeros((n_channels, len(epoch_time)))] * num_objects
+        object_epochs_mean = [
+            np.zeros((n_channels, trimmed_epoch_X.shape[-1]))
+        ] * num_objects
         for i in range(num_objects):
             object_epochs_mean[i] = np.mean(object_epochs[i], axis=0)
 
-        X = np.zeros((num_objects, n_channels, len(epoch_time)))
+        X = np.zeros((num_objects, n_channels, trimmed_epoch_X.shape[-1]))
         for i in range(num_objects):
             X[i, :, :] = object_epochs_mean[i]
         # # object_epochs_mean = np.mean(object_epochs, axis=1)

bci_essentials/paradigm/paradigm.py

@@ -36,23 +36,18 @@ def __init__(self, filters=[5, 30], channel_subset=None):
     def _preprocess(self, eeg, fsample, lowcut, highcut, order=5):
         """
         Preprocess EEG data with the appropriate filter type:
-        - If the data is continuous (i.e., shape is [channels, samples]), a
-        bandpass filter is used.
-
-        - If the data is epoched (i.e., shape is [epoch, channels, samples]),
-        the filter type depends on the signal length relative to the filter's settling time:
-            - If signal length > settling time: use bandpass filter
-            - If signal length ≤ settling time: use lowpass filter
-            - The settling time is calculated by:
-                1. Compute the time constant (tc) of the highpass filter:
-                    tc = 1 / (2 * π * lowcut)
-                2. Compute the settling time with the formula:
-                    settling_time = tc * 5 * order
-                    - In a first-order system, a rule-of-thumb is that the signal settles in approximately 5 time constants (tc * 5)
-                    - In higher-order filters (e.g., a 5th-order filter set as the default), the settling time incrases linearly with the order of the filter (order * tc * 5)
-                    - This is a simplification, but it provides a good approximation for the settling time of the filter.
-                    - For more details, see the reference below:
-                        https://www.analogictips.com/an-overview-of-filters-and-their-parameters-part-4-time-and-phase-issues/
+        - If signal length > settling time: use bandpass filter
+        - If signal length ≤ settling time: use lowpass filter
+        - The settling time is calculated by:
+            1. Compute the time constant (tc) of the highpass filter:
+                tc = 1 / (2 * π * lowcut)
+            2. Compute the settling time with the formula:
+                settling_time = tc * 5 * order
+                - In a first-order system, a rule-of-thumb is that the signal settles in approximately 5 time constants (tc * 5)
+                - In higher-order filters (e.g., a 5th-order filter set as the default), the settling time incrases linearly with the order of the filter (order * tc * 5)
+                - This is a simplification, but it provides a good approximation for the settling time of the filter.
+                - For more details, see the reference below:
+                    https://www.analogictips.com/an-overview-of-filters-and-their-parameters-part-4-time-and-phase-issues/
 
         Parameters
         ----------
@@ -74,31 +69,21 @@ def _preprocess(self, eeg, fsample, lowcut, highcut, order=5):
             Preprocessed EEG. Shape is the same as `eeg`.
 
         """
+        logger.debug("Preprocessing EEG data")
 
-        n_dims = len(eeg.shape)
-        if n_dims == 2:
-            logger.debug("Preprocessing continuous EEG")
+        # Get the length of the signal
+        signal_length = eeg.shape[-1] / fsample
+
+        # Highpass filter settling time
+        tc = 1 / (2 * np.pi * lowcut)
+        settling_time = order * tc * 5
+
+        if signal_length > settling_time:
+            logger.debug("Applied bandpass filter")
             preprocessed_eeg = bandpass(eeg, lowcut, highcut, order, fsample)
-        elif n_dims == 3:
-            logger.debug("Preprocessing epoched EEG")
-
-            # Get the length of the signal
-            signal_length = eeg.shape[2]
-
-            # Highpass filter settling time
-            tc = 1 / (2 * np.pi * lowcut)
-            settling_time = order * tc * 5
-
-            if signal_length > settling_time:
-                logger.debug("Applied bandpass filter to epoched EEG")
-                preprocessed_eeg = bandpass(eeg, lowcut, highcut, order, fsample)
-            else:
-                logger.debug("Applied lowpass filter to epoched EEG")
-                preprocessed_eeg = lowpass(eeg, highcut, order, fsample)
         else:
-            raise ValueError(
-                "Preprocessing failed. EEG must be 2D (continuous) or 3D (epoched)."
-            )
+            logger.debug("Applied lowpass filter")
+            preprocessed_eeg = lowpass(eeg, highcut, order, fsample)
 
         return preprocessed_eeg