Electroencephalogram denoising system integrating multiple artifact recognition and automatic removal methods

Abstract

The embodiment of the invention discloses an electroencephalogram denoising system integrating various artifact recognition and automatic removal methods, which comprises the following steps of: acquiring an original electroencephalogram signal through an electroencephalogram preprocessing step, and preprocessing the original electroencephalogram signal; and identifying electrode peak artifacts, myoelectricity artifacts and electro-oculogram artifacts in the original electroencephalogram signals through an artifact detection and division step, and removing the electrode peak artifacts, the myoelectricity artifacts and the electro-oculogram artifacts through an artifact removal step to obtain target electroencephalogram signals. Therefore, by establishing a unified artifact detection and classification mechanism, the noise is processed in a targeted manner, the loss of effective electroencephalogram information is reduced, the recognition precision is improved, and the calculation overhead is reduced.

Description

Technical Field

[0001] This invention relates to the field of electroencephalogram (EEG) signal monitoring, and more particularly to an EEG denoising method that integrates multiple artifact recognition and automatic removal methods. Background Technology

[0002] EEG (Electroencephalography, scalp electroencephalography) is a low-cost, non-invasive neurophysiological recording technique with high temporal resolution, widely used in fields such as neurological diseases, social psychology, and brain-computer interfaces. EEG signals are weak bioelectrical signals generated by the firing of neurons in the brain, reflecting brain activity. Due to its ease of acquisition, non-invasiveness, and high temporal resolution, EEG plays an increasingly important role in scientific research and disease diagnosis. However, the low signal-to-noise ratio (SNR) of EEG limits its widespread adoption. The main reason for the low SNR is various artifacts, including non-physiological and physiological artifacts. Non-physiological artifacts include those caused by electrode loosening and environmental interference, commonly including transient electrode spike artifacts. Physiological artifacts include electrooculography (EOG) caused by blinking and gliding; electromyography (EMG) caused by facial, neck, and limb muscle movements; and electrocardiography (ECG) caused by heart rhythm. These artifacts are superimposed on normal EEG components, interfering with the analysis of EEG data. Traditional methods include bandpass / notch filtering, wavelet thresholding denoising, and ICA (independent component analysis), but they often lack specificity for different types of noise, easily leading to the loss of effective EEG information. Existing ensemble methods for multiple artifacts lack differentiated modeling for different artifact types, easily mistaking neural signal components as noise to be removed, and incurring high computational costs. Summary of the Invention

[0003] In view of this, embodiments of the present invention provide an EEG denoising method that integrates multiple artifact recognition and automatic removal methods. By establishing a unified artifact detection and classification mechanism, noise is processed in a targeted manner, the loss of effective EEG information is reduced, recognition accuracy is improved, and computational overhead is reduced.

[0004] In a first aspect, embodiments of the present invention provide an EEG denoising method integrating multiple artifact recognition and automatic removal methods, the method comprising: The EEG preprocessing step is used to acquire raw EEG signals and preprocess the raw EEG signals. The artifact detection and classification steps are used to identify electrode spike artifacts, electromyographic artifacts, and electrooculographic artifacts in the original EEG signal; The artifact removal step is used to remove the electrode spike artifacts, electromyographic artifacts, and electrooculographic artifacts to obtain the target electroencephalogram (EEG) signal.

[0005] In some embodiments, the EEG preprocessing step: Low-frequency baseline drift and high-frequency noise in the original EEG signal are removed by zero-phase bandpass filtering, and objective artifact fragments in the original EEG signal are identified and eliminated.

[0006] In some embodiments, the artifact detection and segmentation steps are as follows: The raw EEG signal of a single channel is divided into windows to obtain a first window, and a preset proportion of overlapping area is retained between adjacent first windows; Peak points are obtained by peak detection on the raw EEG signal within the first window, including positive peak points and negative peak points; Peak points with a time interval less than or equal to a preset threshold are grouped together. If the number of peak points in a group exceeds a preset number, the group is deleted. If the number of peak points in a group does not exceed the preset number, the group is determined to be a valid group. The search time range is obtained by weighting the peak points within the effective group based on their amplitudes. Peak detection is performed within the search time range, and the peak point with the largest slope is taken as the peak vertices. The nearest peak point adjacent to the peak is taken as the left and right boundary points of the peak. If the slope between the left and right boundary points and the peak is less than a preset slope threshold, the peak is deleted. If the slope between the left and right boundary points and the peak is greater than or equal to the preset slope threshold, the time position of the peak is recorded to determine the electrode peak artifact.

[0007] In some embodiments, the artifact detection and segmentation steps are as follows: The original electroencephalogram (EEG) signal is divided into segments according to a preset time interval to obtain the first segment; A binary label sequence is obtained based on the first segment using a machine learning model, and the binary label sequence includes the binary label of the first segment; The amplitude envelope is determined based on the binary label sequence; The high-sensitivity artifact interval is determined based on the amplitude envelope using the PELT change point detection algorithm. The candidate change point positions are determined based on the high-sensitivity artifact interval; Electromyographic artifacts are determined based on the candidate variable point locations.

[0008] In some embodiments, the artifact detection and segmentation steps are as follows: Identify segments in the binary label sequence where the consecutive labels are 1, and record their start and end positions; Determine an initial time interval, and then extend the initial time interval by a window with a label of 0 before and after it to determine candidate time intervals; Bandpass filtering is performed on the EEG signals of the candidate time intervals to obtain enhancement signals; The amplitude envelope is determined based on the enhanced signal using the Hilbert transform.

[0009] In some embodiments, the artifact detection and segmentation steps are as follows: The second segment is obtained based on the EEG signals within the sliding window and the high-sensitivity artifact region; Calculate the Euclidean distance between the second segments to obtain the distance matrix; Sort the off-diagonal elements in the distance matrix in ascending order, and extract the median distance of the off-diagonal elements; The kernel function is determined based on the median distance and the radial basis kernel function; The similarity matrix is ​​obtained based on the kernel function and the second segment; The change point detection result is determined based on the similarity matrix using the PELT change point detection algorithm; The candidate change point locations are determined based on the change point check results.

[0010] In some embodiments, the artifact detection and segmentation steps are as follows: If there are two candidate variable points in the test results, the position of the candidate variable point is determined based on the candidate variable point, and electromyographic artifacts are determined based on the position of the candidate variable point; If there are more than two candidate change points in the inspection results, the candidate change points are divided into sub-intervals according to their positions; Calculate the mean, variance, and sum of squares of the EEG signal voltages within the sub-intervals; Electromyographic artifacts are determined based on the mean, variance, and sum of squares of the sub-intervals.

[0011] In some embodiments, the artifact detection and segmentation steps are as follows: The raw EEG signals in the analysis channel and the reference channel are divided into second and third windows according to the preset window length. The overlapping area between adjacent second windows and adjacent third windows is retained by a preset ratio. Zero-phase bandpass filtering is applied to the EEG signals in the second and third windows; Robust standardization is performed on the EEG signals within the second and third windows to obtain second and third standard signals; The average value of the second standard signal is calculated to determine the comprehensive index signal of the forehead channel; The forehead peak value is obtained based on the forehead channel comprehensive index signal and the preset extraction threshold, and the peak time interval is determined based on the forehead peak value. Candidate electrooculogram artifacts are determined based on the third standard signal using the peak time interval and the preset extraction threshold. Determine artifact boundaries based on the candidate electrooculogram artifacts; Based on the artifact boundaries, electrooculography artifacts are identified.

[0012] In some embodiments, the artifact removal step: The artifact-free interval is determined based on the boundary time points of the electrode spike artifacts. The first wavelet coefficients are determined based on the artifact-free interval using the Akima interpolation algorithm, and the electrode spike artifacts are removed based on the first wavelet coefficients.

[0013] In some embodiments, the artifact removal step: The intervals of electromyography artifacts and the intervals of electrooculography artifacts are merged to obtain the physiological artifact intervals. The original EEG signal within the physiological artifact interval is decomposed using wavelet decomposition to obtain the second wavelet coefficients; The quantile distribution is determined based on the second wavelet coefficients; The EEG signals within the physiological artifact intervals are adjusted according to the quantile distribution to remove the electromyographic artifacts and the electrooculographic artifacts.

[0014] The technical solution of this invention obtains raw EEG signals through an EEG preprocessing step, preprocesses the raw EEG signals, identifies electrode spike artifacts, electromyography (EMG) artifacts, and electrooculography (EOG) artifacts in the raw EEG signals through an artifact detection and classification step, and removes these artifacts through an artifact removal step to obtain the target EEG signal. Thus, by establishing a unified artifact detection and classification mechanism, noise can be specifically processed, reducing the loss of effective EEG information, improving recognition accuracy, and reducing computational overhead. Attached Figure Description

[0015] The above and other objects, features and advantages of the present invention will become clearer from the following description of embodiments of the invention with reference to the accompanying drawings, in which: Figure 1 This is a schematic diagram illustrating the steps of the EEG denoising method that integrates multiple artifact recognition and automatic removal methods according to an embodiment of the present invention. Figure 2 This is a flowchart of an EEG preprocessing method according to an embodiment of the present invention; Figure 3 This is a flowchart of an electrode spike artifact detection method according to an embodiment of the present invention; Figure 4 This is a flowchart of a method for identifying and classifying electromyographic artifacts according to an embodiment of the present invention; Figure 5 This is a flowchart of an amplitude envelope determination method according to an embodiment of the present invention; Figure 6 This is a flowchart of a method for determining the position of a candidate change point according to an embodiment of the present invention; Figure 7 This is a flowchart of a method for determining electromyographic artifacts according to an embodiment of the present invention; Figure 8 This is a flowchart of a method for classifying electrooculography artifacts according to an embodiment of the present invention; Figure 9 This is a flowchart of a method for removing electrode spike artifacts according to an embodiment of the present invention; Figure 10 This is a flowchart of a method for removing electromyographic and electrooculographic artifacts according to an embodiment of the present invention; Figure 11 This is a schematic diagram of an electronic device according to an embodiment of the present invention. Detailed Implementation

[0016] The present application is described below based on embodiments, but it is not limited to these embodiments. In the detailed description of the present application below, certain specific details are described in detail. Those skilled in the art can fully understand the present application without these details. To avoid obscuring the substance of the present application, well-known methods, processes, flows, elements, and circuits are not described in detail.

[0017] Furthermore, those skilled in the art should understand that the accompanying drawings provided herein are for illustrative purposes only and are not necessarily drawn to scale.

[0018] Unless the context explicitly requires it, words such as "including" or "contains" throughout the application should be interpreted as including rather than exclusive or exhaustive; that is, meaning "including but not limited to".

[0019] In the description of this application, it should be understood that the terms "first," "second," etc., are used for descriptive purposes only and should not be construed as indicating or implying relative importance. Furthermore, in the description of this application, unless otherwise stated, "a plurality of" means two or more.

[0020] The solutions described in this specification and embodiments, if involving the processing of personal information, will be processed only on the premise of having a legal basis (such as obtaining the consent of the personal information subject, or being necessary for the performance of a contract), and will only be processed within the scope stipulated or agreed upon. A user's refusal to process personal information beyond what is necessary for basic functions will not affect the user's use of basic functions.

[0021] Figure 1 This is a schematic diagram illustrating the steps of an EEG denoising method integrating multiple artifact recognition and automatic removal methods according to an embodiment of the present invention. Figure 1 In the illustrated embodiment, the EEG denoising method includes an EEG preprocessing step 11, an artifact detection and segmentation step 12, and an artifact removal step 13. The EEG preprocessing step 11 is used to acquire raw EEG signals and preprocess them. The artifact detection and segmentation step 12 is used to identify electrode spike artifacts, electromyographic artifacts, and electrooculographic artifacts in the raw EEG signals. The artifact removal step 13 is used to remove the electrode spike artifacts, electromyographic artifacts, and electrooculographic artifacts to obtain the target EEG signal.

[0022] The EEG preprocessing step 11 includes an EEG acquisition sub-step 111 and a preprocessing sub-step 112. The EEG acquisition sub-step 111 acquires the weak voltage signal through the electrode pads, converts the weak voltage signal into a current signal, and then amplifies the current signal to obtain the raw EEG signal.

[0023] The preprocessing sub-step 112 employs a zero-phase bandpass filter of 1~80Hz to remove low-frequency baseline drift and high-frequency noise. The low-frequency baseline drift originates from slow signal shifts caused by breathing, slight body movements, etc., and typically has a frequency less than 1Hz. This can cause an overall shift in the EEG waveform, masking the true EEG components. The high-frequency noise mainly comes from device electromagnetic interference, wire friction, etc., and typically has a frequency greater than 80Hz. It manifests as irregular spikes in the signal, interfering with the identification of high-frequency target signals such as electromyography artifacts.

[0024] Simultaneously, after acquiring the raw EEG signal in the EEG acquisition sub-step, it can be read by dedicated EEG analysis software. The software extracts metadata and voltage data according to the corresponding data format, automatically associates channel names with corresponding voltage sequences, and generates time axis data based on the sampling frequency. Each voltage value is bound to the corresponding time point to form paired data of time and voltage. The parsed time and voltage sequences are reconstructed into the original EEG signal waveform through the software's visualization module. The software defaults to vertically arranging the waveforms by channel (e.g., Fp1 on top, Fp2 in the middle, Cz at the bottom). Each channel independently displays its corresponding time domain waveform, with the vertical axis representing the voltage value and the horizontal axis representing the time axis.

[0025] Preprocessing sub-step 112 can further filter out effective channels with stable signal quality. The total number of effective channels meets the requirement of at least 19 channels for the EEG standard 10-20 system. At the same time, non-critical channels with excessive signal noise or abnormal amplitude fluctuations are removed. If the number of effective channels does not meet the standard, the entire data segment is considered invalid. Preprocessing sub-step 112 is used to remove obvious interference and retain effective data, providing high-quality input for the accurate identification of subsequent electrode spike, electromyography, and electrooculography artifacts, and avoiding invalid data from consuming computing resources or interfering with the detection results.

[0026] The artifact detection and classification step 12 includes an electrode spike artifact detection sub-step 121, an electromyography artifact discrimination and classification sub-step 122, and an electrooculography artifact classification sub-step 123. Specifically, the electrode spike artifact detection sub-step 121 divides the raw EEG signal of a single channel into an overlapping first window, detects positive and negative peak points within the window, groups peak points with time intervals less than or equal to a preset threshold, removes invalid groups with excessive peak points, and calculates a search time range for valid groups using peak point amplitude weighting. Within this range, it locates the peak with the largest slope and its nearest left and right boundary points. If the slope between the boundary point and the peak is greater than or equal to a preset threshold, the peak time position is recorded to determine the electrode spike artifact. Electrode spike artifacts refer to instantaneous non-physiological spike artifacts caused by electrode loosening, environmental interference, and communication interference.

[0027] The electromyography (EMG) artifact identification and segmentation sub-step 122 divides the original EEG signal into a first segment according to a preset time interval. A binary label sequence is output through an ensemble learning model. Candidate time intervals are obtained based on the binary label sequence. After extracting the envelope energy through bandpass filtering at a preset frequency and Hilbert transform, a high-sensitivity artifact interval is determined. The high-sensitivity artifact interval is then divided into a second segment through a sliding window. The Euclidean distance between segments is calculated and the median distance is extracted. A similarity matrix is ​​constructed by combining the radial basis function kernel function. The PELT (Pruned Exact Linear Time) change point detection algorithm is input to obtain candidate change points. If there are two change points, the artifact boundary is directly determined. If there are more than two change points, abnormal intervals are merged by filtering the sub-interval mean, variance, and local energy. Finally, EMG artifacts are determined.

[0028] In step 123 of the electrooculography (EOG) artifact segmentation sub-step, the left frontal channel (Fp1) and right frontal channel (Fp2) are selected as analysis channels. These channels are closest to the eyeball and have extremely high sensitivity to EOG signals, accurately capturing changes in electrical signals generated by blinking and eye movements. Channels unaffected by EOG artifacts are selected as reference channels, such as the central top of the head channel (Cz). This channel is far from the eyes, and the attenuation of EOG signals is almost zero, reflecting only normal brain electrical activity, serving as a benchmark for distinguishing artifacts from normal brain electrical activity. The raw brain electrical signals of the analysis channels (Fp1, Fp2) and the reference channel (Cz) are respectively divided into overlapping second and third windows. Simultaneously, zero-phase bandpass filtering at a preset frequency is performed on all EEG signals within the second and third windows. The filtered signals from the second and third windows are then subjected to robust Z-Score (Standard Score) standardization. The average value of the standardized signal corresponding to the second window for each analysis channel (Fp1, Fp2) is calculated to obtain the comprehensive index signal for the prefrontal channel. Positive and negative peaks with amplitudes exceeding a preset extraction threshold are identified, and the time position and amplitude of each peak are recorded. The preset extraction threshold is manually set, for example, 3 times MAD (Median Absolute Deviation). Using the time position of each peak as the center, preset windows are extended forward and backward to form peak time intervals. For each prefrontal peak time interval, it is checked whether the reference channel has a peak exceeding the preset extraction threshold within that interval. If the reference channel has a peak within that interval, it indicates that the prefrontal peak is normal EEG activity and is removed. If the reference channel has no peak, it is determined to be a candidate EOG artifact peak. In this process, for each candidate electrooculogram (EOG) artifact peak, a boundary expansion operation is performed to determine the complete artifact interval. Centered on the candidate peak, the time range is gradually expanded to the left and right, while simultaneously monitoring signal extreme points. Expansion stops when the expanded interval contains at least three extreme points. The start and end times of the expanded interval are recorded as the boundary of a single candidate EOG artifact. Simultaneously, all candidate EOG artifacts are sorted by start time. If the intervals of two artifacts overlap or the interval interval is less than 0.1 s, they are considered the same EOG event and merged into a single continuous interval.

[0029] The artifact removal step 13 includes an electrode spike artifact removal sub-step 131 and an electromyography (EMG) and electrooculography (EOG) artifact removal sub-step 132. Specifically, the electrode spike artifact removal sub-step 131 first summarizes the start time, end time, peak position, and amplitude of each spike output from the artifact detection and segmentation steps. Based on the obtained spike artifact time range, the sampling points containing these spike artifacts are hollowed out. Then, Akima interpolation is used to fill the hollowed-out signal gaps, thereby eliminating the interference of electrode spike artifacts on the EEG signal.

[0030] In step 132, the removal of electromyographic (EMG) and electrooculography (EOG) artifacts first summarizes the EMG artifact intervals of all channels and the EOG artifact intervals of channels Fp1 and Fp2. Overlapping intervals are merged, intervals with intervals less than a preset time threshold are deleted, and non-overlapping intervals are retained as independent intervals. This results in a set of continuous and non-overlapping physiological artifact intervals (e.g., [0.3s~1.6s, 2.1s~2.5s]), avoiding redundant denoising or artifact omissions. Simultaneously, wavelet decomposition is performed on the signals within the intervals using the db4 wavelet basis to obtain wavelet coefficients in different frequency domains. The reference normal EEG is the clean segment before and after the artifact segment; that is, the reference segment for a certain artifact segment is the segment from the end of the previous artifact segment to the beginning of this artifact segment, and from the end of this artifact segment to the beginning of the next artifact segment. The quantile distribution of the wavelet coefficients of this clean segment is calculated. For wavelet coefficients within the artifact interval, if they exceed a reasonable threshold determined based on the distribution range of normal EEG coefficients, they are judged as abnormal coefficients. Based on the frequency domain distribution characteristics of these abnormal coefficients and their matching degree with known artifact characteristics, they are classified into EEG artifact correlation coefficients or EMG artifact correlation coefficients. The abnormal coefficients are compressed to a reasonable range according to quantile mapping rules while preserving their relative order. Then, all adjusted wavelet coefficients are input into an inverse wavelet transform algorithm to reconstruct the time-domain signal. The reconstructed signal undergoes amplitude calibration to ensure that its mean and standard deviation are consistent with the signal in the artifact-free interval. Finally, the processed signals are integrated channel by channel to obtain the target EEG signal with EMG and EEG artifacts removed.

[0031] This invention acquires raw EEG signals through an EEG preprocessing step, preprocesses the raw EEG signals, identifies electrode spike artifacts, electromyography (EMG) artifacts, and electrooculography (EOG) artifacts in the raw EEG signals through artifact detection and classification steps, and removes these artifacts through an artifact removal step to obtain the target EEG signal. Thus, by establishing a unified artifact detection and classification mechanism, noise can be specifically processed, reducing the loss of effective EEG information, improving recognition accuracy, and reducing computational overhead.

[0032] Figure 2 This is a flowchart of an EEG preprocessing method according to an embodiment of the present invention. Figure 1 The specific details of the steps. For example... Figure 2 As shown, the EEG preprocessing method specifically includes the following steps: Step S110: Remove low-frequency baseline drift and high-frequency noise from the original EEG signal using zero-phase bandpass filtering.

[0033] Specifically, the EEG preprocessing step reads the raw EEG signals from each channel, including channel Fp1 in the left forehead, channel Fp2 in the right forehead, and channel Cz in the center of the top of the head. A zero-phase bandpass filter at a preset frequency is used to remove low-frequency baseline drift and high-frequency noise. Low-frequency baseline drift originates from slow signal shifts caused by breathing, slight body movements, etc. For example, the preset frequency is 1-80Hz, with frequencies typically less than 1Hz, which can cause an overall shift in the EEG waveform, masking the true EEG components. High-frequency noise mainly originates from factors such as electromagnetic interference from equipment and wire friction; its frequency is typically greater than 80Hz, appearing as irregular, spike-like fluctuations in the signal. Filtering out frequencies above 80Hz is primarily due to sampling rate limitations. During EEG signal acquisition, the signal-to-noise ratio is relatively low, making it difficult to extract valuable information from high-frequency signals above 80Hz. Therefore, research is generally conducted only on signals below 80Hz, which also avoids interference from high-frequency noise in the identification of target signals such as electromyography artifacts.

[0034] Step S120: Identify and remove objective artifact fragments from the original EEG signal.

[0035] Specifically, the EEG preprocessing step identifies and removes objective artifact fragments from the original EEG signal based on preset judgment parameters. These judgment parameters can be an amplitude abnormality threshold (±100μV), waveform flattening duration (500ms), high-frequency spike density (50 spikes / second), or impedance threshold (10kΩ). These judgment parameters can be adjusted by the user according to the acquisition device and the scenario. Signals are processed channel by channel using a sliding time window method. For each window, the maximum and minimum amplitudes are calculated. If the maximum value is greater than 100μV or the minimum value is less than -100μV, the window is marked as a suspected artifact due to amplitude anomalies. The standard deviation of the signal within the window is calculated. If the standard deviation is less than or equal to 0.5μV, and this condition is met for five consecutive windows (500ms in total), the time period is marked as a suspected artifact due to waveform flatness. A 100Hz high-pass filter is applied to the signal within the window, and the number of peaks with amplitudes greater than 20μV after filtering is counted. If the number of peaks per second is greater than 50, the window is marked as a suspected artifact due to high-frequency peaks. If the electrode impedance at the corresponding time point of the window is greater than 10kΩ, and the window has already been marked as any of the above suspected artifacts, the marking is upgraded to a high-confidence suspected artifact. Consecutive suspected artifact windows are merged to form complete suspected artifact segments. The start time, end time, channel number, and objective artifact type of each segment are recorded. At the same time, operators can identify artifacts caused by loose electrodes, electrodes not contacting the scalp, amplifier malfunctions, electromagnetic interference from equipment, or poor wire contact, and no further processing is required.

[0036] Figure 3 This is a flowchart of a method for detecting electrode spike artifacts according to an embodiment of the present invention. Figure 3 As shown, Figure 1 The specific steps of the electrode spike artifact detection method include the following: Step S210: Divide the raw EEG signal of a single channel into windows to obtain a first window, and retain a preset proportion of overlapping area between adjacent first windows.

[0037] Specifically, using a fixed window length and a preset overlap ratio, the acquired raw EEG signal is divided into N consecutive overlapping first windows by sliding from the beginning to the end. For example, assuming a 10s EEG signal, the fixed window length is 1s, the overlap ratio is 10%, the signal in the first window is 0~1s, the signal in the second window is 0.9s~1.9s (overlapping with the first window by 0.1s), the signal in the third window is 1.8s~2.8s (overlapping with the second window by 0.1s), and so on, until the EEG signal ends. The number of sampling points per window is obtained by multiplying the sampling frequency by the window length, for example: sampling frequency 256Hz × window length 1s = 256 sampling points (numbered 0~255). The window sliding step size is obtained by subtracting the overlap duration from the window length, for example: window length - overlap duration = 1s - 0.1s = 0.9s.

[0038] Step S220: Perform peak detection on the original EEG signal in the first window to obtain peak points, including positive peak points and negative peak points.

[0039] Specifically, for each sampling point within the first window, the voltage values ​​of adjacent sampling points are compared one by one. If the voltage value of a sampling point is greater than the voltage value of the previous sampling point and greater than the voltage value of the next sampling point, it is determined to be a positive peak point. If the voltage value of a sampling point is less than the voltage value of the previous sampling point and less than the voltage value of the next sampling point, it is determined to be a negative peak point. The time and amplitude of each peak point are recorded.

[0040] Step S230: Group peak points with time intervals less than or equal to a preset threshold into a group. If the number of peak points in a group exceeds a preset number, delete the group. If the number of peak points in a group does not exceed the preset number, it is determined to be a valid group.

[0041] Specifically, the peak points within each first window are sorted chronologically, and the time interval between two adjacent peak points is calculated. If the time interval is less than or equal to a preset threshold, they are grouped together. Based on this, the number of peak points within each group is determined. If the number of peak points in a group exceeds a preset number, it is identified as a dense noise group, and all peak points in that group are deleted and not included in subsequent processes. If the number of peak points in a group does not exceed the preset number, it is identified as a valid group, and the complete information of all peak points in that group is retained.

[0042] Step S240: Obtain the search time range by weighted calculation based on the amplitude of the peak points within the effective group.

[0043] Specifically, the effective group contains n peak points, and the parameters of each peak point are (t1, v1), (t2, v2), ..., (t... n ,v n ), where t is the time axis and v is the amplitude axis.

[0044] The total amplitude is calculated using the following formula: ; The amplitude percentage of each peak point is calculated using the following formula: (i=1~n); The weighted time mean is calculated using the following formula: ; Centered on the weighted time mean t_{core}, expand the window length forward and backward by half (e.g., 0.05s) to form a search time range of [t_{core}-0.05s, t_{core}+0.05s]. If the search range exceeds the current window boundary, it is truncated to the window boundary.

[0045] Step S250: Perform peak detection within the search time range and take the peak point with the largest slope as the peak vertex.

[0046] Specifically, the search time range of the current valid group is extracted. From the sampling point data in the first window, all sampling points with timestamps within the search range are selected to form a list of sampling points within the search range. Peak detection is then re-executed according to the rules in step S220 to ensure that no core peaks are missed. The slope between adjacent peak points Si(ti,vi) and Sj(tj,vj) (j=i+1) within the search range is obtained using the following formula:

[0047] Among them, a positive slope (slope_ij > 0) indicates that the amplitude increases from Si to Sj, and a negative slope (slope_ij < 0) indicates that the amplitude decreases.

[0048] Calculate the absolute value of the slope of all adjacent peak pairs, find the slope with the largest absolute value and denote it as slope_max, and the corresponding peak pair is denoted as (S_prev, S_curr). If slope_max is greater than 0, it is a positive slope, the amplitude of S_curr is greater than the amplitude of S_prev, and S_curr is the extreme value in the peak pair, which is determined to be a peak vertex. If slope_max is less than 0, it is a negative slope, the amplitude of S_prev is less than the amplitude of S_curr, and S_prev is the extreme value in the peak pair, which is determined to be a peak vertex. Each valid group corresponds to one unique peak vertex. The output format of the peak detection can be {peak vertex timestamp: t_vertex, amplitude: v_vertex, polarity: positive / negative, corresponding slope: slope_max}.

[0049] Step S260: Take the nearest peak point adjacent to the peak as the left and right boundary points of the peak. If the slope between the left and right boundary points and the peak is less than a preset slope threshold, delete the peak. If the slope between the left and right boundary points and the peak is greater than or equal to the preset slope threshold, record the time position of the peak to determine the electrode peak artifact.

[0050] Specifically, within the search time range, the two closest peak points adjacent to the peak vertex are found. The closest one on the left is the left boundary point, and the closest one on the right is the right boundary point. If the peak vertex is the first peak in the search range, only the right boundary point is taken; if it is the last peak, only the left boundary point is taken. If the search range contains only one peak vertex, one sampling point before and after the vertex is taken as virtual boundary points. If the absolute value of the slope of either boundary is greater than or equal to a preset threshold, it is determined to be an electrode peak artifact, and the precise time position of the peak vertex is recorded. If the absolute values ​​of the slopes of both boundaries are less than the preset threshold, it is determined to be a normal EEG peak, and only that peak vertex is deleted without affecting other peak points in the same valid group. At the same time, for peak vertices identified by adjacent windows, if the time interval between the peak vertices is less than or equal to a preset overlapping time interval, they are determined to be the same peak artifact, and only the vertex with the middle time position is recorded.

[0051] Figure 4 This is a flowchart of a method for identifying and classifying electromyographic artifacts according to an embodiment of the present invention. Figure 4 As shown, Figure 1 The specific steps of the electromyography artifact discrimination and classification method include the following: Step S310: Divide the original EEG signal into segments according to a preset time interval to obtain the first segment.

[0052] Specifically, the electromyography (EMG) artifact detection and segmentation steps begin at the start of the original EEG signal, and the original EEG signal is divided into first segments at preset time intervals, with no overlap between the first segments. In this non-overlapping segmentation method, the preset time interval is equal to the duration of each first segment. For example, if the window length is set to 2 seconds, then the step size should also be 2 seconds, so that the sequentially segmented first segments will not overlap and can completely and non-repeatingly cover the original EEG signal. Furthermore, the preset time interval can also be greater than the duration of each first segment to achieve non-overlapping segmentation.

[0053] Step S320: Obtain a binary label sequence based on the first segment using a machine learning model. The binary label sequence includes the binary labels of the first segment.

[0054] Specifically, the EMG artifact detection and segmentation step labels each first segment. Segments containing EMG artifacts, where the frequency band energy percentage is greater than or equal to a preset proportion or the amplitude fluctuation range is greater than or equal to a preset EEG amplitude range, are identified as positive samples and labeled as 1. Segments without EMG artifacts, conforming to the normal EEG amplitude range (within ±50μV), are identified as negative samples and labeled as 0. The 0 and 1 samples are obtained by manually labeling real EEG data to obtain positive samples containing EMG artifacts (labeled 1) and negative samples of normal EEG without EMG artifacts (labeled 0). Alternatively, Gaussian white noise with different signal-to-noise ratios can be superimposed on publicly available clean EEG datasets to generate corresponding sample data, which are then labeled according to the above judgment rules, serving as another source of samples for training the model. For each first segment, 16-dimensional features are extracted to form a feature vector. These 16 features include Shannon entropy, second-order Galley entropy, delta-band Shannon entropy, beta-band Shannon entropy, γ-band Shannon entropy, Lyapunov exponent, fractal dimension, Hjorth activity, Hjorth complexity, median frequency, kurtosis, first-order autocorrelation coefficient, second-order autocorrelation coefficient, third-order autocorrelation coefficient, fourth-order autocorrelation coefficient, and coefficient of variation. Simultaneously, the feature vector is divided into training, validation, and test sets according to a certain ratio.

[0055] The Shannon entropy is obtained based on the amplitude histogram within the window to determine the probability. probability distribution satisfy The calculation formula is: Where i is the index of the amplitude interval, i = 1, 2, 3, ..., N, and N is the total number of intervals in the amplitude histogram. A more dispersed amplitude distribution indicates stronger randomness and higher entropy. Electromyography artifacts with more dispersed amplitude distributions often exhibit higher Shannon entropy.

[0056] Galley entropy is the Shannon entropy used to calculate non-additive systems, and it better characterizes long-tailed distributions or non-stationary signals. When the non-Gaussian distribution characteristics of electromyographic artifacts are more pronounced, Tsallis entropy is more sensitive than Shannon entropy. Its formula is... Where i is the index of the amplitude interval, i = 1, 2, ..., B, B is the total number of intervals in the amplitude histogram, and q is the order. Here, second-order Galley entropy is used, and the formula is... .

[0057] Spectral entropy refers to the Shannon entropy calculated for energy in a specific frequency band. Here, the spectral entropy of three frequency bands, δ (1~4Hz), β (13~30Hz), and γ (30~80Hz), is calculated to reveal baseline drift and high-frequency interference caused by electromyography artifacts, which can be distinguished from clean electroencephalogram (EEG) signals.

[0058] The Lyapunov exponent describes the exponential divergence rate of adjacent trajectories in phase space, and is used to measure the degree of chaos and short-term instability. Its calculation formula is as follows: Where t is the sequence number of the signal sampling point, t=1, 2, ..., n, and n is the total number of valid sampling points within the calculation window. Let be the instantaneous derivative of the signal at the t-th sampling point. Electromyographic artifacts introduce rapid irregular components, which are generally more chaotic, highly nonlinear, and have a larger Lyapunov exponent.

[0059] Fractal dimension with different scaling The length of the curve below and The dimension of the log-linear relationship is estimated, and the slope is fitted using the least squares method: , where k is the step size parameter for signal segmentation, with a maximum value of 10. This is used to measure signal complexity, while electromyographic artifacts introduce more irregular components, resulting in a higher fractal dimension.

[0060] Hjorth Parameters include two characteristics: activity and complexity. Activity reflects frequency characteristics, and is calculated using the following formula: , wherein Original signal The first derivative with respect to time, the original signal This refers to unprocessed EEG signals, reflecting the instantaneous rate of change of the signal at each moment. Complexity reflects changes in the signal's morphology, and is calculated using the following formula: , wherein Original signal The second derivative with respect to time, i.e., the first derivative The derivative of the signal reflects the rate of change of the instantaneous change rate. Electromyographic artifacts have many high-frequency components, thus exhibiting high activity and correspondingly increased complexity.

[0061] The median frequency is the frequency at which the cumulative power in the power spectral density reaches half of the total power. The calculation formula is: , wherein The median frequency, the For power spectral density, This is the maximum analysis frequency. Electromyography artifacts enhance high-frequency energy in EEG, shifting the median frequency towards higher frequencies.

[0062] Kurtosis measures the kurtosis of a distribution, and is calculated using the following formula: Where x is the amplitude of a single sampling point, the The mean of the amplitudes of all sampling points within the window. This represents the standard deviation of the amplitude values ​​at all sampling points within the window. Normal EEG is close to a Gaussian distribution (kurtosis ≈ 3), but EMG artifacts contain transient high-amplitude fluctuations, which significantly increase the kurtosis.

[0063] The autocorrelation coefficient reflects the temporal correlation of a signal. Normal EEG has strong rhythmicity (high correlation), while EMG artifacts approximate noise, causing its autocorrelation value to decay rapidly. Here, only the first to fourth order autocorrelation coefficients are calculated: Where t is the sampling point index of the original signal, Let be the amplitude of the original signal at the t-th sampling point. The amplitude of the original signal at the (t+k)th sampling point, Let N be the square of the amplitude at the t-th sampling point, and N be the number of sampling points within the window.

[0064] The coefficient of variation measures the magnitude of the standard deviation relative to the mean, and is calculated using the following formula: ,in Indicates standard deviation, The value represents the mean. Electromyography artifacts show dramatic fluctuations and large values ​​(CVs), while normal EEG shows more stable fluctuations and smaller values ​​(CVs).

[0065] The standardized training set feature vectors, labels, and validation set data are imported into the model. Based on the initial hyperparameters and a custom loss function, the XG Boost (eXtreme Gradient Boosting) model is initialized, and an early stopping mechanism is set, for example, stopping if there is no improvement in AUC after 10 consecutive rounds. The custom loss function is defined using the `objective` parameter of XG Boost, and the formula is as follows:

[0066] The gradient calculation formula is:

[0067] The formula for calculating the second derivative is:

[0068] Where y_i is the true label (1 = positive sample, 0 = sample), y^i is the positive sample probability predicted by the model (0~1), and wfp is the initial false positive weight factor.

[0069] The trained weak learners are combined to predict the training set, outputting the positive sample probability y^i for each sample. This is then substituted into the loss function with wfp to calculate the loss value for each sample. Simultaneously, the negative gradient of the loss function is calculated, and a new weak learner is trained using this negative gradient as the target value. The trained and newly trained weak learners are then used to predict the validation set, and two core metrics are calculated: False positive rate (FPR): (FP represents the number of negative samples that were misclassified as positive samples, and TN represents the number of negative samples that were correctly classified.) True positive rate (TPR): (TP is the number of positive samples correctly identified, and FN is the number of positive samples misidentified as negative samples).

[0070] Simultaneously, the false positive rate (wfp) is adjusted based on the validation set FPR results. Assuming wfp is 1.5, if the FPR is greater than the first threshold (0.05), the false positive rate is too high, and wfp is increased by 0.2 (e.g., from 1.5 to 1.7) to strengthen the penalty for false positive errors. If the FPR is less than or equal to the first threshold (0.05) but greater than the second threshold (0.03), the false positive rate is acceptable, and wfp remains unchanged. If the FPR is less than or equal to the second threshold (0.03), the false positive rate is too low, which may lead to missed positive samples, and wfp is decreased by 0.1 (e.g., from 1.5 to 1.4) to balance false positive control and recall. At the same time, the weights of the new weak learners are calculated using a greedy algorithm, and the new weak learners are integrated into the ensemble model according to the calculated weights to form the updated model.

[0071] Simultaneously, the validation set AUC (Area Under the Curve) value is calculated to evaluate the model's overall discriminative ability. If the AUC improvement is less than 0.001 over 10 consecutive training rounds, the model is considered converged, and training is stopped. Furthermore, the core hyperparameters can be adjusted using Bayesian optimization with the goal of maximizing the validation set AUC, until the core requirements for EMG artifact recognition are met.

[0072] Step S330: Determine the amplitude envelope based on the binary label sequence.

[0073] Specifically, the model identifies the leading and trailing edges of potentially missed EMG artifacts by first obtaining candidate time intervals based on a binary label sequence. Then, by expanding the window with a label of 0, it strives to ensure complete coverage of the EMG artifacts. A bandpass filter at a preset frequency is applied to the signal; this operation preserves the high-frequency characteristics of the EMG artifacts while filtering out low-frequency normal EEG signals and high-frequency device noise, thereby enhancing the energy difference between the EMG artifacts and normal EEG signals. The filtered signal is then subjected to a Hilbert transform to obtain an analytic signal, and its amplitude is calculated as the amplitude envelope.

[0074] Figure 5 This is a flowchart of an amplitude envelope determination method according to an embodiment of the present invention. Figure 5 As shown, Figure 1 The specific steps of the amplitude envelope determination method include the following: Step S331: Identify segments in the binary label sequence where the consecutive labels are 1, and record their start and end positions.

[0075] Specifically, the binary label sequence is traversed to identify all segments with consecutive labels of 1. The start and end segment numbers of each segment are recorded. Based on the time range corresponding to the segment number, the actual time range of the consecutive 1 segments is calculated. For example, if the start segment is (1.0~1.5s) and the end segment is (2.0~2.5s), then the time range is [1.0s, 2.5s].

[0076] Step S332: Determine the initial time interval, and expand the initial time interval with a window labeled 0 before and after it to determine the candidate time interval.

[0077] Specifically, for each initial time interval, a segment labeled 0 is extended before and after it. If the extended segment exceeds the range of the original signal, the signal boundary is used as the limit. For example, for the initial time interval [1.0s, 2.5s], a segment is extended forward (0.5~1.0s) and then backward (2.5~3.0s) to obtain the candidate time interval [0.5s, 3.0s].

[0078] Step S333: Bandpass filtering is performed on the EEG signals of the candidate time intervals to obtain enhancement signals.

[0079] Specifically, for the EEG signals within each candidate time interval, a preset frequency zero-phase bandpass filter is performed to obtain the enhanced signal, preserving the high-frequency characteristics of EMG artifacts, filtering out low-frequency normal EEG signals and high-frequency device noise, and enhancing the energy difference of the artifact signals.

[0080] Step S334: Determine the amplitude envelope based on the enhanced signal using Hilbert transform.

[0081] Specifically, the Hilbert transform is performed on the enhanced signal to obtain the analytic signal, and its amplitude is calculated as the amplitude envelope.

[0082] Step S340: Determine the high-sensitivity artifact interval based on the amplitude envelope using the PELT change point detection algorithm.

[0083] Specifically, the PELT change point detection algorithm can automatically identify points in the amplitude envelope where the energy changes significantly, dividing the amplitude envelope into different intervals. The intervals with significant energy changes identified by the PELT change point detection algorithm are designated as high-energy intervals, and the timestamps corresponding to these high-energy intervals are recorded. Adjacent high-energy intervals are merged into continuous segments. For example, high-energy intervals with a time interval less than or equal to 0.02 s are considered adjacent high-energy intervals, and each continuous segment is a candidate high-sensitivity interval. If the interval between candidate high-sensitivity intervals is less than or equal to 0.1 s, these candidate high-sensitivity intervals are merged into a continuous interval, and the high-sensitivity artifact interval is finally determined based on the merged continuous interval.

[0084] Step S350: Determine the candidate change point position based on the high-sensitivity artifact interval.

[0085] Specifically, a second segment is acquired from the enhanced signal in the high-sensitivity artifact region using a specific sliding window and step size. The Euclidean distance between segments is calculated to obtain a distance matrix. The median of the off-diagonal elements is then taken, and the kernel function is determined based on this and the radial basis function kernel function formula. This leads to the similarity matrix, which is then used in conjunction with the PELT change point detection algorithm and the similarity matrix to determine the change point detection results. Finally, change points exceeding the interval are removed, and candidate change point locations are determined based on change point distance or penalty term verification. This provides a systematic process for accurately determining candidate change point locations from the high-sensitivity artifact region of EEG signals, facilitating the precise capture of key signal changes and providing a reliable basis for subsequent EEG signal analysis and processing.

[0086] Figure 6 This is a flowchart of a method for determining the location of candidate change points according to an embodiment of the present invention. Figure 6 As shown, Figure 1 The specific steps of the method for determining candidate change point locations include the following: Step S351: Obtain the second segment based on the EEG signals within the sliding window and the high-sensitivity artifact region.

[0087] Specifically, the enhancement signal within the high-sensitivity artifact region is extracted, and the second segment is obtained by dividing it into segments with a 0.2s sliding window and a 0.1s step size (50% overlap).

[0088] Step S352: Calculate the Euclidean distance between the second segments to obtain the distance matrix.

[0089] Specifically, for all second segments within the same high-sensitivity artifact interval, calculate the Euclidean distance pairwise. Let the number of second segments be M, then the distance matrix is ​​an M×M matrix, with diagonal elements being 0 (distance between a segment and itself), and off-diagonal elements representing the Euclidean distance dij between corresponding segments. For example, when M=3, the distance matrix is... , where dij is the Euclidean distance between segments i and j.

[0090] Step S353: Sort the off-diagonal elements in the distance matrix in ascending order, and extract the median distance of the off-diagonal elements.

[0091] Specifically, the off-diagonal elements of the distance matrix are extracted, sorted in ascending order, and the median value of the sorted sequence is taken as the median distance. For example, d 12 = 2.3 =4.1, d21=2.3, d23=1.8, d31=4.1, d32=1.8. After sorting the off-diagonal elements, the result is [1.8, 2.3, 2.3, 4.1, 4.1, 4.5]. The median is... .

[0092] Step S354: Determine the kernel function based on the median distance and the radial basis kernel function.

[0093] Specifically, the kernel function is determined using the RBF kernel function based on the median distance and the radial basis function, and the specific formula is as follows:

[0094] in equal , For radial basis kernel functions, .

[0095] Step S345: Obtain the similarity matrix based on the kernel function and the second segment.

[0096] Specifically, for any two second segments x_i and x_j, the Euclidean distance dij is substituted into the RBF kernel function to calculate the kernel function value k(x_i,x_j). The similarity matrix is ​​an M×M symmetric matrix, and the element value is the kernel function value k(x_i,x_j), which takes values ​​in the range [0,1]. The closer the value is to 1, the more similar the segments are.

[0097] For example, M=3, =0.0488, d_{12}=2.3; , The similarity matrix is ​​as follows:

[0098] Step S356: Determine the change point detection result based on the similarity matrix using the PELT change point detection algorithm.

[0099] Specifically, the change point detection result is determined by the PELT change point detection algorithm similarity matrix, the preset penalty term, and the timestamp information of the second segment.

[0100] Step S357: Determine the candidate change point positions based on the change point check results.

[0101] Specifically, change points whose timestamps exceed the high-sensitivity artifact range are removed. The remaining change points are then validated, and key parameters are introduced to further filter change points and determine candidate change point locations. These key parameters are the distance between change points or a penalty term. Specifically, if the timestamp interval between two change points is less than or equal to a preset change point distance threshold (this threshold can be set according to actual signal characteristics and requirements, for example, 0.05s), the change point with the smaller timestamp is retained, while the change point with the larger timestamp is removed. This avoids too many densely packed change points in areas where local signal changes are small. The penalty term limits the total number of change points, preventing excessive unnecessary change points due to overfitting. In the change point detection algorithm, a larger penalty term coefficient results in a stricter limit on the number of change points and fewer output change points. Conversely, a smaller penalty term coefficient allows for a larger number of change points. Based on actual needs and data characteristics, the penalty term coefficient is set reasonably to ensure that key signal changes are captured while controlling the number of change points within a reasonable range.

[0102] Step S360: Determine electromyographic artifacts based on the candidate variable point locations.

[0103] Specifically, for scenarios with different numbers of candidate variable points, electromyographic artifacts are determined based on the positions of the candidate variable points. Two variable points correspond to the start and end boundaries of an artifact, directly forming a complete artifact interval. Multiple variable points split a highly sensitive artifact interval into multiple sub-intervals. The true artifact sub-intervals need to be screened out using statistical features and then merged into a complete artifact.

[0104] Figure 7 This is a flowchart of a method for determining electromyographic artifacts according to an embodiment of the present invention. Figure 7 As shown, Figure 1 The specific steps of the electromyography artifact detection method include the following: Step S361: If there are two candidate change points in the inspection results, determine the position of the candidate change point based on the candidate change point, and determine the electromyographic artifacts based on the position of the candidate change point.

[0105] Specifically, if there are two candidate change points in the inspection results, the two candidate change points are sorted in ascending order of timestamps. The first change point (timestamp t1) is the artifact start boundary, and the second change point (timestamp t2) is the artifact end boundary. Electromyographic artifacts are determined based on the candidate change point positions [t1, t2] of the two change points. At the same time, it is ensured that the interval duration is greater than or equal to 0.05s (the minimum duration of electromyographic artifacts); if it is less than 0.1s, it is discarded.

[0106] Step S362: If there are more than two candidate change points in the inspection results, divide the candidate change points into sub-intervals according to the change point positions.

[0107] Specifically, if there are more than two candidate change points in the inspection results, they are sorted in ascending order by timestamp, and two adjacent change points form a sub-interval.

[0108] Step S363: Calculate the average value, variance, and sum of squares of the EEG signal voltage within the sub-interval.

[0109] Specifically, for each sub-interval, extract the enhanced signal sampling points within that interval and calculate the mean, variance, and sum of squares.

[0110] Wherein, the average value is: K is the number of sampling points in the sub-interval. This represents the voltage value at the k-th sampling point; variance:

[0111] Sum of squares: .

[0112] Step S364: Determine electromyographic artifacts based on the mean, variance, and sum of squares of the sub-intervals.

[0113] Specifically, the characteristics of each sub-interval are compared with the mean, variance, and sum of squares of the artifact-free interval. A sub-interval is identified as an artifact sub-interval if its variance is greater than or equal to a preset multiple and its local energy is greater than or equal to a preset multiple of its sum of squares. Based on these artifact sub-intervals, electromyographic artifacts are determined. Simultaneously, adjacent artifact sub-intervals with a time interval less than or equal to a preset time threshold are merged into a continuous interval, and isolated artifact sub-intervals with a duration less than the preset time threshold are removed.

[0114] Figure 8 This is a flowchart of a method for classifying artifacts related to electrooculography according to an embodiment of the present invention. Figure 8 As shown, Figure 1 The specific steps of the electrooculography artifact classification method include the following: Step S410: Divide the raw EEG signals in the analysis channel and the reference channel into a second window and a third window according to the preset window length. The second window and the third window retain a preset ratio of overlapping area.

[0115] Specifically, the analysis channel is divided into second windows according to a preset window length, with a preset overlap ratio maintained between adjacent second windows. Similarly, the reference channel is divided into third windows according to the preset window length, with a preset overlap ratio maintained between adjacent third windows. The analysis channels are located at the left forehead (Fp1) and right forehead (Fp2) positions. These channels are at least 2 cm from the eyeball and have a sensitivity of at least 0.1 μV to voltage changes in electrooculography (EOG) activity, accurately capturing electrical signals generated by blinking and eye movements. The reference channel is located at the center of the top of the head (Cz) position. This channel is at least 10 cm from the eyeball, has an EOG signal attenuation rate at least 90%, and only reflects normal brain activity, serving as a benchmark for distinguishing artifacts from normal brain activity.

[0116] Step S420: Perform zero-phase bandpass filtering on the EEG signals in the second and third windows.

[0117] Specifically, FIR zero-phase bandpass filtering is performed independently on the original signals within each second and third window.

[0118] Step S430: Robustly standardize the EEG signals in the second and third windows to obtain the second and third standard signals.

[0119] Specifically, a second standard signal and a third standard signal are obtained based on the EEG signals within the second and third windows using a robust standardization formula, which is as follows:

[0120] in, The median of the signal within the window. This represents the absolute deviation of the median. The calibration coefficients are used to make the MAD of normally distributed data equivalent to the standard deviation, ensuring that the signal scale is uniform after standardization.

[0121] Step S440: Calculate the comprehensive index signal of the forehead channel of the second standard signal.

[0122] Specifically, the second standard signals Fp1 and Fp2 are matched according to the window number, the signal value of each timestamp is averaged, and the forehead channel comprehensive index signal is determined based on the average value, where the formula is:

[0123] and The second standard signal values ​​for Fp1 and Fp2 at the same timestamp.

[0124] Step S450: Obtain the forehead peak value based on the forehead channel comprehensive index signal and the preset extraction threshold, and determine the peak time interval based on the forehead peak value.

[0125] Specifically, based on the statistical characteristics of the forehead channel comprehensive index signal, a preset extraction threshold is established, positive and negative peaks in the comprehensive index signal whose amplitudes exceed the preset extraction threshold are identified, the time position of each peak is recorded, and the peak time interval is determined by extending 0.05s forward and backward from the time position of each peak as the center.

[0126] Step S460: Determine candidate electrooculogram artifacts based on the third standard signal using the peak time interval and the preset extraction threshold.

[0127] Specifically, for each peak time interval, the third standard signal in the reference channel within that interval is extracted. If the amplitude of the third standard signal within that interval does not exceed the preset extraction threshold, it is determined to be a candidate electrooculogram artifact. If it exceeds the threshold, it indicates normal electroencephalogram activity and is removed.

[0128] Step S470: Determine the artifact boundaries based on the candidate electrooculogram artifacts.

[0129] Specifically, for each candidate electrooculogram (EOG) artifact, the time range is gradually expanded to the left and right based on the peak time interval, while simultaneously monitoring signal extreme points. The process stops when the expanded interval contains at least three extreme points, ensuring coverage of the complete EOG event including the rising edge, peak, and falling edge. Simultaneously, the start and end times of the expanded interval are recorded as the boundaries of a single candidate EOG artifact.

[0130] Step S480: Determine the electrooculography artifacts based on the artifact boundaries.

[0131] Specifically, all candidate electrooculogram (EOG) artifacts are sorted by start time. If the intervals of two artifacts overlap or the interval is less than a preset merging threshold, they are merged into a continuous interval. At the same time, extremely short intervals with a duration less than a preset time threshold are removed, and valid intervals are retained, thereby identifying EOG artifacts.

[0132] Figure 9 This is a flowchart of a method for removing electrode spike artifacts according to an embodiment of the present invention. Figure 9 As shown, Figure 1 The specific steps of the electrode spike artifact removal method include the following: Step S510: Determine the artifact-free interval based on the boundary time points of the electrode spike artifacts.

[0133] Specifically, the electrode spike artifact detection step outputs the boundary time points of the electrode spike artifacts, with the start time t1 and the end time t2. Using the boundary time points of the electrode spike artifacts as the boundary, the signal segment outside the spike t1~t2 is locked, which is the artifact-free interval.

[0134] Step S520: Determine the first wavelet coefficients based on the artifact-free interval using the Akima interpolation algorithm, and remove the electrode spike artifacts based on the first wavelet coefficients.

[0135] Specifically, sampling points within the artifact-free interval are extracted and input into the Akima interpolation algorithm. The Akima interpolation algorithm calculates the reconstructed voltage value (i.e., the first wavelet coefficient) for each sampling point within the peak interval [t1, t2] based on the voltage value and timestamp of the sampling point. The logic for the reconstructed voltage value is as follows: ensuring that the amplitude of the reconstructed signal at t1 equals the amplitude v_0 at the end of the left artifact-free interval, and the amplitude at t2 equals the amplitude v_1 at the beginning of the right artifact-free interval, and that the slope of the reconstructed curve smoothly connects with the slopes of the signals in the left and right artifact-free intervals. The calculated first wavelet coefficients are then used to replace all sampling point values ​​within the peak interval [t1, t2] of the original signal, thereby removing the electrode peak artifacts.

[0136] Figure 10 This is a flowchart of a method for removing electromyographic and electrooculographic artifacts according to an embodiment of the present invention. Figure 10 As shown, Figure 1 The specific steps of the electromyography and electrooculography artifact removal methods include the following: Step S610: Merge the intervals of the electromyography artifacts and the intervals of the electrooculography artifacts to obtain the physiological artifact intervals.

[0137] Specifically, the system receives the intervals of electromyographic artifacts output from the electromyographic artifact discrimination and segmentation step and the intervals of electrooculography artifacts output from the electrooculography artifact segmentation step. The electromyographic artifact intervals and electrooculography artifact intervals in the same channel are grouped together, and physiological artifact intervals are obtained by merging them channel by channel.

[0138] Step S620: Perform wavelet decomposition on the original EEG signal within the physiological artifact interval to obtain the second wavelet coefficients.

[0139] Specifically, the original EEG signals within each physiological artifact interval are extracted independently by channel to form interval signal sequences. A wavelet decomposition function is then called to perform a 5-level db4 wavelet decomposition on the interval signals to obtain the second wavelet coefficients. Simultaneously, the decomposed coefficients can be stored by frequency band classification, with each group of coefficients containing sampling point information in the time dimension.

[0140] Step S630: Determine the quantile distribution based on the second wavelet coefficients.

[0141] Specifically, from the raw EEG signals acquired from the current channel, it is first necessary to identify the physiological artifact regions. Physiological artifact regions are signal segments containing interference such as electromyographic (EMG) artifacts and electrooculography (EOG) artifacts. Normal EEG signals are the portion of the raw signal excluding these physiological artifact regions. For example, in a continuously recorded EEG data stream, specific detection algorithms or manual labeling are used to identify time periods with obvious artifacts; the signals corresponding to the remaining time periods are considered normal EEG signals. For each artifact segment, the reference normal EEG (clean segment) is selected from the portions before and after the artifact segment. Specifically, the reference segments for a given artifact include the time period from the end of the previous artifact segment to the beginning of this artifact segment, and the time period from the end of this artifact segment to the beginning of the next artifact segment. The db4 wavelet basis is used to perform wavelet decomposition on the normal EEG signals and the signals within the artifact regions to be processed. The signals are decomposed into different frequency domains, yielding wavelet coefficients at different levels. These wavelet coefficients contain information about the signal at different frequencies and times. The wavelet coefficients of normal EEG are divided into frequency bands and correspond one-to-one with the second wavelet coefficients. The wavelet coefficients of each group of normal EEG are sorted in ascending order, and the quantile distribution of the wavelet coefficients of each group of normal EEG is calculated.

[0142] Step S640: Adjust the EEG signal within the physiological artifact interval according to the quantile distribution to remove the electromyographic artifacts and the electrooculographic artifacts.

[0143] Specifically, for wavelet coefficients within the artifact interval, a reasonable threshold is determined based on the previously calculated distribution range of normal EEG coefficients. This threshold can be set according to quantile distribution, for example, selecting the value corresponding to a specific quantile as the threshold limit. If the wavelet coefficients within the artifact interval exceed this reasonable threshold, they are judged as abnormal coefficients. These abnormal coefficients are likely caused by interference from EMG or EMG artifacts, and show significant differences from the wavelet coefficient characteristics of normal EEG. Based on the frequency domain distribution characteristics of the abnormal coefficients and their matching degree with known artifact characteristics, they are classified as EMG artifact correlation coefficients or EMG artifact correlation coefficients. The classified abnormal coefficients are then compressed to a reasonable range according to the quantile mapping rule. The quantile mapping rule is a method that adjusts abnormal coefficients to a distribution range similar to normal coefficients based on the quantile distribution of normal EEG wavelet coefficients. During the adjustment of abnormal coefficients, the relative order of the coefficients must be preserved. All adjusted wavelet coefficients (including normal coefficients and adjusted abnormal coefficients) are input into the inverse wavelet transform algorithm. Through the inverse transform, the decomposed wavelet coefficients can be reconstructed back into the time domain signal, obtaining the EEG signal with preliminary artifact removal. Simultaneously, the reconstructed signal is amplitude-calibrated to ensure that the target EEG signal after artifact removal is consistent with the normal EEG signal in amplitude, while preserving the time-frequency characteristics of normal EEG and leaving no obvious artifact residue.

[0144] This invention acquires raw EEG signals through an EEG preprocessing step, preprocesses the raw EEG signals, identifies electrode spike artifacts, electromyography (EMG) artifacts, and electrooculography (EOG) artifacts in the raw EEG signals through artifact detection and classification steps, and removes these artifacts through an artifact removal step to obtain the target EEG signal. Thus, by establishing a unified artifact detection and classification mechanism, noise can be specifically processed, reducing the loss of effective EEG information, improving recognition accuracy, and reducing computational overhead.

[0145] Figure 11 This is a schematic diagram of an electronic device according to an embodiment of the present invention. (For example...) Figure 11 As shown, Figure 11 The illustrated electronic device is an EEG noise reduction device, which includes a general computer hardware architecture, comprising at least a processor 91 and a memory 92. The processor 91 and memory 92 are connected via a bus 93. The memory 92 is adapted to store instructions or programs executable by the processor 91. The processor 91 can be a standalone microprocessor or a collection of one or more microprocessors. Thus, the processor 91 executes the instructions stored in the memory 92, thereby performing the method flow of the embodiments of the present invention as described above to process data and control other devices. The bus 93 connects the aforementioned components together, and also connects these components to a display controller 94, a display device, and an input / output (I / O) device 95. The input / output (I / O) device 95 can be a mouse, keyboard, modem, network interface, touch input device, motion-sensing input device, printer, and other devices known in the art.

[0146] Those skilled in the art will understand that embodiments of this application can be provided as methods, apparatus (devices), or computer program products. Therefore, this application can take the form of a completely hardware embodiment, a completely software embodiment, or an embodiment combining software and hardware aspects. Furthermore, this application can take the form of a computer program product embodied on one or more computer-readable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0147] This application is described with reference to flowchart illustrations of methods, apparatus (devices), and computer program products according to embodiments of this application. It should be understood that each step in the flowchart can be implemented by computer program instructions.

[0148] These computer program instructions may be stored in a computer-readable storage medium that can direct a computer or other programmable data processing device to function in a particular manner, such that the instructions stored in the computer-readable storage medium produce an article of manufacture including an instruction means, the implementation process of which is described in the instruction means. Figure 1 The function specified in one or more processes.

[0149] These computer program instructions may also be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing device to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing device, produce instructions for implementing processes. Figure 1 A device for a function specified in one or more processes.

[0150] Another embodiment of the present invention relates to a non-volatile storage medium for storing a computer-readable program for use by a computer to execute some or all of the above-described method embodiments.

[0151] That is, those skilled in the art will understand that all or part of the steps in the methods of the above embodiments can be implemented by a program specifying the relevant hardware. This program is stored in a storage medium and includes several instructions to cause a device (which may be a microcontroller, chip, etc.) or processor to execute all or part of the steps of the methods described in the embodiments of this application. The aforementioned storage medium includes various media capable of storing program code, such as a USB flash drive, a portable hard drive, a read-only memory (ROM), a random access memory (RAM), a magnetic disk, or an optical disk.

[0152] The above description is merely a preferred embodiment of this application and is not intended to limit this application. Various modifications and variations can be made to this application by those skilled in the art. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of this application should be included within the protection scope of this application.

Claims

1. A brainwave denoising method integrating multiple artifact recognition and automatic removal methods, characterized in that, The method includes: The EEG preprocessing step is used to acquire raw EEG signals and preprocess the raw EEG signals. The artifact detection and classification steps are used to identify electrode spike artifacts, electromyographic artifacts, and electrooculographic artifacts in the original EEG signal; The artifact removal step is used to remove the electrode spike artifacts, electromyographic artifacts, and electrooculographic artifacts to obtain the target electroencephalogram (EEG) signal.

2. The method according to claim 1, characterized in that, The EEG preprocessing steps: Low-frequency baseline drift and high-frequency noise in the original EEG signal are removed by zero-phase bandpass filtering, and objective artifact fragments in the original EEG signal are identified and eliminated.

3. The method according to claim 1, characterized in that, The artifact detection and segmentation steps are as follows: The raw EEG signal of a single channel is divided into windows to obtain a first window, and a preset proportion of overlapping area is retained between adjacent first windows; Peak points are obtained by peak detection on the raw EEG signal within the first window, including positive peak points and negative peak points; Peak points with a time interval less than or equal to a preset threshold are grouped together. If the number of peak points in a group exceeds a preset number, the group is deleted. If the number of peak points in a group does not exceed the preset number, the group is determined to be a valid group. The search time range is obtained by weighting the peak points within the effective group based on their amplitudes. Peak detection is performed within the search time range, and the peak point with the largest slope is taken as the peak vertices. The nearest peak point adjacent to the peak is taken as the left and right boundary points of the peak. If the slope between the left and right boundary points and the peak is less than a preset slope threshold, the peak is deleted. If the slope between the left and right boundary points and the peak is greater than or equal to the preset slope threshold, the time position of the peak is recorded to determine the electrode peak artifact.

4. The method according to claim 1, characterized in that, The artifact detection and segmentation steps are as follows: The original electroencephalogram (EEG) signal is divided into segments according to a preset time interval to obtain the first segment; A binary label sequence is obtained based on the first segment using a machine learning model, and the binary label sequence includes the binary label of the first segment; The amplitude envelope is determined based on the binary label sequence; The high-sensitivity artifact interval is determined based on the amplitude envelope using the PELT change point detection algorithm. The candidate change point positions are determined based on the high-sensitivity artifact interval; Electromyographic artifacts are determined based on the candidate variable point locations.

5. The method according to claim 4, characterized in that, The artifact detection and segmentation steps are as follows: Identify segments in the binary label sequence where the consecutive labels are 1, and record their start and end positions; Determine an initial time interval, and then extend the initial time interval by a window with a label of 0 before and after it to determine candidate time intervals; Bandpass filtering is performed on the EEG signals of the candidate time intervals to obtain enhancement signals; The amplitude envelope is determined based on the enhanced signal using the Hilbert transform.

6. The method according to claim 4, characterized in that, The artifact detection and segmentation steps are as follows: The second segment is obtained based on the EEG signals within the sliding window and the high-sensitivity artifact region; Calculate the Euclidean distance between the second segments to obtain the distance matrix; Sort the off-diagonal elements in the distance matrix in ascending order, and extract the median distance of the off-diagonal elements; The kernel function is determined based on the median distance and the radial basis kernel function; The similarity matrix is ​​obtained based on the kernel function and the second segment; The change point detection result is determined based on the similarity matrix using the PELT change point detection algorithm; The candidate change point locations are determined based on the change point check results.

7. The method according to claim 6, characterized in that, The artifact detection and segmentation steps are as follows: If there are two candidate variable points in the test results, the position of the candidate variable point is determined based on the candidate variable point, and electromyographic artifacts are determined based on the position of the candidate variable point; If there are more than two candidate change points in the inspection results, the candidate change points are divided into sub-intervals according to their positions; Calculate the mean, variance, and sum of squares of the EEG signal voltages within the sub-intervals; Electromyographic artifacts are determined based on the mean, variance, and sum of squares of the sub-intervals.

8. The method according to claim 1, characterized in that, The artifact detection and segmentation steps are as follows: The raw EEG signals in the analysis channel and the reference channel are divided into second and third windows according to the preset window length. The overlapping area between adjacent second windows and adjacent third windows is retained by a preset ratio. Zero-phase bandpass filtering is applied to the EEG signals in the second and third windows; Robust standardization is performed on the EEG signals within the second and third windows to obtain second and third standard signals; The average value of the second standard signal is calculated to determine the comprehensive index signal of the forehead channel; The forehead peak value is obtained based on the forehead channel comprehensive index signal and the preset extraction threshold, and the peak time interval is determined based on the forehead peak value. Candidate electrooculogram artifacts are determined based on the third standard signal using the peak time interval and the preset extraction threshold. Determine artifact boundaries based on the candidate electrooculogram artifacts; Based on the artifact boundaries, electrooculography artifacts are identified.

9. The method according to claim 1, characterized in that, The artifact removal steps are as follows: The artifact-free interval is determined based on the boundary time points of the electrode spike artifacts. The first wavelet coefficients are determined based on the artifact-free interval using the Akima interpolation algorithm, and the electrode spike artifacts are removed based on the first wavelet coefficients.

10. The method according to claim 1, characterized in that, The artifact removal steps are as follows: The intervals of electromyography artifacts and the intervals of electrooculography artifacts are merged to obtain the physiological artifact intervals. The original EEG signal within the physiological artifact interval is decomposed using wavelet decomposition to obtain the second wavelet coefficients; The quantile distribution is determined based on the second wavelet coefficients; The EEG signals within the physiological artifact intervals are adjusted according to the quantile distribution to remove the electromyographic artifacts and the electrooculographic artifacts.

Images