PhotometryExperiment

This class is used to process and slice photometry experiments, with support for both dual and single channel experiments.

The two main APIs:

• preprocess_signal() applies a low-pass filter, fits a reference trace, and applies correction methods, with optionally functionality to detect and correct movement artifacts.

• extract_trial_data() uses the events timestamps in the .events attribute to slice, align, and pass down trial-relative event timesteps to a 2D PhotometryData object.

---

Example Usage

Dual channel preprocessing

exp.preprocess_signal(
    # lowpass butterworth params
    cutoff_frequency=3,
    order=4,

    # correction method and isosbestic fit params
    correction_method='dF/F',
    fit_using='IRLS',
    maxiter=1000,
    c=2,
)

Single channel preprocessing

# import artifact handlers
from pyFiberPhotometry.analysis.artifact import ODS_Detector, Spline_Corrector

# instantiate artifact detector and corrector
detector = ODS_Detector(
    score_threshold=5,
    jump_score_threshold=10,
    expand_sec=(0.5, 2),
    buffer_sec=1.5,
    n_chunks=20,
)

corrector = Spline_Corrector(
    anchor_sec=(0.2, 0.2),
    correct_spikes=True,
    correct_jumps=True,
)

exp.preprocess_signal(
    cutoff_frequency=3,
    order=4,
    correction_method='dB/B',
    fit_using='IRLS',
    maxiter=1000,
    c=2,

    # pass in artifact handlers
    artifact_detector=detector,
    artifact_corrector=corrector,
)

Trial extraction

exp.extract_trial_data(
    # what event should we consider the "start" of a trial
    align_to='event',

    # what events do we want to center on
    center_on=['lever1', 'lever2'],

    # how long in seconds should our trials be relative to "center_on"
    trial_bounds=(-10, 10),

    # expected range of our events relative to "align_to"
    event_tolerences={
        'lever1':(2, 10),
        'lever2':(2, 10),
        'loud_noise':(2, 12),
    },

    # which trial-wise normalization should be preformed
    trial_normalization='none',

    # if multiple of the same event are within tolerences which should be picked
    event_conflict_logic='first',

    # should an error be thrown if multiple "center_on" event are present in the same trial
    check_overlap=True,
)

---

PhotometryExperiment(raw_signal, raw_isosbestic, time, events={}, metadata={}, frequency=None)

Handle processing of raw photometry data.

Initialize a photometry experiment.

Parameters:

• raw_signal (ndarray) –

  Raw signal channel values.

• raw_isosbestic (ndarray) –

  Raw isosbestic channel values. If None, experiment is assumed to be single channel.

• time (ndarray) –

  Time points corresponding to the raw signals.

• events (dict, default: {} ) –

  Mapping of event labels to timestamp arrays. Defaults to {}.

• metadata (dict, default: {} ) –

  Additional experiment metadata. Defaults to {}.

• frequency (float | None, default: None ) –

  Sampling frequency in Hz. If None (default), it is estimated from raw_signal and time.

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

def __init__( 
        self,
        raw_signal: np.ndarray,
        raw_isosbestic: np.ndarray | None,
        time: np.ndarray,
        events: dict = {},
        metadata: dict = {},
        frequency: float | None = None,
        ):
    """Initialize a photometry experiment.

    Args:
        raw_signal (np.ndarray): Raw signal channel values.
        raw_isosbestic (np.ndarray): Raw isosbestic channel values. If
            ``None``, experiment is assumed to be single channel.
        time (np.ndarray): Time points corresponding to the raw signals.
        events (dict, optional): Mapping of event labels to timestamp
            arrays. Defaults to ``{}``.
        metadata (dict, optional): Additional experiment metadata. Defaults
            to ``{}``.
        frequency (float | None, optional): Sampling frequency in Hz. If
            ``None`` (default), it is estimated from ``raw_signal`` and
            ``time``.
    """
    # save inputs
    self.raw_signal = raw_signal
    self.raw_isosbestic = raw_isosbestic
    self.time = time
    self.events = events
    self.metadata = metadata
    self.frequency = raw_signal.size / (self.time.max() - self.time.min() + 1) if frequency is None else frequency

    self.signal = None
    self.fitted_ref = None
    self.trial_data = None

    # add relevant metadata
    self.metadata['frequency'] = self.frequency
    self.metadata['channel_mode'] = self.channel_mode

load_TDT(data_folder, box, event_labels, signal_label, isosbestic_label, downsample=10) classmethod

Load photometry data from TDT format.

Parameters:

• data_folder (str) –

  Path to the TDT block folder.

• box (str) –

  TDT box identifier used in stream and epoc labels.

• event_labels (list[str]) –

  Event labels to extract from epocs.

• signal_label (str) –

  Base label for the signal channel.

• isosbestic_label (str) –

  Base label for the isosbestic channel.

• downsample (int, default: 10 ) –

  Downsampling factor for the raw streams (mean pooling). Defaults to 10.

Returns:

• PhotometryExperiment ( PhotometryExperiment ) –

  Loaded experiment instance.

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

@classmethod
def load_TDT(
        cls,
        data_folder: str,
        box: str,
        event_labels: list[str],
        signal_label: str,
        isosbestic_label: str,
        downsample: int = 10,
        ) -> PhotometryExperiment:
    """Load photometry data from TDT format.

    Args:
        data_folder (str): Path to the TDT block folder.
        box (str): TDT box identifier used in stream and epoc labels.
        event_labels (list[str]): Event labels to extract from epocs.
        signal_label (str): Base label for the signal channel.
        isosbestic_label (str): Base label for the isosbestic channel.
        downsample (int, optional): Downsampling factor for the raw streams
            (mean pooling). Defaults to ``10``.

    Returns:
        PhotometryExperiment: Loaded experiment instance.
    """
    from .PhotometryLoader import TDTLoader
    loader = TDTLoader(
        data_folder=data_folder,
        box=box,
        event_labels=event_labels,
        signal_label=signal_label,
        isosbestic_label=isosbestic_label,
        downsample=downsample,
    )
    return loader.load()

run_pipeline()

Run the full processing pipeline for one experiment.

This method is intended to be implemented by child classes.

Returns:

• None –

  None

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

def run_pipeline(self) -> None:
    """Run the full processing pipeline for one experiment.

    This method is intended to be implemented by child classes.

    Returns:
        None
    """
    return

preprocess_signal(cutoff_frequency=3.0, order=4, signal_normalization='none', correction_method='dF/F', fit_using='IRLS', maxiter=1000, c=3, artifact_detector=None, artifact_corrector=None)

Low-pass filter and preprocess the signal using isosbestic fitting.

Parameters:

• cutoff_frequency (float, default: 3.0 ) –

  Low-pass cutoff frequency in Hz. Values between 1 and 5 recommended. Defaults to 3.0.

• order (int, default: 4 ) –

  Butterworth filter order. Values >3 are recommended. Defaults to 4.

• correction_method (Literal['dF/F', 'dF', 'dB/B', 'dB', 'none'] | Callable, default: 'dF/F' ) –

  Reference trace correction method:
  * 'dF/F' and 'dB/B = (signal - fitted reference) / fitted reference * 'dF' and 'dB' = (signal - fitted reference) For dual channel experiments the reference is the isosbestic, for single channel it is a fit photobleaching curve. Use 'dF/F' or 'dF' for dual channel and 'dB/B' or 'dB' for single channel. A custom function that takes in the positional arguements signal, fitted_reference and returns a 1D np.ndarray can also be passed. Defaults to 'dF/F'.

• signal_normalization (Literal['zscore', 'nullZ', 'none'] | Callable, default: 'none' ) –

  Method for whole-signal normalization; 'none' for dF/F and 'zscore' for dF is recommended. A custom function that takes in the positional arguements signal and returns a 1D np.ndarray can also be passed. Defaults to 'none'.

• fit_using (Literal['OLS', 'IRLS', 'IRLS_no_intercept', 'OLS_no_intercept'] | Callable, default: 'IRLS' ) –

  Model used to fit isosbestic to experimental signal. IRLS methods recommended. Use a no intercept type model if large global change is present in the experimental signal. A custom function that takes in the positional arguements signal, isosbestic and returns a 1D np.ndarray and a sequence of params can also be passed. Note custom functions will only apply to isosbestic fits (i.e. dual channel experiments). Defaults to 'IRLS'.

• maxiter (int, default: 1000 ) –

  Maximum iterations of the IRLS isosbestic fit. Defaults to 1000.

• c (float | None, default: 3 ) –

  Constant for IRLS fits; smaller values mean more agressive downweighting. 1.4 <= c <= 3 is recommended unless there is large global drift in the experimental signal, in which case large values (>5) are better. Defaults to 3.

• artifact_detector (ArtifactDetector | None, default: None ) –

  Detector object with method .detect(signal, reference, time) used to detect artifacts. If None, detection is skipped. Defaults to None.

• artifact_corrector (ArtifactCorrector | None, default: None ) –

  Corrector object with method .correct(signal, time, artifacts) used to correct artifacts. If None, correction is skipped. Defaults to None.

Returns:

• None –

  None

Raises:

• ValueError –
  • If correction_method is incompatible with self.channel_mode
  • If artifact_corrector specified but not artifact_detector

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

def preprocess_signal(
        self,
        cutoff_frequency: float = 3.0, 
        order: int = 4,
        signal_normalization: Literal['zscore', 'nullZ', 'none'] | Callable = 'none',
        correction_method: Literal['dF/F', 'dF', 'dB/B', 'dB', 'none'] | Callable = 'dF/F',
        fit_using: Literal['OLS', 'IRLS', 'IRLS_no_intercept', 'OLS_no_intercept'] | Callable = 'IRLS',
        maxiter: int = 1000,
        c: float | None = 3,
        artifact_detector: ArtifactDetector | None = None,
        artifact_corrector: ArtifactCorrector | None = None,
        ) -> None:
    """Low-pass filter and preprocess the signal using isosbestic fitting.

    Args:
        cutoff_frequency (float, optional): Low-pass cutoff frequency in Hz.
            Values between 1 and 5 recommended. Defaults to ``3.0``.

        order (int, optional): Butterworth filter order. 
            Values >3 are recommended. Defaults to ``4``.

        correction_method (Literal['dF/F', 'dF', 'dB/B', 'dB', 'none'] | Callable, optional):
            Reference trace correction method:  
            * ``'dF/F'`` and ``'dB/B`` = (signal - fitted reference) / fitted reference
            * ``'dF'`` and ``'dB'`` = (signal - fitted reference) 
            For dual channel experiments the reference is the isosbestic, for single channel
            it is a fit photobleaching curve. Use ``'dF/F'`` or ``'dF'`` for dual channel
            and ``'dB/B'`` or ``'dB'`` for single channel.
            A custom function that takes in the positional arguements ``signal``, ``fitted_reference``
            and returns a 1D ``np.ndarray`` can also be passed. Defaults to ``'dF/F'``.

        signal_normalization (Literal['zscore', 'nullZ', 'none'] | Callable, optional): 
            Method for whole-signal normalization; ``'none'``
            for dF/F and ``'zscore'`` for dF is recommended. 
            A custom function that takes in the positional arguements ``signal``
            and returns a 1D ``np.ndarray`` can also be passed.
            Defaults to ``'none'``.

        fit_using (Literal['OLS', 'IRLS', 'IRLS_no_intercept', 'OLS_no_intercept'] | Callable, optional): 
            Model used to fit isosbestic to experimental signal. IRLS methods recommended. 
            Use a no intercept type model if large global change is present in the experimental signal.
            A custom function that takes in the positional arguements ``signal``, ``isosbestic``
            and returns a 1D ``np.ndarray`` and a sequence of params can also be passed.
            Note custom functions will only apply to isosbestic fits (i.e. dual channel experiments).
            Defaults to ``'IRLS'``.

        maxiter (int, optional): Maximum iterations of the IRLS isosbestic
            fit. Defaults to ``1000``.

        c (float | None, optional): Constant for IRLS fits; smaller values
            mean more agressive downweighting. ``1.4 <= c <= 3`` is
            recommended unless there is large global drift in the
            experimental signal, in which case large values (>5) are better.
            Defaults to ``3``.

        artifact_detector (ArtifactDetector | None, optional): Detector object with method
            ``.detect(signal, reference, time)`` used to detect artifacts.
            If ``None``, detection is skipped. Defaults to ``None``.

        artifact_corrector (ArtifactCorrector | None, optional): Corrector object with method
            ``.correct(signal, time, artifacts)`` used to correct artifacts.
            If ``None``, correction is skipped. Defaults to ``None``.

    Returns:
        None

    Raises:
        ValueError: 
            * If ``correction_method`` is incompatible with ``self.channel_mode``
            * If ``artifact_corrector`` specified but not ``artifact_detector``
    """
    # validate inputs
    dual_channel_methods = ['dF/F', 'dF']
    single_channel_methods = ['dB/B', 'dB']

    if self.channel_mode == 'single' and correction_method in dual_channel_methods:
        raise ValueError(f'Correction methods {",".join(dual_channel_methods)} are for dual channel experiments.')

    if self.channel_mode == 'dual' and correction_method in single_channel_methods:
        raise ValueError(f'Correction methods {",".join(single_channel_methods)} are for single channel experiments.')

    if artifact_corrector is not None and artifact_detector is None:
        raise ValueError(f'artifact_detector not specified but artifact_corrector is.')

    # apply lowpass butterworth filter
    filt_sig = self.low_frequency_pass_butter(
        self.raw_signal, 
        self.frequency,
        cutoff_frequency=cutoff_frequency,
        order=order
    )

    if self.channel_mode == 'dual':
        # apply lowpass to isosbestic
        filt_iso = self.low_frequency_pass_butter(
            self.raw_isosbestic, 
            self.frequency,
            cutoff_frequency=cutoff_frequency,
            order=order
        )

        # trim to minimum length
        min_len = min(filt_sig.size, filt_iso.size)
        filt_sig = filt_sig[:min_len]
        filt_iso = filt_iso[:min_len]

        # fit isosbestic to signal
        fitted_ref, r2_val, coeffs = self.fit_isosbestic_to_signal(
            filt_sig, 
            filt_iso,
            maxiter=maxiter,
            fit_using=fit_using,
            c=c,
        )
        reference_type = 'isosbestic'

    elif self.channel_mode == 'single':
        # fit photobleaching curve
        fitted_ref, r2_val, coeffs = self.fit_photobleaching_curve(
            signal=filt_sig
        )
        reference_type = 'photobleaching'

    else:
        raise ValueError(f"Channel mode, {self.channel_mode}, not recognized.")

    # save relevant data
    self.metadata['reference_fit'] = dict(
        type = reference_type,
        r2_val = r2_val,
        coeffs = coeffs,
    )
    self.filt_sig = filt_sig
    self.fitted_ref = fitted_ref

    # apply reference curve correction
    signal = self._apply_correction_method(
        correction_method=correction_method,
        signal=filt_sig,
        fitted_ref=fitted_ref,
    )

    # apply signal normalization
    signal = self._apply_signal_normalization(
        signal_normalization=signal_normalization,
        signal=signal,
    )

    # detect artifacts
    artifacts = None
    if artifact_detector is not None:
        reference = None if self.channel_mode == 'single' else self.fitted_ref

        artifacts = self._detect_artifacts(
            detector=artifact_detector,
            signal=signal,
            reference=reference,
            time=self.time,
        )

        # correct artifacts
        if artifact_corrector is not None:
            signal = self._correct_artifacts(
                corrector=artifact_corrector,
                artifacts=artifacts,
                signal=signal,
                time=self.time,
            )

    # save results and end pipeline
    self.signal = signal
    self.artifacts = artifacts

    if callable(correction_method): self.metadata['correction_method'] = correction_method.__name__
    else: self.metadata['correction_method'] = correction_method

    return

extract_trial_data(align_to, center_on, trial_bounds, baseline_bounds=None, event_tolerences={}, trial_normalization='none', check_overlap=True, time_error_threshold=0.1, window_alignment='nearest', invalid_window_policy='drop', event_conflict_logic='first')

Build trial-wise windows, normalize, and store trial data.

Parameters:

• align_to (str) –

  Event label used to align and identify trial, should be one per trial.

• center_on (list[str]) –

  Event labels to center trial windows on. Events should be mutually exclusive (i.e. lever press choice). If no center_on events are present within an identified trial, align_to will be centered on.

• trial_bounds (tuple[float, float]) –

  Trial window bounds relative to center_on events.

• baseline_bounds (tuple[float, float] | None, default: None ) –

  Baseline window bounds relative to align_to event used for per-trial normalizations. Defaults to None.

• event_tolerences (dict[str, tuple[float, float]], default: {} ) –

  Time tolerances for event annotation, relative to align_to. Defaults to {}.

• trial_normalization (Literal['zscore', 'zero', 'mad', 'amp', 'none'] | Callable, default: 'none' ) –

  Normalization method for trial signals based on baselines. A custom function that takes in the positional arguements trial_signals, baseline_signals and returns a 2D np.ndarray of shape (n_trials, n_times) can also be passed. Defaults to 'none'.

• check_overlap (bool, default: True ) –

  Whether to throw an error when multiple center_on events are found in the same trial. Defaults to True.

• time_error_threshold (float, default: 0.1 ) –

  Maximum allowed mean timing error. Defaults to 0.1.

• window_alignment (Literal['nearest', 'interp'], default: 'nearest' ) –

  Stategy for aligning trial times. In 'nearest' mode, events times are rounded to the nearest sampling times, giving a maximum event alignment error of +/- 0.5/frequency. In 'interp' mode, signals are linearly interpolated to an exact event-centered time grid, removing event alignment error but introduction signal interpolation error. Use 'nearest' if event times are already locked to time sampling points.

• invalid_window_policy (Literal['drop', 'error'], default: 'drop' ) –

  Policy for handling trials whose windows extend outside the signal bounds. 'drop' drops the invalid windows while 'error' raises an error if any invalid windows are present.

• event_conflict_logic (Literal['first', 'last', 'mean'], default: 'first' ) –

  Logic for choosing center-on event timestamps if multiple of the same event are present. Defaults to 'first'.

Returns:

• None –

  None

Raises:

• ValueError –

  If baseline_bounds is None and trial_normalization requires baselines.

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

def extract_trial_data(
        self,
        align_to: str,
        center_on: list[str],
        trial_bounds: tuple[float, float],
        baseline_bounds: tuple[float, float] | None = None,
        event_tolerences: dict[str, tuple[float, float]] = {},
        trial_normalization: Literal['zscore', 'zero', 'mad', 'amp', 'none'] | Callable = 'none',
        check_overlap: bool = True,
        time_error_threshold: float = 0.1,
        window_alignment: Literal['nearest', 'interp'] = 'nearest',
        invalid_window_policy: Literal['drop', 'error'] = 'drop',
        event_conflict_logic: Literal['first', 'last', 'mean'] = 'first',
        ) -> None:
    """Build trial-wise windows, normalize, and store trial data.

    Args:
        align_to (str): Event label used to align and identify trial, should be one per trial.

        center_on (list[str]): Event labels to center trial windows on.
            Events should be mutually exclusive (i.e. lever press choice).
            If no ``center_on`` events are present within an identified
            trial, ``align_to`` will be centered on.

        trial_bounds (tuple[float, float]): Trial window bounds relative to
            ``center_on`` events.

        baseline_bounds (tuple[float, float] | None, optional): Baseline
            window bounds relative to ``align_to`` event used for per-trial
            normalizations. Defaults to ``None``.

        event_tolerences (dict[str, tuple[float, float]], optional): Time
            tolerances for event annotation, relative to ``align_to``.
            Defaults to ``{}``.

        trial_normalization (Literal['zscore', 'zero', 'mad', 'amp', 'none'] | Callable, optional):
            Normalization method for trial signals based on baselines.
            A custom function that takes in the positional arguements 
            ``trial_signals``, ``baseline_signals`` and returns a 2D 
            ``np.ndarray`` of shape (n_trials, n_times) can also be passed.
            Defaults to ``'none'``.

        check_overlap (bool, optional): Whether to throw an error when
            multiple ``center_on`` events are found in the same trial.
            Defaults to ``True``.

        time_error_threshold (float, optional): Maximum allowed mean timing
            error. Defaults to ``0.1``.

        window_alignment (Literal['nearest', 'interp'], optional):
            Stategy for aligning trial times. In ``'nearest'`` mode, events times
            are rounded to the nearest sampling times, giving a maximum
            event alignment error of +/- 0.5/frequency. In ``'interp'`` mode,
            signals are linearly interpolated to an exact event-centered time
            grid, removing event alignment error but introduction signal
            interpolation error. Use ``'nearest'`` if event times are already
            locked to time sampling points.

        invalid_window_policy (Literal['drop', 'error'], optional):
            Policy for handling trials whose windows extend outside the signal
            bounds. ``'drop'`` drops the invalid windows while ``'error'`` raises
            an error if any invalid windows are present.

        event_conflict_logic (Literal['first', 'last', 'mean'], optional):
            Logic for choosing center-on event timestamps if multiple of the
            same event are present. Defaults to ``'first'``.

    Returns:
        None

    Raises:
        ValueError: If ``baseline_bounds`` is ``None`` and ``trial_normalization``
            requires baselines.
    """
    # validate inputs
    if align_to not in self.events: 
        raise KeyError(f"align_to '{align_to}' not found in events: {list(self.events)}")

    missing = [lab for lab in center_on if lab not in self.events]
    if missing: raise KeyError(f"center_on labels not found in events: {missing}")

    if baseline_bounds is None:
        calc_baselines = False
        if trial_normalization in ['zscore', 'zero']:
            raise ValueError(f'Baseline bounds have to be specified to for normalization method {trial_normalization}')
    else:
        calc_baselines = True

    # build trials around align_to event
    align_events = self.events[align_to].copy()
    if align_events.size == 0: raise ValueError(f"No '{align_to}' events found.")

    # annotate events based on tolerences
    events_selected = self._annotate_intervals(
        align_to=align_to,
        series=self.time, 
        centers=align_events,
        events=self.events, 
        tolorences=event_tolerences,
        logic=event_conflict_logic
        )

    # find window centers using the selected events
    trial_window_centers = self._find_window_centers(
        center_on=center_on, 
        align_on=align_to, 
        events=events_selected,
        check_overlap=check_overlap,
        )

    # construct trial windows
    trial_windows = self._create_windows(
        signal=self.signal,
        time=self.time,
        events=events_selected,
        centers=trial_window_centers,
        bounds=trial_bounds,
        strategy=window_alignment,
    )

    # do the same for baseline
    if calc_baselines:
        baseline_windows = self._create_windows(
            signal=self.signal,
            time=self.time,
            events=events_selected,
            centers=align_events,
            bounds=baseline_bounds,
            strategy=window_alignment,
        )
    else:
        baseline_windows = WindowResult.empty()

    # handle invalid windows
    self._handle_invalid_windows(
        trial_windows=trial_windows,
        baseline_windows=baseline_windows,
        policy=invalid_window_policy,
    )

    # apply trial-wise normalization method
    processed_trial_signals = self._apply_trial_normalization(
        trial_normalization=trial_normalization, 
        raw_trial_signals=trial_windows.signals, 
        baseline_signals=baseline_windows.signals,
    )

    # save trial data
    self.raw_trial_signal = trial_windows.signals
    self.trial_data = PhotometryData.from_arrays(
        obs=pd.DataFrame(trial_windows.events),
        data=processed_trial_signals,
        time_points=trial_windows.time_grid,
        metadata=self.metadata.copy(),
    )

    # assign trial numbers
    self.trial_data.obs.insert(0, 'trial_num', np.arange(self.trial_data.n_trials, dtype=int) + 1)

    # save baseline data if desired
    if calc_baselines:
        self.baseline_data = PhotometryData.from_arrays(
            obs=pd.DataFrame(baseline_windows.events),
            data=baseline_windows.signals,
            time_points=baseline_windows.time_grid,
            metadata=self.metadata.copy(),
        )

    # check error in times
    time_err = trial_windows.times.std(axis=0).mean()
    if time_err > time_error_threshold:
        raise ValueError(f'Error in rounded times is too high: {time_err} > {time_error_threshold}')

low_frequency_pass_butter(signal, sample_frequency, cutoff_frequency=30.0, order=4, axis=0)

Apply a low-pass Butterworth filter to a signal.

Parameters:

• signal (ndarray) –

  Input signal array.

• sample_frequency (float) –

  Sampling frequency in Hz.

• cutoff_frequency (float, default: 30.0 ) –

  Low-pass cutoff frequency in Hz. Defaults to 30.0.

• order (int, default: 4 ) –

  Butterworth filter order. Defaults to 4.

• axis (int, default: 0 ) –

  Axis of the array to perform a low-pass on. Defaults to 0.

Returns:

• ndarray –

  np.ndarray: Filtered signal.

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

def low_frequency_pass_butter(
        self,
        signal: np.ndarray, 
        sample_frequency: float, 
        cutoff_frequency: float = 30.0, 
        order: int = 4,
        axis: int = 0,
    ) -> np.ndarray:
    """Apply a low-pass Butterworth filter to a signal.

    Args:
        signal (np.ndarray): Input signal array.
        sample_frequency (float): Sampling frequency in Hz.
        cutoff_frequency (float, optional): Low-pass cutoff frequency in Hz.
            Defaults to ``30.0``.
        order (int, optional): Butterworth filter order. Defaults to ``4``.
        axis (int, optional): Axis of the array to perform a low-pass on.
            Defaults to ``0``.

    Returns:
        np.ndarray: Filtered signal.
    """
    normalized_frequency = cutoff_frequency / (sample_frequency / 2)
    sos = butter(order, normalized_frequency, btype='low', output='sos') 
    return sosfiltfilt(sos, signal, axis=axis, padtype='odd', padlen=None)

fit_isosbestic_to_signal(signal, isosbestic, fit_using='IRLS', maxiter=1000, c=None)

Fit the isosbestic channel to the signal.

Parameters:

• signal (ndarray) –

  Filtered signal trace.

• isosbestic (ndarray) –

  Filtered isosbestic trace.

• fit_using (Literal['OLS', 'IRLS', 'IRLS_no_intercept', 'OLS_no_intercept'] | Callable, default: 'IRLS' ) –

  Model used to fit isosbestic. Defaults to 'IRLS'.

• maxiter (int, default: 1000 ) –

  Maximum iterations of IRLS isosbestic fit. Defaults to 1000.

• c (float | None, default: None ) –

  Constant for IRLS fits; smaller values mean more agressive downweighting. 1.4 <= c <= 3 is recommended. Defaults to None.

Returns:

• tuple[ndarray, float, Any] –

  tuple[np.ndarray, float, Any]: Fitted isosbestic, R-squared value, and fit coefficients.

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

def fit_isosbestic_to_signal(
        self, 
        signal: np.ndarray, 
        isosbestic: np.ndarray, 
        fit_using: Literal['OLS', 'IRLS', 'IRLS_no_intercept', 'OLS_no_intercept'] | Callable = 'IRLS',
        maxiter: int = 1000,
        c: float | None = None,
        ) -> tuple[np.ndarray, float, Any]:
    """Fit the isosbestic channel to the signal.

    Args:
        signal (np.ndarray): Filtered signal trace.
        isosbestic (np.ndarray): Filtered isosbestic trace.
        fit_using (Literal['OLS', 'IRLS', 'IRLS_no_intercept', 'OLS_no_intercept'] | Callable, optional): 
            Model used to fit isosbestic. Defaults to ``'IRLS'``.
        maxiter (int, optional): Maximum iterations of IRLS isosbestic fit.
            Defaults to ``1000``.
        c (float | None, optional): Constant for IRLS fits; smaller values
            mean more agressive downweighting. ``1.4 <= c <= 3`` is
            recommended. Defaults to ``None``.

    Returns:
        tuple[np.ndarray, float, Any]: Fitted isosbestic, R-squared value,
            and fit coefficients.
    """
    # fit using specified model
    match fit_using:
        case func if callable(func):
            fitted_iso, params = fit_using(signal, isosbestic)
        case 'OLS':
            fitted_iso, params = OLS_fit(signal, isosbestic, add_intercept=True)
        case 'IRLS':
            fitted_iso, params = IRLS_fit(signal, isosbestic, maxiter=maxiter, c=c, add_intercept=True)
        case 'OLS_no_intercept':
            fitted_iso, params = OLS_fit(signal, isosbestic, add_intercept=False)
        case 'IRLS_no_intercept':
            fitted_iso, params = IRLS_fit(signal, isosbestic, maxiter=maxiter, c=c, add_intercept=False)
        case _:
            raise ValueError(f'{fit_using} fitting method not recognized!')

    # cheack output
    if np.isnan(fitted_iso).any():
        raise ValueError(f'NaN values in fitted isosbestic.')
    r2_val = r2_score(signal, fitted_iso)
    return fitted_iso, r2_val, params

fit_photobleaching_curve(signal, window_dur=5)

Fit a curve with a negative bi-exponential photobleaching model. Uses a soft_l1 least squares to fit to the sliding-window median downsampled signal.

Parameters:

• signal (ndarray) –

  Signal array to fit the photobleaching curve to.

• window_dur (float, default: 5 ) –

  Length of the window in seconds used for sliding-window median downsampling. Defaults to 5.

Returns:

• ndarray –

  tuple[np.ndarray, list[float]]: Fitted curve and fitted parameter

• float –

  values.

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

def fit_photobleaching_curve(
        self,
        signal: np.ndarray, 
        window_dur: float =  5,
        ) -> tuple[np.ndarray, float, list[float]]:
    """Fit a curve with a negative bi-exponential photobleaching model.
    Uses a soft_l1 least squares to fit to the sliding-window median 
    downsampled signal.

    Args:
        signal (np.ndarray): Signal array to fit the photobleaching curve to.
        window_dur (float, optional): Length of the window in seconds used
            for sliding-window median downsampling. Defaults to ``5``.

    Returns:
        tuple[np.ndarray, list[float]]: Fitted curve and fitted parameter
        values.
    """
    window_len = int(window_dur * self.frequency)
    fitted_curve, params = fit_photobleaching(signal, self.time, window=window_len)
    r2_val = r2_score(signal, fitted_curve)

    return fitted_curve, r2_val, params

median_centered_abs_max_check(trial_signal, threshold=0.075)

Check whether a trial signal should be flagged as poor quality

Tests a threshold for minimum MAD (robust Z-score). Threshold needs to be tuned for specific experiments. Generally not recommended.

Parameters:

• trial_signal (ndarray) –

  Trial-by-time signal array.

• threshold (float, default: 0.075 ) –

  Mean absolute-max threshold applied after median centering. Defaults to 0.075.

Returns:

• is_poor_signal ( bool ) –

  True if the signal is classified as poor quality.

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

def median_centered_abs_max_check(
        self, 
        trial_signal: np.ndarray, 
        threshold: float = 0.075
        ):
    """Check whether a trial signal should be flagged as poor quality

    Tests a threshold for minimum MAD (robust Z-score).
    Threshold needs to be tuned for specific experiments.
    Generally not recommended.

    Args:
        trial_signal (np.ndarray): Trial-by-time signal array.
        threshold (float, optional): Mean absolute-max threshold applied
            after median centering. Defaults to ``0.075``.

    Returns:
        is_poor_signal (bool): ``True`` if the signal is classified as poor quality.
    """
    median_centered = trial_signal - np.median(trial_signal, axis=1, keepdims=True)
    abs_max = np.abs(median_centered).max(axis=1)
    is_poor_signal = np.mean(abs_max) < threshold
    return is_poor_signal

dashboard(save=None, downsample=20)

Plot a quick dashboard for the experiment.

Plots the raw, fitted, and, if .preprocess_signal() has been run, processed signal, isosbestic trace, and the fitted photobleaching curve if available.

Parameters:

• save (str | None, default: None ) –

  Path to save the figure. If None, the figure is not saved. Defaults to None.

• downsample (int | optional, default: 20 ) –

  Downsample factor for signals before plotting. Defaults to 20.

Returns:

• None –

  None

Source code in pyFiberPhotometry/core/PhotometryExperiment.py

def dashboard(self, save: str | None = None, downsample: int = 20) -> None:
    """Plot a quick dashboard for the experiment.

    Plots the raw, fitted, and, if ``.preprocess_signal()`` has been run, 
    processed signal, isosbestic trace, and the fitted photobleaching curve 
    if available.

    Args:
        save (str | None, optional): Path to save the figure. If ``None``,
            the figure is not saved. Defaults to ``None``.
        downsample (int | optional): Downsample factor for signals before
            plotting. Defaults to 20.

    Returns:
        None
    """
    # determine is self.process_singal has been run
    if self.fitted_ref is not None:
        fig = self.dashboard_full()
    else:
        fig = self.dashboard_raw()

    if save is not None:
        plt.savefig(save, bbox_inches='tight')