Aluode/PerceptionLabPortable

1# Authors: The MNE-Python contributors.2# License: BSD-3-Clause3# Copyright the MNE-Python contributors.4 5import numpy as np6 7from .._fiff.pick import _picks_to_idx8from ..evoked import Evoked9from ..utils import _ensure_int, _validate_type, verbose10 11 12def _temp_proj(ref_2, ref_1, raw_data, n_proj=6):13    # Orthonormalize gradiometer and magnetometer data by a QR decomposition14    ref_1_orth = np.linalg.qr(ref_1.T)[0]15    ref_2_orth = np.linalg.qr(ref_2.T)[0]16 17    # Calculate cross-correlation18    cross_corr = np.dot(ref_1_orth.T, ref_2_orth)19 20    # Channel weights for common temporal subspace by SVD of cross-correlation21    ref_1_ch_weights, _, _ = np.linalg.svd(cross_corr)22 23    # Get temporal signals from channel weights24    proj_mat = ref_1_orth @ ref_1_ch_weights25 26    # Project out common subspace27    filtered_data = raw_data28    proj_vec = proj_mat[:, :n_proj]29    weights = filtered_data @ proj_vec30    filtered_data -= weights @ proj_vec.T31 32 33@verbose34def cortical_signal_suppression(35    evoked, picks=None, mag_picks=None, grad_picks=None, n_proj=6, *, verbose=None36):37    """Apply cortical signal suppression (CSS) to evoked data.38 39    Parameters40    ----------41    evoked : instance of Evoked42        The evoked object to use for CSS. Must contain magnetometer,43        gradiometer, and EEG channels.44    %(picks_good_data)s45    mag_picks : array-like of int46        Array of the first set of channel indices that will be used to find47        the common temporal subspace. If None (default), all magnetometers will48        be used.49    grad_picks : array-like of int50        Array of the second set of channel indices that will be used to find51        the common temporal subspace. If None (default), all gradiometers will52        be used.53    n_proj : int54        The number of projection vectors.55    %(verbose)s56 57    Returns58    -------59    evoked_subcortical : instance of Evoked60        The evoked object with contributions from the ``mag_picks`` and ``grad_picks``61        channels removed from the ``picks`` channels.62 63    Notes64    -----65    This method removes the common signal subspace between two sets of66    channels (``mag_picks`` and ``grad_picks``) from a set of channels67    (``picks``) via a temporal projection using ``n_proj`` number of68    projection vectors. In the reference publication :footcite:`Samuelsson2019`,69    the joint subspace between magnetometers and gradiometers is used to70    suppress the cortical signal in the EEG data. In principle, other71    combinations of sensor types (or channels) could be used to suppress72    signals from other sources.73 74    References75    ----------76    .. footbibliography::77    """78    _validate_type(evoked, Evoked, "evoked")79    n_proj = _ensure_int(n_proj, "n_proj")80    picks = _picks_to_idx(evoked.info, picks, none="data", exclude="bads")81    mag_picks = _picks_to_idx(evoked.info, mag_picks, none="mag", exclude="bads")82    grad_picks = _picks_to_idx(evoked.info, grad_picks, none="grad", exclude="bads")83    evoked_subcortical = evoked.copy()84 85    # Get data86    all_data = evoked.data87    mag_data = all_data[mag_picks]88    grad_data = all_data[grad_picks]89 90    # Process data with temporal projection algorithm91    data = all_data[picks]92    _temp_proj(mag_data, grad_data, data, n_proj=n_proj)93    evoked_subcortical.data[picks, :] = data94 95    return evoked_subcortical96