Harmonic connectivity#

biotuner.harmonic_connectivity — cross-frequency / cross-channel coupling pipeline.

Module type: Object

Computes phase-amplitude coupling (PAC), cross-frequency coupling (CFC), endogenous intermodulation, and harmonicity matrices between biosignal channels or frequency bands.

class harmonic_connectivity(sf=None, data=None, peaks_function='EMD', precision=0.1, n_harm=10, harm_function='mult', min_freq=2, max_freq=80, n_peaks=5)[source]#

Bases: object

Class used to compute harmonicity metrics between pairs of sensors.

Parameters:
  • sf (int) – sampling frequency (in Hz)

  • data (2D array (elec, numDataPoints)) – Electrodes x Time series to analyse.

  • peaks_function (str, default=’EMD’) – See compute_biotuner class for details.

  • precision (float, default=0.1) – Precision of the peaks (in Hz) When HH1D_max is used, bins are in log scale.

  • n_harm (int, default=10) – Set the number of harmonics to compute in harmonic_fit function

  • harm_function (str, default=’mult’) – Computes harmonics from iterative multiplication (x, 2x, 3x, …nx) or division (x, x/2, x/3, …x/n). Choose between ‘mult’ and ‘div’

  • min_freq (float, default=2) – Minimum frequency (in Hz) to consider for peak extraction.

  • max_freq (float, default=80) – Maximum frequency (in Hz) to consider for peak extraction.

  • n_peaks (int, default=5) – Number of peaks to extract per frequency band.

compute_harm_connectivity(metric='harmsim', delta_lim=20, save=False, savename='_', graph=True, FREQ_BANDS=None, max_denom_rrci=16)[source]#

Computes the harmonic connectivity matrix between electrodes.

Parameters:
  • metric (str, optional) – The metric to use for computing harmonic connectivity. Default is ‘harmsim’.

    Possible values are:

    • ‘harmsim’:

      computes the harmonic similarity between each pair of peaks from the two electrodes. It calculates the ratio between each pair of peaks and computes the mean harmonic similarity.

    • ‘euler’:

      computes the Euler’s Gradus Suavitatis on the concatenated peaks of the two electrodes.

    • ‘harm_fit’:

      computes the number of common harmonics between each pair of peaks from the two electrodes. It evaluates the harmonic fit between each peak pair and counts the number of common harmonics.

    • ‘subharm_tension’:

      computes the tension between subharmonics of two electrodes. It evaluates the tension between subharmonics of the two electrodes by comparing the subharmonics and their ratios.

    • ‘RRCi’:

      computes the Rhythmic Ratio Coupling with Imaginary Component (RRCi) metric between each pair of peaks from the two electrodes, using a bandwidth of 2 Hz and a max_denom of 16. This metric calculates the imaginary part of the complex phase differences between two filtered signals.

    • ‘wPLI_crossfreq’:

      computes the weighted Phase Lag Index (wPLI) for cross-frequency coupling between each pair of peaks from the two electrodes. The wPLI measures the phase synchronization between two signals, with a value close to 0 indicating no synchronization and a value close to 1 indicating perfect synchronization.

    • ‘wPLI_multiband’:

      computes the weighted Phase Lag Index (wPLI) for multiple frequency bands between the two electrodes. It calculates wPLI for each frequency band and returns an array of wPLI values for the defined frequency bands.

  • delta_lim (int, optional) – The delta limit for the subharmonic tension metric. Default is 20.

  • save (bool, optional) – Whether to save the connectivity matrix. Default is False.

  • savename (str, optional) – The name to use when saving the connectivity matrix. Default is ‘_’.

  • graph (bool, optional) – Whether to display a heatmap of the connectivity matrix. Default is True.

  • FREQ_BANDS (list, optional) – The frequency bands to use for the computation of the wPLI_multiband metric. Default is None. If None, the following frequency bands will be used: [2, 3.55], [3.55, 7.15], [7.15, 14.3], [14.3, 28.55], [28.55, 49.4].

  • max_denom_rrci (int, optional) – The maximum denominator to use for the computation of the RRCi metric. Default is 16.

Returns:

conn_matrix (numpy.ndarray) – The harmonic connectivity matrix between electrodes.

compute_IMF_correlation(nIMFs=5, freq_range=(1, 60), precision=0.5, delta_lim=50)[source]#

Compute the correlation, coherence, and peak frequency between each pair of IMFs for each pair of electrodes.

Returns:

df (pd.DataFrame) – DataFrame with columns: elec1, elec2, imf1, imf2, pearson, coherence, and peak_freq.

compute_time_resolved_harm_connectivity(sf, nIMFs, metric='harmsim', delta_lim=50)[source]#

Computes the time-resolved harmonic connectivity matrix between electrodes, which is a harmonic connectivity matrix for each intrinsic mode function (IMF), and each time point.

Parameters:
  • data (numpy.ndarray) – Input data with shape (num_electrodes, numDataPoints)

  • sf (int) – Sampling frequency

  • nIMFs (int) – Number of intrinsic mode functions (IMFs) to consider.

  • metric (str, optional) – The metric to use for computing harmonic connectivity. Default is ‘harmsim’.

  • delta_lim (int, optional) – The delta limit for the subharmonic tension metric. Default is 20.

Returns:

connectivity_matrices (numpy.ndarray) – Time-resolved harmonic connectivity matrices with shape (IMFs, numDataPoints, electrodes, electrodes).

Notes

!!! This method is very computationally expensive and can take a long time to run. !!!

transitional_connectivity(data=None, sf=None, mode='win_overlap', overlap=10, delta_lim=20, graph=False, n_trans_harm=3)[source]#

This function calculates the transitional connectivity among electrodes, utilizing concepts of transitional harmony and temporal correlation. It does so by employing a two-step approach:

  • Transitional Harmony Calculation: For every electrode, the transitional harmony is determined first.

    This measure reflects the patterns or sequences of signal transitions over time that each electrode experiences.

  • Temporal Correlation with FDR Correction: After obtaining transitional harmonies, the function evaluates the temporal correlation

    between these harmonies for each pair of electrodes. Temporal correlation quantifies how similar the timing and pattern of transitional harmonies are between each pair of electrodes.

The entire process takes into account multiple comparisons, utilizing False Discovery Rate (FDR) correction to reduce the likelihood of false positives. This correction is crucial in maintaining the validity and robustness of the results, especially when dealing with a large number of comparisons, as is the case with numerous electrodes.

Parameters:
  • data (numpy.ndarray) – Multichannel EEG data with shape (n_electrodes, n_timepoints).

  • sf (float) – Sampling frequency of the EEG data in Hz.

  • mode (str, optional, default=’win_overlap’) – The mode to compute the transitional harmony. ‘win_overlap’ computes the transitional harmony using a sliding window with overlap.

  • overlap (int, optional, default=10) – The percentage of overlap between consecutive windows when computing the transitional harmony. Default is 10.

  • delta_lim (int, optional, default=20) – The minimum delta value for the computation of subharmonic tension, in ms. Default is 20.

  • graph (bool, optional, default=False) – If True, it will plot the graph of the transitional harmony. Default is False.

  • n_trans_harm (int, optional, default=3) – Number of transitional harmonics to compute. Default is 3.

Returns:

  • conn_mat (numpy.ndarray) – Connectivity matrix of shape (n_electrodes, n_electrodes) representing the transitional connectivity between electrodes.

  • pval_mat (numpy.ndarray) – P-value matrix of shape (n_electrodes, n_electrodes) with FDR-corrected p-values for the computed connectivity values.

plot_conn_matrix(conn_matrix=None, node_names=None, n_lines=50)[source]#

Plots a connectivity matrix in a circle plot.

Parameters:
  • conn_matrix (ndarray, optional) – The connectivity matrix to plot. If None, uses the object’s own attribute conn_matrix. If conn_matrix is still None, a ValueError will be raised.

  • node_names (list, optional) – The labels for the nodes. If None, a range object will be converted to a list of string values for node names.

  • n_lines (int, default=50) – The number of lines to draw in the plot. Default is 50.

Returns:

fig (matplotlib.figure.Figure) – A matplotlib figure object.

Raises:

ValueError – Raised when no connectivity matrix is found.

Notes

This method uses the function plot_connectivity_circle to plot the connectivity matrix.

compute_harmonic_spectrum_connectivity(sf=None, data=None, precision=0.5, fmin=None, fmax=None, noverlap=1, power_law_remove=False, n_peaks=5, metric='harmsim', n_harms=10, delta_lim=0.1, min_notes=2, plot=False, smoothness=1, smoothness_harm=1, phase_mode=None, save_fig=False, savename='harmonic_spectrum_connectivity.png')[source]#

Computes the harmonic spectrum connectivity between pairs of electrodes. For more details on the harmonic spectrum, see the function compute_cross_spectrum_harmonicity().

Parameters:
  • sf (float, optional) – Sampling frequency of the data. If not provided, uses the object’s own attribute sf.

  • data (ndarray, optional) – 2D data array of shape (n_electrodes, n_datapoints). If not provided, uses the object’s own attribute data.

  • precision (float, default=0.5) – Precision of the frequency axis of the cross-spectrum.

  • fmin (float, optional) – Minimum frequency for computation. If not provided, uses the smallest possible frequency (i.e., zero).

  • fmax (float, optional) – Maximum frequency for computation. If not provided, uses the Nyquist frequency.

  • noverlap (int, default=1) – Argument passed to the scipy.signal.stft() function.

  • power_law_remove (bool, default=False) – If True, removes the power-law noise from the cross-spectrum.

  • n_peaks (int, default=5) – Number of peaks to derive from the harmonic spectrum.

  • metric (str, default=’harmsim’) – Name of the metric to be used in the computation of harmonicity. Choose between:

    • ‘harmsim’

    • ‘subharm_tension’

  • n_harms (int, default=10) – Number of harmonics to consider in the computation.

  • delta_lim (float, default=0.1) – Threshold for the delta limit used in subharmonic tension computation.

  • min_notes (int, default=2) – Minimum number of notes required for a harmonic pattern to be considered valid, when subharmonic tension is used as the metric.

  • plot (bool, default=False) – If True, plots the results.

  • smoothness (int, default=1) – Smoothness factor of the power spectrum. When smoothness=1, no smoothing is applied.

  • smoothness_harm (int, default=1) – Smoothness factor of the harmonic spectrum. When smoothness_harm=1, no smoothing is applied.

  • phase_mode (str, optional) – If set to ‘weighted’, the phase coupling is weighted by the power of associated frequencies.

  • save_fig (bool, default=False) – If True, saves the resulting plot as a .png file.

  • savename (str, default=’harmonic_spectrum_connectivity.png’) – Name of the .png file to save if save_fig is True.

Returns:

output (pandas.DataFrame) – DataFrame containing the results of the harmonic spectrum connectivity computation.

Notes

The harmonic spectrum connectivity is computed between each pair of electrodes. This is achieved by computing the cross-spectrum harmonicity for each pair of electrodes, and storing the results in a DataFrame. The DataFrame includes indices of the electrode pairs, alongside the results of the harmonicity computation.

get_harm_spectrum_metric_matrix(metric)[source]#

Method to retrieve a matrix of harmonic spectrum metric values between pairs of electrodes.

Parameters:

metric (str) – The specific metric from the harmonic spectrum connectivity data to create the matrix from.

Returns:

matrix (pandas.DataFrame or None) – A pivot table with ‘elec1’ and ‘elec2’ as indices and the provided ‘metric’ as values. If ‘harmonic_spectrum_connectivity’ is None or ‘metric’ is not found in ‘harmonic_spectrum_connectivity’, this method will print an error message and return None.

Notes

This method checks if ‘harmonic_spectrum_connectivity’ exists in the object. If it does, it further checks if ‘metric’ exists in the DataFrame’s columns. If both checks pass, a pivot table is created using ‘elec1’ and ‘elec2’ as indices and ‘metric’ as values.

plot_IMF_correlation_matrix(df, elec1, elec2, variable='pearson', savepath=None)[source]#

Plot a heatmap of correlations between IMFs for specified electrodes.

Parameters:
  • df (pd.DataFrame) – DataFrame containing the correlation data with columns: elec1, elec2, imf1, imf2, and pearson.

  • elec1 (int) – The first electrode of interest.

  • elec2 (int) – The second electrode of interest.

Returns:

None

wPLI_crossfreq(signal1, signal2, peak1, peak2, sf)[source]#

Computes the weighted phase lag index (wPLI) between two signals in specified frequency bands, centered around provided peak frequencies.

Parameters:
  • signal1 (ndarray) – First signal in time series.

  • signal2 (ndarray) – Second signal in time series.

  • peak1 (float) – Center of the frequency band for the first signal.

  • peak2 (float) – Center of the frequency band for the second signal.

  • sf (float) – Sampling frequency.

Returns:

wPLI (float) – Weighted phase lag index between two signals in specified frequency bands.

Examples

>>> signal1 = np.random.normal(0, 1, 5000)
>>> signal2 = np.random.normal(0, 1, 5000)
>>> peak1 = 10.0
>>> peak2 = 20.0
>>> sf = 100.0
>>> wPLI = wPLI_crossfreq(signal1, signal2, peak1, peak2, sf)
wPLI_multiband(signal1, signal2, freq_bands, sf)[source]#

Computes the weighted phase lag index (wPLI) between two signals for multiple frequency bands.

Parameters:
  • signal1 (ndarray) – First signal in time series.

  • signal2 (ndarray) – Second signal in time series.

  • freq_bands (list of tuple) – List of frequency bands. Each band is represented as a tuple (lowcut, highcut).

  • sf (float) – Sampling frequency.

Returns:

wPLI_values (list of float) – List of wPLI values for each frequency band.

Examples

>>> signal1 = np.random.normal(0, 1, 5000)
>>> signal2 = np.random.normal(0, 1, 5000)
>>> freq_bands = [(8, 12), (13, 30), (30, 70)]
>>> sf = 100.0
>>> wPLI_values = wPLI_multiband(signal1, signal2, freq_bands, sf)
n_m_phase_locking(signal1, signal2, n, m, sf)[source]#

Calculate n:m phase locking between two signals.

Parameters:
  • signal1 – First time series.

  • signal2 – Second time series.

  • n – Multiplicative factor for the first signal.

  • m – Multiplicative factor for the second signal.

  • sf – Sampling frequency of the signals.

Returns:

Mean resultant length (Rn:m) indicating the phase locking value.

cross_frequency_rrci(signal1, signal2, sfreq, freq_peak1, freq_peak2, bandwidth=1, max_denom=50)[source]#

Computes the Rhythmic Ratio Coupling Index (RRCI) between two signals for cross-frequency.

The function first calculates the rhythmic ratio between two peak frequencies. It then filters the two input signals around a frequency band centered on their respective peak frequencies. Finally, it calculates the rhythmic ratio coupling index (RRCI) for these filtered signals. The RRCI is a measure of how much the rhythms of the two signals, in terms of their phase information, are coupled across different frequencies. In other words, it provides a measure of phase-to-phase coupling across these frequencies.

Parameters:
  • signal1 (ndarray) – First signal in time series.

  • signal2 (ndarray) – Second signal in time series.

  • sfreq (float) – Sampling frequency.

  • freq_peak1 (float) – Peak frequency for the first signal.

  • freq_peak2 (float) – Peak frequency for the second signal.

  • bandwidth (float, default=1) – Frequency bandwidth for filtering the signals.

  • max_denom (int, default=50) – The maximum denominator for the rhythmic ratio.

Returns:

rrci (float) – Rhythmic Ratio Coupling Index between the two signals.

Examples

>>> signal1 = np.random.normal(0, 1, 5000)
>>> signal2 = np.random.normal(0, 1, 5000)
>>> sfreq = 100.0
>>> rrci = cross_frequency_rrci(signal1, signal2, sfreq, 10, 20, 2, 4)
compute_rhythmic_ratio(freq1, freq2, max_denom)[source]#
rhythmic_ratio_coupling_imaginary(signal1, signal2, rhythmic_ratio, sfreq, freq_band1, freq_band2)[source]#

Computes the Imaginary part of Rhythmic Ratio Coupling between two signals. From the paper : “On cross-frequency phase-phase coupling between theta and gamma oscillations in the hippocampus”

Parameters:
  • signal1 (ndarray) – First signal in time series.

  • signal2 (ndarray) – Second signal in time series.

  • rhythmic_ratio (tuple) – Rhythmic ratio (numerator, denominator).

  • sfreq (float) – Sampling frequency.

  • freq_band1 (tuple) – Frequency band for the first signal (lowcut, highcut).

  • freq_band2 (tuple) – Frequency band for the second signal (lowcut, highcut).

Returns:

imaginary_part (float) – Imaginary part of the Rhythmic Ratio Coupling.

HilbertHuang1D_nopeaks(data, sf, graph=False, nIMFs=5, min_freq=1, max_freq=45, precision=0.5, bin_spread='log', smooth_sigma=None)[source]#

The Hilbert-Huang transform provides a description of how the energy or power within a signal is distributed across frequency. The distributions are based on the instantaneous frequency and amplitude of a signal.

Parameters:
  • data (array (numDataPoints,)) – Single time series.

  • sf (int) – Sampling frequency.

  • graph (bool, default=False) – Defines if a graph is generated.

  • nIMFs (int, default=5) – Number of intrinsic mode functions (IMFs) to keep when Empirical Mode Decomposition (EMD) is computed.

  • min_freq (float, default=1) – Minimum frequency to consider.

  • max_freq (float, default=80) – Maximum frequency to consider.

  • precision (float, default=0.1) – Value in Hertz corresponding to the minimal step between two frequency bins.

  • bin_spread (str, default=’log’) –

    • ‘linear’

    • ‘log’

  • smooth_sigma (float, default=None) – Sigma value for gaussian smoothing.

Returns:

  • IF (array (numDataPoints,nIMFs)) – instantaneous frequencies associated with each IMF.

  • IP (array (numDataPoints,nIMFs)) – instantaneous power associated with each IMF.

  • IA (array (numDataPoints,nIMFs)) – instantaneous amplitude associated with each IMF.

  • spec (array (nIMFs, nbins)) – Power associated with all bins for each IMF

  • bins (array (nIMFs, nbins)) – Frequency bins for each IMF

EMD_time_resolved_harmonicity(time_series1, time_series2, sf, nIMFs=5, method='harmsim')[source]#

Computes the harmonicity between the instantaneous frequencies (IF) for each point in time between all pairs of corresponding intrinsic mode functions (IMFs).

Parameters:
  • time_series1 (array (numDataPoints,)) – First time series.

  • time_series2 (array (numDataPoints,)) – Second time series.

  • sf (int) – Sampling frequency.

  • nIMFs (int, default=5) – Number of intrinsic mode functions (IMFs) to consider.

Returns:

harmonicity (array (numDataPoints, nIMFs)) – Harmonicity values for each pair of corresponding IMFs.

temporal_correlation_fdr(data)[source]#

Compute the temporal correlation for each pair of electrodes and output a connectivity matrix and a matrix of FDR-corrected p-values. This function is used in the computation of the transitional harmony connectivity.

Parameters:

data (array-like) – An array of shape (electrodes, samples) containing the electrode recordings.

Returns:

  • connectivity_matrix (ndarray) – A connectivity matrix of shape (electrodes, electrodes) with the temporal correlation for each pair of electrodes.

  • fdr_corrected_pvals (ndarray) – A matrix of FDR-corrected p-values of shape (electrodes, electrodes).

Notes

This function calculates the temporal correlation for each pair of electrodes, creating a connectivity matrix. Simultaneously, it calculates a matrix of p-values and corrects for multiple comparisons using the False Discovery Rate (FDR) method.

compute_cross_spectrum_harmonicity(signal1, signal2, precision_hz, fmin=None, fmax=None, noverlap=1, fs=44100, power_law_remove=False, n_peaks=5, metric='harmsim', n_harms=10, delta_lim=0.1, min_notes=2, plot=False, smoothness=1, smoothness_harm=1, phase_mode=None, save_fig=False, save_name='harmonic_spectrum_connectivity.png')[source]#

Compute the cross-spectrum harmonicity between two signals. This function is useful for analyzing the interaction between frequency components of two different signals.

Parameters:
  • signal1, signal2 (array_like) – Input signals.

  • precision_hz (float) – The precision in Hz for the short time Fourier transform.

  • fmin, fmax (float or None) – Minimum and maximum frequency for computing the power spectral density.

  • noverlap (int) – Number of points to overlap between segments for computing the power spectral density.

  • fs (int) – The sampling frequency of the signals.

  • power_law_remove (bool) – If True, power-law noise is removed from the power spectral densities of the signals.

  • n_peaks (int) – The number of peaks to identify in the spectra.

  • metric ({“harmsim”, “subharm_tension”}) – The metric to use for computing the harmonicity between the signals.

  • n_harms (int) – Number of harmonics to consider for computing subharmonic tension.

  • delta_lim (float) – Delta limit for computing subharmonic tension.

  • min_notes (int) – Minimum number of notes for computing subharmonic tension.

  • plot (bool) – If True, the function will plot the cross harmonic and phase coupling spectrum.

  • smoothness (int) – Smoothness of the Fourier transform. Higher values result in a smoother transform.

  • smoothness_harm (int) – Smoothness of the harmonic spectrum. Higher values result in a smoother spectrum.

  • phase_mode ({None, “weighted”}) – If “weighted”, the phase coupling is computed as a weighted sum. If None, it is computed as a simple average.

  • save_fig (bool) – If True, the plot is saved as a png file.

  • save_name (str) – Name of the png file if save_fig is True.

Returns:

DataFrame – A pandas DataFrame with columns representing various computed metrics of the signals, such as spectral flatness, spectral entropy, etc. The DataFrame also contains columns for peak frequencies in the spectra and their harmonic similarities.