Variation Minimization Based Electrocardiogram Artifacts Removal for Local Field Potentials From Neurostimulator

Local field potential (LFP) recorded by sensing-enabled neurostimulators provided chronic observation of deep brain activities for the research of brain disorders. However, the contamination from the electrocardiogram (ECG) deteriorated the extraction of effective information from LFP. This study proposed a novel algorithm based on minimizing the variance combining template subtraction to improve the performance of ECG artifact removal for LFP. Four patients with implanted electrodes were recruited, and eight real LFP records were collected from their left and right hemispheres, respectively. The results showed that the proposed method improved the accuracy of artifact peak detection in LFP, and the subsequent signal quality after template subtraction compared to the traditional Pan-Tompkins (PT) method. The outcome of this study benefited the LFP-based brain research, promoting the application of sensing-enabled neurostimulators in more areas.


I. INTRODUCTION
D EEP brain stimulation (DBS) is a well-established func- tional neurosurgical technique that is used to treat a variety of neurological and psychiatric diseases [1].This neurosurgical procedure involves precisely implanting electrodes in targeted brain areas by stereotactic techniques.The mechanisms of DBS efficacy are still not clear, and current views mainly attribute it to modulating pathological neuronal network activities [2].The implanted DBS electrode provides a rare opportunity to record deep brain activities directly in the form of local field potential (LFP) from the electrodes stimulating subcortical structures [3].The recorded LFP represented the summed and synchronized activities of a neuronal population near the recording contact [4].The abnormalities of oscillatory synchronization are thought to be part of the development of signs and symptoms of many brain disorders, such as Parkinson's disease (PD) [5], obsessivecompulsive disorder [6], dystonia [7], etc.The recording of the LFP signals opened up the possibility of obtaining the biomarkers of the diseases for diagnosis and treatment.For example, in Parkinson's disease (PD), excessive synchronized beta band (13-30 Hz) activity within the subthalamic nucleus (STN) has been thought to be responsible for Parkinsonian motor symptoms such as rigidity and akinesia [8].It could be used as an electrophysiologic biomarker for guiding DBS lead implantation [9] and stimulation parameters adjusting [10] in PD.
In earlier studies, the LFP recordings in patients who underwent DBS surgery were implemented by intraoperative leads or temporarily externalized leads in the postoperative period [11].These recordings were limited by short duration and an unnatural environment.To address these limitations, the development of neurostimulation devices, such as the sensing-enabled impulse generator (IPG), offered the chance of recording long-term LFP [12].Benefiting the technology of sensing-enabled IPG, large 'real-world' neural LFP datasets were recorded for the discovery of electrophysiologic biomarkers.It improved the efficiency of traditional open-loop DBS.However, open-loop DBS delivered stimulation with fixed parameters, incapable of automatically modifying parameters in real-time to respond to fluctuations in disease states [13].The efficacy of the treatment might be reduced.Recently, the development of adaptive DBS (aDBS), or closed-loop DBS, enabled the automatic adjustment of stimulation parameters in response to the dynamic characteristics of neural oscillations in LFP signals.It was believed to be a more effective therapy technique than the open-loop DBS [14], [15].For closedloop DBS, high quality of LFP signals was essential for the extraction of the accurate biomarker.However, previous studies showed that electrocardiographic (ECG) was commonly presented in the LFP recordings of monopolar mode [16], [17].It was mainly attributed to the fact that the IPG simultaneously served as the reference for the LFP recordings and the anode for the stimulation settings [16].The contaminated LFP signals made it difficult to implement precise parameter adjustments for closed-loop DBS, degrading its efficacy.
Recoding LFP in bipolar mode, i.e., the differential of two adjacent contacts on one electrode, could reject ECG artifacts for regarding them as common-mode signals.However, other common-mode information could also be rejected or reduced in bipolar mode, making the interpretation of neural oscillations not sufficiently accurate [11], [18].Recently, several signal processing algorithms were proposed to eliminate or decrease the influence of ECG artifacts on LFP signals, such as QRS interpolation, adaptive filtering, template subtraction, singular value decomposition (SVD), etc.Among these algorithms, template subtraction was recommended for the advantages including high suppression performance of ECG artifacts [19], [20], no additional reference signal, simple parameter settings, and subjective research input, etc.
The main steps of template subtraction consisted of R peak detection, template extraction, and subtraction.The accuracy of R peak detection was important for the performance of the following steps.The state-of-the-art method of R peak detection from ECG signals was the Pan-Tompkins (PT) algorithm, which employed a series of filters to highlight the rapid increase of the amplitude from heart depolarization [21].Chen et al. employed the PT algorithm to detect R peaks and successfully extracted a template from the simulated LFP signals with ECG artifacts [22].However, for the real LFP signals, the waveforms of ECG artifacts were not as clear as the simulated signals, and the contamination degree could be varied with different electrode and IPG positions [19].For practical applications, the effectiveness of the algorithm needs to be evaluated with real LFP signals.
The performance of template extraction and subtraction depended on the results of R peak detection.As the amplitude of ECG artifacts was much larger than that of LFP signals, the signal variation would be reduced drastically if the ECG artifacts were subtracted correctly.Conversely, incorrect subtraction would increase signal variation.As such, signal variation could be employed as a measurement of peak detection performance from contaminated LFP signals for template subtraction.With this basis, this study proposed a novel algorithm based on variation minimization, termed TSVM, to increase the peak detection accuracy from contaminated LFP signals and improve the subsequent template subtraction performance to suppress the influence of ECG artifacts.Real LFP signals with ECG artifacts were collected from four patients with electrodes implanted in the subthalamic nucleus (STN) after surgery.The performance was compared with PT in peak detection, and noise removal of time and spectral domain after template subtraction.The results of this study would facilitate the development of aDBS with accurate feedback.

II. METHOD A. Subject
Four patients with PD who underwent bilateral implantation of DBS electrodes in the STN (electrode lead model 3389, IPG model SR 1181, SceneRay Inc., Suzhou, China) were recruited in the study.The implantation surgeries were all performed at West China Hospital in 2022.The demographic data is presented in Table I.The patients were instructed to abstain from any dopaminergic medication for a minimum of 12 hours before recording.The locations of two implanted electrodes for one representative patient are displayed in Fig. 1.The recordings were performed within six months after initial DBS while the DBS efficacy was stable.Informed consent forms were signed and received before the experiment.This experimental protocol was in accordance with the Declaration of Helsinki, and approved by the ethics committee of West China Hospital (#.2022923).

B. Data Collection and Preprocessing
The LFP signals were recorded by a sensing-enabled IPG which had been successfully used in our previous work [23].The IPG can synchronically and wirelessly transmit the LFP data to an external computer equipped with a telemetry head at a sampling rate of up to 1,000 Hz.Before the experiment, the subject was instructed to sit comfortably in a chair with the arms supported and asked to relax with eyes open (Fig. 2).
Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

TABLE II INFORMATION OF LFP RECORDS
DBS was ceased for approximately five minutes before the start of LFP recording.Once the LFP sensing mode of the IPG is activated, the number of IPG can subsequently be searched by the software in an external computer.The monopolar mode of the electrodes was selected, and the LFP signals were recorded.There were two recordings for each patient, which were from the electrodes of the left and right hemispheres, respectively.Each recording was 100 s.The data information is displayed in Table II.The number of peaks was visually inspected for each recording.
Data preprocessing was performed before the removal of the artifacts.All the data was first filtered between 1 and 100 Hz using a third-order Butterworth filter to suppress the baseline fluctuations and the noise in the frequency band.Each recording was then normalized to the range between −1 and 1 by the respective maximum of the absolute value for the threshold optimization in the subsequent artifact removal.

C. Artifact Removal
The proposed peak detection method TSVM was based on the minimization of the signal variance after artifact removal.As the amplitude of the peaks from the ECG artifact was larger than that of LFP, the variance of the signal was increased after the introduction of the artifact.Conversely, with artifact removal, the variance would be decreased, which could be used to evaluate the fidelity of the signal.The peak detection and template subtraction were integrated to obtain optimal performance.
The implementation of the proposed peak detection method is as follows (See Algorithm 1).The initial peaks were detected by a method [24] searching the points larger than the surroundings by a predefined threshold.These local maximums were then refined by removing the points of which: 1) the interval from its precedent was smaller than 100 ms; and 2) the amplitude was smaller than the difference between the average and three standard deviations.The procedures were used to remove the abnormal points with respect to interval and amplitude.The remaining points were used to extract the template where the length was half of the peak-to-peak interval.The averaged template was subtracted from the signal at each peak point, and the variance was computed and compared to that before subtraction.With the analysis of the signal variance change after peak removal, the threshold with the minimum variance was adopted for the final peak detection to remove the artifacts.As an example, the relationship between the threshold, signal variance after template subtraction, and the number of correctly detected peaks for Data L1 are displayed in Fig. 3, which is consistent with the analysis.
For the evaluation of the proposed method, the performance of the traditional method, i.e., PT, was calculated.PT employed a series of filters to remove the noise and highlight the content of the interest, with an adaptive threshold for peak detection [21].The comparison was made for the performance of peak detection, as well as the signal quality after template subtraction.

D. Performance Evaluation
For the evaluation of peak detection performance, sensitivity (Sen), Positive predictive value (PPV), and error rate (Err) were employed.They are defined as: Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

Algorithm 1 Peak Detection
Input: LFP data S = {s i } N i=1 Output: the sequence of the detected peaks X Initialization: T ← S 1. for delta = 0.05 to 1 step 0.05 do 2.
for i = 1 to N do 4.
LocalMax ← s i if LocalMax > s i

14.
Perform subtraction with averaged template to obtain the data D if Var(T) > Var(D) 16.
X ← Y and T ← D end 17.end where TP is true positive, representing the number of correctly detected peaks in one record, and FN is false negative, representing the number of missed peaks in one record.A peak is correctly detected if the time provided by the algorithm is within 50 ms of the occurrence, otherwise, it is considered missed.FP is false positive, representing the number of the peaks falsely detected by the algorithm in one record.Sen measures the ability to detect the peaks, and a result of no peak is reliable if the value of Sen is high.PPV measures the reliability of a positive result, i.e., a peak is detected.Err measures the error rate of the peak detection events.
To evaluate the performance of artifact removal after template subtraction, the signal-to-noise ratio (SNR) and root mean square error of the power spectral density (RMSEP) were calculated, quantifying the performance in the time and frequency domain, respectively.The SNR was defined as the power ratio of the signal to the noise for the signal after the template subtraction.As the true LFP signal was unavailable, the sections with strong LFP power and low ECG interference were used to estimate the power of the signal, which discarded the segments half RR interval long and centered at the detected peaks (S −peak ).The sections with low LFP and high ECG interference were used to estimate the power of the noise, which consists of the segments discarded in constructing the signal (S peak ).The SNR is calculated as SNR = 10log 10 P S −peak P S peak (4) where P (S) = E S 2 is the power of the signal S, and E (•) is the expectation.The larger the SNR is, the better the artifact removal is.
The RMSEP evaluates the spectrum difference between the signal and the noise defined above.The RMSEP is calculated as: RMSEP = E F S −peak − F S peak 2 (5) where F (•) is the power spectral density (PSD), which is estimated by Welch's method.The frequency range is set between 1 and 40 Hz, as well as the beta band (13-30 Hz), which is related to Parkinsonian motor symptoms.The smaller the RMSEP is, the better the artifact removal is.

A. Peak Detection
Both TSMV and PT were employed to detect peaks of eight records in Table II.The performance comparison is displayed in Table III.TSVM performed better than PT on TP and FN, and they both had zero FP.It indicated that both TSVM and PT did not have wrongly detected peaks, but PT missed more peaks than TSVM.This resulted in higher sensitivity and lower error rates for TSVM compared to PT. Fig. 4 shows the sample segments of eight records with detected peaks by both methods.The positions detected by TSVM were exactly on the peaks of the signal, while the positions detected by PT Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.were not.As such, the TP of TSVM was higher than that of PT.However, if the distances between the true peaks and the positions detected by PT were within the predefined tolerance, these positions were still regarded as true positives, not false positives.As such, the FP values of both algorithms were zero.

B. ECG Artifact Removal
Template subtraction was performed after peak detection to remove the ECG artifacts from LFP signals.Sample records of the raw signals and the signals after template subtraction are displayed in Fig. 5 (See Supplement for the details between 5 and 10 s).The results after template subtraction were different among the data of the four patients.In general, the number of the remaining peaks from the raw signal was less for TSVM compared to PT.For PT, all the data of the four patients had the remaining peak artifacts.For TSVM, the peak artifacts of the second patient (L2, R2), the third patient (L3, R3), and the fourth patient (L4, R4), were almost removed.However, there were still many peaks with relatively low amplitude in L2 and R2, which might come from the T waves.The average time between the R-peak and the following suspected T-peak was 0.2395 s and 0.2402 s for L2 and R2, respectively, which were in the normal range [25], [26].The deep S waves were also observed in R2 and L2.Besides, it seemed there was a small number of T-waves in R1.These factors affected the noise suppression.The signal quality of L3, R3, L4, and R4 was higher than that of the others after TSVM, and their waveforms after TSVM were similar to the cleaned LFP signals, from which the artifacts cannot be visually identified.
The quantitative measurements of the artifact removal performance are displayed in Table IV.As the beta band from 13 to 30 Hz was important for PD research, the RMSEP was calculated in both the full band and beta band.The performance of PT was inferior to that of TSVM in both SNR and RMSEP for all the data.Besides, the differences of the two metrics between TSVM and PT in L1 and R1 of the first patient were smaller than the differences in the recordings of other patients, which was consistent with the performance of the peak removal in Fig. 5.For RMSEP, the difference between TSVM and PT was increased when the spectrum was reduced from full band to the beta band.The spectrum comparison of the beta band is displayed in Fig. 6.Compared to the raw data, the spectrum became smooth after template subtraction.Further, compared to PT, the spectrum of the data after TSVM was much closer to the spectrum of the data discarding the artifacts.The observations and measurements indicated the signal quality of LPF after TSVM was better than PT.

IV. DISCUSSION
This study proposed a novel algorithm for ECG artifact removal from LFP for accurate observation of deep brain activities.The method was based on the variance minimization Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.The number of remaining peaks after TSVM is smaller than that after PT.The artifacts cannot be identified visually from the waveforms of L3, R3, L4 and R4 after TSVM.
of the signal after template subtraction for the magnitude of the ECG artifacts much higher than that of LFP.Compared to the traditional PT method, the proposed method TSVM demon-strated higher TP and lower FN in peak detection of ECG artifacts from real LFP of the left and right hemispheres of four patients.The signal after TSVM-based template subtraction Fig. 6.Spectrum comparison of beta band (13-30 Hz) among the raw data (Raw), the data after template subtraction (TSVM and PT), and the data discarding the segments of the artifacts, which are centered at the peak and half RR interval long (Clean).The spectrum of TSVM is much closer to that of Clean compared to PT, indicating the signal quality after TSVM is better.

TABLE IV
COMPARISON OF SIGNAL QUALITY AFTER TEMPLATE SUBTRACTION also had higher SNR of the time domain and lower RMSEP of the frequency domain, indicating that the signal quality was improved after artifact removal.As the parameters of template were kept the same for both methods, it indicated that of detection performance enhanced the suppression of the noise.
Except for template subtraction, there were several algorithms proposed in previous studies for ECG artifact removal from LFP [20], such as adaptive filtering, SVD, etc.As the system of the sensing-enabled neurostimulator was implanted in the patient, simplicity would be favored to keep the system stable and save the power for long-term work.The steps of template subtraction were simple compared to other algorithms.No additional reference was required, and the performance was not sensitive to the values of the hyperparameters, which would be hard to select for healthcare professionals.Besides, its high suppression performance of ECG artifacts was demonstrated in previous studies [19], [20].Towards the real-world application, template subtraction was selected in this study.As there were no hyperparameters in TSVM, the entire signal processing was simple and easy to implement for future medical applications.
There were significant individual differences in the degree of ECG artifacts on LFP among these four patients.Several factors could account for these differences, including the specific locations of the DBS electrodes and the IPG, changes in impedance at the tissue-electrode interface, and the potential leakage of fluid into the IPG [16], [20].These factors collectively influence the signal chain of the device and, consequently, contribute to variations in ECG contamination.It's worth noting that the dissimilarities in ECG contamination can significantly impact the performance of artifact removal algorithms.For example, in this study, almost all the peaks were removed for the data from the second to the fourth patient, but not the first.And for the second patient, there were still small peaks after removal, which might be the T waves.As the TSVM focused on the noise around R peaks, the disturbances from other components of ECG remained.For the signal quality of the third and fourth patient after TSVM, the artifacts were visually removed and the power density of beta Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
frequency was close to the clean (Fig. 6).The varied efficiency of artifacts removal indicated the importance of testing with real LFP signals with ECG artifacts from different patients, and the necessity of tailoring the algorithms for the specific patients.
Except for signal quality, there were still some challenges for the applications of aDBS in clinics, including the accurate detection of physiological markers [27], the target population of aDBS, the long-term stability, and the power consumption of complex algorithms for IPG.For the limitation of the study, due to the difficulties of signal collection, only eight 100-second records were obtained from four patients.The performance of the proposed method should be further verified by more patients with longer signals.Besides, though PT was the state-of-the-art in R peaks detection, there were many other effective methods [16], [20].The comparison of TSVM with them needed to be performed in future.Another limitation is that only offline analysis was performed.As the closed-loop DBS was a potential application, the online performance of the artifact removal method was necessary with the extraction of the essential information from the restored signals for accurate feedback.
V. CONCLUSION DBS-LFP recordings are often contaminated with ECG artifacts, leading to LFP oscillations unusable for therapeutic algorithms in closed-loop DBS.In this study, a novel algorithm was proposed to remove ECG artifacts from real LFP in the DBS target of STN in four patients.The proposed method TSVM demonstrated higher TP and lower FN in peak detection of ECG artifacts than the traditional PT method.Further, the LFP signal quality after TSVM was improved with the SNR of the time domain and lower RMSEP of the frequency domain.As there were no hyper-parameters, the whole signal processing of TSVM would be easy to implement in the LFPsensing IPG, facilitating the development of the closed-loop DBS and improving the efficacy.

Fig. 1 .
Fig. 1.MRI image with preoperative surgical plan of Deep Brain Stimulation (DBS) for electrode locations of one representative subject.The sloid white lines represent the electrode locations, and the dotted lines are the electrode path.The square represents the location of the Anterior Commissure (AC) and Posterior Commissure (PC), which serve as the reference points in establishing the coordinate system for the DBS target.The rings represent the entrance to the skull and end of the electrode path, respectively.

Fig. 2 .
Fig. 2. One representative patient sitting in chair and relaxing before the experiment.The patient is asked to sit still during the data collection.

Fig. 3 .
Fig.3.The relationship between the threshold, signal variance after template subtraction, and the number of correctly detected peaks for Data L1 with the proposed TSVM.The signal value is normalized with respect to the maximum.The minimum variance corresponds to the largest number of correctly detected peaks.

Fig. 4 .
Fig.4.Comparison of peak detection performance on the data from 10 to 20 s for the eight records.All the peaks detected by TSVM are exactly on the top of the true peaks, while some peaks detected by PT are on the beginning, such as L2, and the middle, such as R3, of the true peaks, a little away from the local maximal point.

Fig. 5 .
Fig. 5. Waveforms of the signals from 10 to 60 s between the raw and the results after template subtraction based on TSVM and PT, respectively.The number of remaining peaks after TSVM is smaller than that after PT.The artifacts cannot be identified visually from the waveforms of L3, R3, L4 and R4 after TSVM.