A graph fourier domain 256 lead electroencephalogram time offset synchronization calibration method

Abstract

The application discloses a kind of graph fourier domain 256 guide brain electricity time offset synchronous calibration method, belong to biomedical signal processing technical field.The present application constructs 256 lead spatial topology graph and its graph Laplace matrix for the problem that the source positioning accuracy of high-density brain electricity system between each channel sub-millisecond level hardware time offset leads to decline, obtains graph fourier transform base by characteristic decomposition;The multi-channel signal of each time is transformed to graph frequency domain, and the phase difference of each channel relative to virtual reference is calculated;Linear relationship of phase difference and graph frequency is fitted using weighted least squares regression, and the sub-sampling level time offset of each channel is estimated;Finally, signal resampling calibration is realized using fractional delay interpolation.The present application does not need external training data, and the calibration accuracy is better than 0.02 times sampling interval, and offline calculation graph base is processed online with less than 5ms, which can significantly improve the synchronization quality and source positioning accuracy of 256 guide brain electricity.

Description

Technical Field

[0001] This invention relates to the field of biomedical signal processing technology, and in particular to a method for synchronous calibration of EEG time shift in the 256-channel graph Fourier domain. Background Technology

[0002] Electroencephalography (EEG), as a non-invasive brain imaging technique, is widely used in fields such as epileptic focus localization, cognitive science, and brain-computer interfaces. High-density EEG (256 leads or more) has become a cutting-edge tool in neuroimaging research due to its ability to significantly increase spatial sampling density. However, high-density EEG systems face an inherent engineering challenge: although all channels share the same master clock, differences in the delay of the front-end amplifier group, the aperture jitter of the analog-to-digital converter, and the inconsistent trace lengths on the printed circuit board result in sub-millisecond (approximately 0.1ms to 1ms) relative time offsets between the actual acquired signals. For a 256-lead system, with a lead spacing of only about 2-3cm, the time offset is equivalent to a source localization error of tens of millimeters, severely undermining the spatial advantages of high-density EEG.

[0003] Existing solutions mainly include hardware improvements (high cost and difficult to modify) and software post-processing methods (cross-correlation methods rely on stationarity assumptions, phase-locked loop methods are sensitive to noise, and deep learning methods require extensive annotation and have weak interpretability). Therefore, there is an urgent need for a synchronous calibration method that does not rely on hardware modifications, has clear physical meaning, and fully utilizes the high-density topological information of the 256 leads. Summary of the Invention

[0004] The purpose of this invention is to provide a synchronous calibration method for time offset of 256-channel EEG in the graphical Fourier domain, which solves the problem that it is difficult to accurately and robustly estimate the sub-millisecond time offset of multi-channel 256-channel EEG systems in the prior art.

[0005] To achieve the above objectives, the present invention provides a method for synchronous calibration of EEG time shift in the 256-channel graph Fourier domain, comprising the following steps: Step S1: Simultaneously acquire raw EEG signals using 256 acquisition channels and preprocess them to obtain a multi-channel time-domain signal matrix; Step S2: Construct the graph Laplace matrix based on the lead spatial location ; Step S3: For Perform eigenvalue decomposition The graph Fourier transform basis is obtained. sum graph frequency ; Step S4: For the graph signal at each time point Perform graph Fourier transform The spectrum of the graph is obtained; Step S5: Calculate the phase difference of each channel relative to the virtual reference. ; Step S6: For each channel, perform a weighted least squares regression of the phase difference against the graph frequency to fit a straight line. Solve for the time offset The weight Take the spectrum amplitude energy; Step S7: According to The original signal is calibrated by fractional delay interpolation, and the synchronized 256-channel EEG signal is output.

[0006] Preferably, in step S2, the adjacency matrix is ​​constructed using a threshold Gaussian kernel function: (When the distance of the geodesic line) The Graph Laplace matrix is .

[0007] Preferably, in step S4, before performing the graph Fourier transform, the graph signal is first weighted using a graph domain Hanning window, and the window function is... .

[0008] Preferably, in step S5, the virtual reference signal is taken as the complex average of the spectral coefficients of all channels. Phase difference .

[0009] Preferably, in step S6, the Huber loss function is used for iterative reweighting to suppress outliers, and convergence is achieved after 3 iterations.

[0010] Preferably, in step S7, the fractional delay filter uses fourth-order Lagrange interpolation, with interpolation coefficients of... Output signal .

[0011] Preferably, the preprocessing in step S1 includes mean removal, energy normalization, and electrooculography / electromyography artifact removal based on independent component analysis.

[0012] Preferably, the method further includes step S8, which outputs the synchronized signal and the estimated time offset to the neuroimaging system for lead field delay compensation in source localization.

[0013] Preferably, the entire method does not require external training data during calibration, and the time offset estimation accuracy is better than 0.02 times the sampling interval.

[0014] Therefore, the present invention, employing the above-described structure, provides a graphical Fourier domain 256-channel EEG time shift synchronization calibration method, which has the following beneficial effects: (1) This invention applies the phase linearity property of the graph Fourier transform to the synchronous calibration of EEG multichannels. By utilizing the high-density topological information of 256 leads, the time-domain delay problem is transformed into a graph frequency domain phase regression problem. The physical meaning is clear, and the time offset has a closed analytical solution. .

[0015] (2) This invention fully utilizes the phase information of high-energy components of the graph spectrum through weighted least squares and Huber loss, resulting in high estimation accuracy and robustness.

[0016] (3) The eigenvalue decomposition of the graph Laplacian matrix of the present invention only needs to be calculated offline once, and the online calibration only involves matrix multiplication and linear regression. The calibration time is less than 5ms, and the computation efficiency is high.

[0017] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. Attached Figure Description

[0018] Figure 1 This is a schematic diagram of the flowchart of a method for synchronizing the temporal offset of 256-channel EEG in the Fourier domain according to the present invention. Figure 2 This invention relates to a 256-lead scalp geodesic distance and topology map for a 256-lead EEG time-shift synchronization calibration method in the graphical Fourier domain. Detailed Implementation

[0019] The technical solution of the present invention will be further described below with reference to the accompanying drawings and embodiments.

[0020] Unless otherwise defined, the technical or scientific terms used in this invention shall have the ordinary meaning understood by one of ordinary skill in the art to which this invention pertains. The terms "first," "second," and similar terms used in this invention do not indicate any order, quantity, or importance, but are merely used to distinguish different components. Terms such as "comprising" or "including" mean that the element or object preceding the word encompasses the elements or objects listed following the word and their equivalents, without excluding other elements or objects. Terms such as "connected" or "linked" are not limited to physical or mechanical connections, but can include electrical connections, whether direct or indirect. Terms such as "upper," "lower," "left," and "right" are used only to indicate relative positional relationships; when the absolute position of the described object changes, the relative positional relationship may also change accordingly.

[0021] Example like Figure 1 As shown, this invention provides a method for synchronizing and calibrating the time shift of 256-channel EEG in the Fourier domain, as detailed below: Step S1: Data Acquisition and Preprocessing Raw EEG signals were acquired simultaneously using 256 acquisition channels. The sampling rate was set at... Sampling interval The signal of each channel is subjected to mean-reduction processing: ; Then, energy normalization is performed: ; To reduce interference from electrooculography (EOG) and electromyography (EMG) artifacts, independent component analysis (ICA) was used to remove the three components most correlated with the EOG template before reconstructing the signal. The final preprocessed signal matrix was obtained. .

[0022] Step S2: Constructing the lead space map like Figure 2 As shown, the geodesic distance between any two leads is calculated based on the three-dimensional coordinates of each lead provided by the international 10-10 system. (Using the Vincenty formula for ellipsoidal surfaces). The adjacency matrix is ​​constructed using a threshold Gaussian kernel function. : ; Among them, the cutoff distance Take 1.2 times the median distance of all geodesic lines. Calculate the degree matrix. (on the diagonal) Then construct the normalized graphical Laplace matrix: ; Step S3: Solving for the graphical Fourier transform basis Perform eigenvalue decomposition on the graph Laplacian matrix: ; in The eigenvalues ​​are arranged in ascending order. ; This is the eigenvector matrix. Eigenvectors Constructing the Fourier transform basis of the graph, eigenvalues Defined as graph frequency. This decomposition only needs to be calculated offline once.

[0023] Step S4: Graphical Fourier Transform For each sampling time (Corresponding to one sampling point), extract the signal values ​​of all channels at that moment to form a graph signal vector. To suppress spectral leakage, first... Perform Hanning window weighting on the graph domain: ; in For nodes The degree, To maximize the value, then perform a graphical Fourier transform: ; Obtain the graph spectrum vector , its first Each component Corresponding frequency The spectral coefficients below.

[0024] Step S5: Phase difference calculation Construct a virtual reference signal and take the complex average of the spectral coefficients of all channels: ; For each channel i, its phase difference relative to the reference is: ; in Returns the principal argument value (range) of a complex number. To avoid phase entanglement, when the absolute value of the phase difference between adjacent frequency bands exceeds... At that time, perform untangling: Select integer This ensures that the phase difference between adjacent points is continuous.

[0025] Step S6: Weighted least squares regression to estimate time shift For each channel Obtain data points According to graph signal processing theory, time delay (In units of sampling interval) and phase difference satisfy a linear relationship: ; We use weighted least squares to fit a straight line passing through the origin, with the weights taken as the amplitude energy of the corresponding frequency components in the graph: ; The regression solution is: ; To improve robustness, Huber loss is used for iterative reweighting. The Huber loss function is defined as follows: ; Among residuals threshold (Twice the absolute deviation of the median). The final result is obtained after iteratively updating the weights three times. .

[0026] Step S7: Subsampling level resampling calibration For channels Original discrete signal Delayed reconstruction required Subsequent signals: ; Implemented using a fractional delay filter with 4th-order Lagrange interpolation: ; The interpolation coefficients are: After performing the above operations on all 256 channels, a time-synchronized EEG signal matrix is ​​obtained. .

[0027] Step S8 (optional): Output and Application synchronized signal and the estimated time offset for each channel Output to the neuroimaging reconstruction system. In the source localization calculation, time delay compensation is performed on the lead field matrix of the forward model: .

[0028] Working principle: Based on graph signal processing theory: If a signal defined on a graph has a small time-domain shift... The phase of the Fourier transform spectrum of the graph changes linearly with the frequency of the graph, and the slope is... By constructing a Laplace graph of geodesic distance, eigenvalue decomposition yields the graph Fourier basis. Phase difference is extracted after signal transformation at each time step, and weighted regression is used to estimate the offset. Finally, fractional delay interpolation achieves subsampling-level calibration. This method requires no stationarity assumption, is independent of training data, and possesses analytical solutions and high robustness.

[0029] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit them. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the technical solutions of the present invention, and these modifications or equivalent substitutions cannot cause the modified technical solutions to deviate from the spirit and scope of the technical solutions of the present invention.

Claims

1. A method for synchronous calibration of time shift in 256-channel EEG in the graphical Fourier domain, characterized in that, Includes the following steps: Step S1: Simultaneously acquire raw EEG signals using 256 acquisition channels and preprocess them to obtain a multi-channel time-domain signal matrix; Step S2: Construct the graph Laplace matrix based on the lead spatial location ; Step S3: For Perform eigenvalue decomposition The graph Fourier transform basis is obtained. sum graph frequency ; Step S4: For the graph signal at each time point Perform graph Fourier transform The spectrum of the graph is obtained; Step S5: Calculate the phase difference of each channel relative to the virtual reference. ; Step S6: For each channel, perform a weighted least squares regression of the phase difference against the graph frequency to fit a straight line. Solve for the time offset The weight Take the spectrum amplitude energy; Step S7: According to The original signal is calibrated by fractional delay interpolation, and the synchronized 256-channel EEG signal is output.

2. The method for synchronous calibration of EEG time shift in the 256-channel Fourier domain according to claim 1, characterized in that: In step S2, the adjacency matrix is ​​constructed using a threshold Gaussian kernel function: The Graph Laplace matrix is .

3. The method for synchronous calibration of EEG time shift in the 256-channel Fourier domain according to claim 1, characterized in that: In step S4, before performing the graph Fourier transform, the graph signal is first weighted using a graph domain Hanning window, and the window function is... .

4. The method for synchronous calibration of EEG time shift in the 256-channel Fourier domain according to claim 1, characterized in that: In step S5, the virtual reference signal is taken as the complex average of the spectral coefficients of all channels. Phase difference .

5. The method for synchronous calibration of EEG time shift in the 256-channel Fourier domain according to claim 1, characterized in that: In step S6, the Huber loss function is used for iterative reweighting to suppress outliers, and convergence is achieved after 3 iterations.

6. The method for synchronous calibration of EEG time shift in the 256-channel Fourier domain according to claim 1, characterized in that: In step S7, the fractional delay filter uses fourth-order Lagrange interpolation, with interpolation coefficients of... Output signal .

7. The method for synchronous calibration of EEG time shift in the 256-channel Fourier domain according to claim 1, characterized in that: The preprocessing in step S1 includes mean removal, energy normalization, and electrooculography / electromyography artifact removal based on independent component analysis.

8. The method for synchronous calibration of EEG time shift in the 256-channel Fourier domain according to claim 1, characterized in that: In step S3, the graph Laplacian matrix is ​​normalized, and its eigenvalues ​​satisfy... .

9. The method for synchronous calibration of EEG time shift in the 256-channel Fourier domain according to claim 1, characterized in that: It also includes step S8, which outputs the synchronized signal and the estimated time offset to the neuroimaging system for lead field time delay compensation in source localization.

10. The method for synchronous calibration of EEG time shift in the 256-channel Fourier domain according to claim 1, characterized in that: The entire method requires no external training data during calibration, and the time offset estimation accuracy is better than 0.02 times the sampling interval.

Images