Powerline noise from recording biomedical devices have been known to introduce distortion to recorded signals and, as a result, compromise their integrity and negatively affect their interpretations. Consequently, fields such as neural engineering, neurosurgery, cardiology, and drug discovery need to refine raw data by eliminating powerline noise. An adulteration of neural signals by powerline noise makes it difficult to understand the properties of ion channels and how ensembles of neurons interact to perform specific computations for observed behaviors. In effect, clinical applications of neural interfaces such as brain-machine and brain-computer interfaces become difficult to implement. Neurosurgical procedures often involve recording neural activity for intraoperative monitoring, localization of brain regions, and proper placement of stimulating electrodes for deep brain stimulation (Finnis et al., 2003; Guridi et al., 1999; Schramm et al., 1990). Special rhythms of the brain such as high synchronous intra-cortical activity in the gamma frequency band (25 Hz to 100 Hz) often lead to dyskinesia (Brown, 2006). Since powerline noise exists within the frequency range of some of these important neural oscillations, it is essential to extinguish it to make valid neuro-pathological conclusions and, in succession, develop appropriate approaches for treatment. Details about the condition of the heart and, in turn, the flow of blood in the body is vital in disease diagnosis. By the same token, electrocardiography has become an indispensable clinical tool in cardiology. The electrocardiogram (ECG) is a representation of the electrical activity of cardiac tissue during systole and diastole. This signal is typically obtained by measuring the electrical potential between a specific spatial combination of recording electrodes. As the cardiac depolarization-mediated mean electrical vector approaches a positive electrode, the ECG signal increases positively in amplitude and vice versa. The P wave, QRS complex and T wave in ECG signals are manifestations of atrial depolarization, ventricular depolarization and ventricular repolarization respectively. Although its origin is still under debate, an inversion of the U wave has been shown to be an indicator of myocardial ischemia (Miwa et al., 1993; Gerson et al., 1979; Van Eck, Kors & van Herpen, 2005). In the presence of powerline noise, a malformation in the morphology of the ECG waveform—and by extension, the U wave—is very difficult to detect with the naked eye. Lastly, microelectrodes have been demonstrated to be very useful in drug discovery (Cui et al., 2006; Jones et al., 2011; Jain & Muthuswamy, 2008; Stett et al., 2003). Signal perturbations influenced by powerline noise make it difficult to understand drug-tissue interactions for the development of clinically-viable pharmaceutical products prior to clinical trials. For these reasons, the elimination of powerline noise has been a vigorous field of investigation (Keshtkaran & Yang, 2014; Agrawal & Gupta, 2013; Van Alste & Schilder, 1985; Levkov et al., 2005; Piskorowski, 2012; Poungponsri & Yu, 2013; Keshtkaran & Yang, 2012; Poornachandra & Kumaravel, 2008; Clancy, Morin & Merletti, 2002; Philips, 1996).
Powerline noise is characterized by a chronic sinusoidal 50/60 Hz element which can be observed in raw recordings of biomedical data. The sinusoidal component is usually a result of the use of devices that employ alternating current as a source of power. Alternating current has been used in the design of biomedical devices because it has been demonstrated to possess the quality of being relatively stable, especially over long distances, as opposed to direct current. In some cases, powerline noise is removed by low-pass filters with cut-off frequencies below 50/60 Hz. Although this approach solves the problem of extinguishing powerline noise and has its applications, it is challenging to implement on some forms of biomedical data because of the importance of broadband signals. These include, but are not limited to, electroencephalogram signals and extracellular neural recordings. For some purposes, such as the extraction of action potentials from broadband neural tissue recordings, the pitfall of obtaining noisy recordings can be potentially avoided by employing a high-pass filter with a cut-off frequency above 250 Hz (Oweiss, 2010). Band-stop filters have been used to attenuate powerline noise to inconspicuous levels. Nevertheless, due to the instability of biomedical signals, band-stop filters sometimes fail in reducing noise with 50/60 Hz center frequency and thus may have to rely on correction methods (Ferdjallah & Barr, 1994; Ferdjallah & Barr, 1990; Hamilton, 1996). Although band-stop Butterworth filters obtain better powerline noise removal results with increasing filter orders, their step responses are usually characterized by ringing and overshoot. In the same vein, depending on the filter order, the rate of attenuation can be low. Although the rate of attenuation in Chebyshev filters are higher than Butterworth filters, their step responses are marked by higher levels of ringing. Bessel filters do have significantly lower levels of ringing and overshoot, however their slower attenuation rates make it possible for powerline noise to leak into the signal.
Many signal processing algorithms have been proposed to solve the issue of powerline noise interference. Some of the notable procedures involve blind source separation. Empirical mode decomposition has been proposed as a potent approach for eliminating powerline noise (Blanco-Velasco, Weng & Barner, 2008; De-xiang, Xiao-pei & Xiao-jing, 2008; Li et al., 2012; Naji, Firoozabadi & Kahrizi, 2011; Nimunkar & Tompkins, 2007; Zhang & Zhou, 2013; Chang, 2010; Chang & Liu, 2011; Zivanovic & González-Izal, 2013). Independent component analysis has also been explored by many researchers as a potential approach for removing powerline noise (Xue et al., 2006; Iriarte et al., 2003; Castellanos & Makarov, 2006; Kuzilek et al., 2014; Chawla, 2011; Kuzilek, Kremen & Lhotska, 2014). A combination of empirical mode decomposition and independent component analysis has also been looked into as a potential approach for eliminating powerline noise (Mariyappa et al., 2015).
This paper proposes an algorithm which uses blind source separation and wavelet analysis to detect and remove powerline noise in biomedical signals. The approach is subsequently compared with a band-stop 4th order Butterworth filter, empirical mode decomposition, independent component analysis and a combination of empirical mode decomposition with independent component analysis. This unsupervised machine learning approach is fully automatable and void of the need to apply adaptive self-correction mechanisms.
The proposed approach is motivated by potency of the ensemble empirical mode decomposition algorithm in decomposing a signal into amplitude–frequency modulations—whose linear combination is results in the reconstruction of the signal. In this framework, these modulations are assumed to be a linearly mixed representation of source signals. Some of the modulations are then selected based on their frequency properties and un-mixed into their statistically independent sources via independent component analysis. The Morlet wavelet was used to describe the time-frequency properties of the decompositions and independent component analysis aided in extracting an optimal data-driven model of powerline noise. An inversion of the wavelet transform at 50/60 Hz of the powerline noise model was used to recover the powerline noise that exists within the signal. The result was subtracted from the original signal to obtain its denoised version.
Current algorithms that employ a combination of empirical mode decomposition and independent component analysis require the user to manually select the amplitude–frequency modulations and/or independent components of interest. An innovation of the proposed framework is that the signals of interest are automatically selected based on a pre-defined frequency threshold around their wavelet transformations. Another innovation of the proposed technique is that the inverse wavelet transform is applied on the powerline noise model—which was obtained via ensemble empirical mode decomposition and independent component analysis—to extract the 50/60 Hz component.
Pseudo-code for proposed approach
The data used for the evaluation of this approach consisted of broadband neural activity from ensembles of hippocampal neurons, electrocardiographs and electroencephalogram signals (EEG signals).
One set of data was obtained from the Collaborative Research in Computational Neuroscience initiative (Mizuseki et al., 2009). The recordings were made in the Cornu Ammonis 1 layer of the right dorsal hippocampus of Long Evans rats during open field foraging. This data was sampled at 20 kHz.
Three randomly chosen sets of data with ECG signals from seemingly healthy volunteers were also employed (Goldberger et al., 2000; Garcia-Gonzalez et al., 2013; Citi, Poli & Cinel, 2010). These were sampled at 5 kHz.
Six data sets of EEG signals from humans used for brain-computer interface applications were used (Goldberger et al., 2000). Three of them, which were previously used to elucidate the limitations of using event-related potentials for brain-computer interface applications, were sampled at 2,048 Hz (Citi, Poli & Cinel, 2010). The other half of the data sets were aimed at general purpose brain-computer interface design; they were sampled at 160 Hz (Schalk et al., 2008).
Sliding time window
For the purpose of computational speed, a sliding window with no overlap was used.
Decomposition of raw signal into amplitude-frequency modulations
Empirical mode decomposition
Suppose we want to eliminate a chronic sinusoidal xHz noise from a signal y(t) in a single channel. We can employ empirical mode decomposition (EMD) to split the time series signal into narrow-band amplitude–frequency sections (Huang et al., 1998; Wang et al., 2014). These sections, called intrinsic mode functions (IMFs), are obtained by an iterative sifting procedure. Initially, the local maxima and minima of the signal are detected and are connected by cubic splines to form the upper and lower envelopes respectively. The mean of the envelopes is then subtracted from the original signal to give rise to the first IMF. Subsequently, the first IMF is subtracted from the original signal to result in the residue. This process is repeated k times with the residue obtained at the end of each iteration serving as the input for the next. Ultimately, k IMFs and one residue will be obtained. The signal y(t) can be reconstructed by a summation of the IMFs ci(t) and the residue r(t): (1)
The method by which the local extrema are detected is described in the Appendix (Ensemble empirical mode decomposition: detecting the local extrema).
Ensemble empirical mode decomposition
The ensemble empirical mode decomposition algorithm (EEMD) is an extension of EMD which is aimed at making the extraction of IMFs a robust process (Wu & Huang, 2009). In this algorithm, Gaussian noise of zero mean and a specified standard deviation is added to the input signal before EMD is applied. For the extraction of each IMF, the addition of Gaussian noise to the input is done over a specified number of ensembles and the mean extract in all the ensembles of putative IMFs is selected as the true IMF. This is a truly noise-assisted method of blind source separation because Gaussian noise forces the EMD algorithm to consider all options when sifting. With a high enough number of ensembles for each EEMD iteration, it can be inferred by central limit theorem that the mean of the ensembles is representative of the most likely IMF; thus, the fittest survive. The input signal for each ensemble can be summarized as follows: (2) where nl is the desired inverse signal-to-noise ratio of the input signal and ehin is the input signal for each ensemble h.
Powerline noise detection
Due to the fact that IMFs obtained via EEMD have the attribute of being amplitude–frequency modulations (Fig. 1A), it is essential to un-mix the cocktail into their statistically independent sources. In principle, independent constituents of y(t) with specific frequencies—such as a chronic sinusoidal xHz noise—can be extracted with the aid of ci(t) and r(t). However, due to the high level of variability in the frequency spectrum for each IMF and the occasional overlap in their frequency spectra (Fig. 1B), selecting which IMFs will undergo independent component analysis will help reduce the probability of having noisy estimates of source components. Each IMF selected should be a putative powerline noise in accordance with the properties of its frequency.
IMF selection via wavelet analysis
A wavelet φ(a,b)(t) is square-integrable function that can be dilated (or constricted) along a and translated along b; it is written in the following form: (3) For this framework, the real component of the consecutive Morlet wavelet was employed: (4) where f and j represent the center frequency (scale) and translation respectively. In order to describe the frequency properties of each independent component, a projection of the Morlet wavelet unto ci(t) = [c1(t), …, ck(t), r(t)]T was used. This projection p(F, i) was by accomplished by finding the frequency f in the set of frequencies that maximizes a pseudo-convolution between φ(f,j)(t) and ci(t) and a summation of the pseudo-convolution values at all temporal locations t: (5) where jk, tk and fk are integers representing the kth translation, time and frequency evaluated respectively. A summation of the pseudo-convolution values at all temporal locations exposes the overall frequency morphology that describes a particular IMF. As mentioned previously, although each consecutive IMF exist within a lower frequency range, there is a likelihood of obtaining IMFs with partially overlapping frequency spectra (Fig. 1B). Thus, by choosing IMFs with dominant frequencies within a specified range, there is an appreciable degree of assurance that the xHz signal is in either of the chosen IMFs. Essentially, with Eq. (5), if xHz − bHz ≤ p(F, i) ≤ xHz + bHz (with bHz serving as a threshold for IMF selection), then ∃i such that i is the index that defines ci(t) as a putative xHz powerline noise. In the following subsections, any ci(t) such that xHz − bHz ≤ p(F, i) ≤ xHz + bHz is denoted as nq(t). This method is further elucidated in the Appendix (IMF selection).
Independent component analysis
Independent component analysis (ICA) is a data-driven approach to blind source separation of mixed input signals into their independent components. The premise of ICA is the assumption that each vectorial input signal nq(t) in n = [n1(t), …, nk(t)]T is an observed signal which is a linearly mixed representation of its intrinsic source component si(t) and components of statistically independent origin (Hyvärinen, Karhunen & Oja, 2004). This assumption infers the following: (6) where A is the mixing matrix and the source component s = [s1(t), …, sk(t)]. Inevitably, finding the un-mixing matrix R to obtain a matrix of output signals u which are estimates of si(t) is the essence of ICA: (7)
In this implementation of ICA, R was determined by an approximate simultaneous diagonalization of the third and fourth order cumulant tensors. Independent component analysis: finding R of the Appendix provides an explanation of how R was computed.
Powerline noise recognition
The powerline noise in u is recognized by extending the method outlined in Eq. (5) to search for the ui whose frequency properties resemble that of xHz powerline noise the most: (8) where x represents the expected frequency of the powerline noise (50 Hz or 60 Hz). In essence, the ui(t) that satisfies Eq. (8) is the most appropriate approximation of powerline noise. The motivation behind this approach has been explicated in Powerline noise recognition of the Appendix.
Powerline noise is modeled by inverting a convolution between dcomp(t) and the real portion of the Morlet wavelet, and normalizing the root-mean-square amplitude—via learning and estimating the least square distance aggregate—to compensate for the loss of amplitude information due to ICA (Appendix, Signal reconstruction): (9)
In the same light as Eq. (1), the denoised signal can be reconstructed in the following manner: (10)
Results and Discussion
In this section, properties of the proposed algorithm are explored. The performance of the proposed denoising algorithm is evaluated in the presence of artificial and natural powerline noise. Subsequently, the performance with varying SNRs is compared with that of a band-stop 4th order Butterworth filter, empirical mode decomposition, independent component analysis and the combination of empirical mode decomposition with independent component analysis.
Noise extraction and signal reconstruction
Characteristics of the proposed approach
To demonstrate the viability of the process of powerline noise removal, artificial 60 Hz noise was added to local field potentials (SNR = 4.2816dB; amplitude = 700) and was subsequently reconstructed. A 200 ms time window (without overlap) was used to denoise 250 ms of the mentioned neural data. As shown in Figs. 2A and 2C, the procedure was able to recover a sufficient amount of the original signal with little Manhattan distances. The morphology of the Manhattan distance suggests that the procedure extracted the appropriate frequency, but obtained relatively erroneous amplitude information. As shown in Figs. 2B and 2D, the pseudo-convolution procedure for extraction of powerline noise chose the 60 Hz noise—which is the right frequency for the noise model in this context. Since the pseudo-convolution between the extracted alternating current noise peaks at about 60 Hz on all translations, it implies that the summation of the pseudo-convolution values on every translation should result in a function that also peaks at about 60 Hz. Thus, whichever independent component that has its summed pseudo-convolution peak closest to 60 Hz is indeed the desired powerline noise. Although the final round of blind-source separation in this algorithm requires prior knowledge about the statistical properties of the source data being approximated, it is plausible that this constraint is not a hindrance but, instead, facilitates the process. With the assumption of statistical independence, at least one of the resulting approximations of the source data is forced to look undeniably unique. Signals of this form are usually some fluctuations that are chronically present in the mixed data. For this reason, ICA effectively serves as a helper for the extraction and identification of the desired powerline noise to be removed. The Manhattan distance between the original signal and the denoised signal increased during the final 50 ms of the neural data. This suggests the approach does not provide good results with very small time windows. Essentially, the best approach will be to run the algorithm without a time window. In conclusion, it is most appropriate for offline signal analysis.
Learning the noise amplitude
In accordance with the fact that the problem of deconvolution is ill-posed, the solutions are not unique; this issue is tackled by learning the appropriate amplitude for the noise model. To exhibit the amplitude learning procedure, an EEG with natural 50 Hz noise was denoised using the proposed framework. As presented in Fig. 3C, the learning procedure aimed at minimizing the cost function has a unique solution. This indicates that there is some level of consistency in the process. In view of the fact that there is virtually no ground truth with natural AC noise, the power spectra of the original signal and the denoised signal were compared. The traces in Fig. 3A show that there is little difference between the original signal and the recovered signal. However, Fig. 3B reveals that this manifestation is primarily due to the low amplitude of the 50 Hz AC noise. Although the amplitude is very small, the algorithm was able to extinguish an appreciable amount of powerline noise. The same figure indicates that if the amplitude of the AC noise is small enough, it might go unnoticed after Fourier transformation. Further, in order to obtain an accurate view of how a biomedical signal behaves under various frequency specifications, it is better to employ the wavelet transform (which is a convolution between a wavelet and a signal) rather than the Fourier transform. This is due to the fact that the Fourier transform assumes the signal being transformed is of a sinusoidal origin (Bloomfield, 2004), which is not true in most cases. Although the short-term Fourier transform is used to avoid this pitfall, the wavelet transform (which handles frequencies in a logarithmic fashion) adapts better to highly variable signals (Daubechies, 1990). The suggested algorithm employed a transformation analytically similar to the wavelet transform. By the same token, the frequency properties revealed by a summation of the pseudo-convolution between the real component of the Morlet wavelet and the signals aided in reliably in extracting powerline noise; this, in turn, provided a solid model foundation for learning the most appropriate amplitude.
Artificial corruption of biomedical signals with powerline noise
Sub-cortical local field potentials and extra-cranial EEGs were resampled to have equal sampling rate, concatenated and adulterated with artificial powerline noise. A band-stop 4th order Butterworth filter was used to filter out signals between 59.5 Hz and 60.5 Hz, and the proposed framework was also applied to the adulterated signal. The introduction of powerline noise introduced a slight amplitude decorrelation between the original signal and the adulterated signal. The aim of Fig. 4 was to characterize the level of re-correlation after removing the powerline noise using the proposed approach and the mentioned infinite impulse response filter (IIR filter). It was noted that although the frequency spectra indicated that both methods expelled powerline noise, their correlations differed minutely. The proposed approach recovered the signals with a slightly higher correlation than the named IIR filter (proposed = 1.0000, Butterworth filter = 0.9999).
In order to understand the subtle differences between the proposed method and the mentioned IIR filter, the Manhattan distance between the Fourier transformation of the original signal which had no powerline noise and the denoised signal using both approaches were evaluated (Fig. 5). It is shown that the difference between the power spectrum of the original signal and the one reconstructed using the proposed approach was smaller than that of the band-stop Butterworth filter. This was further validated via a Wilcoxon rank-sum test; the p-value obtained by comparing the original power spectrum with that of the one reconstructed using the proposed approach was 0.9903, while that of the band-stop Butterworth filter was 0.7627. In spite of the fact that there was no statistically significant difference between the reconstructed power spectra and that of the original, there was a statistically significant difference between the power spectrum of the signals reconstructed using the suggested algorithm and that of the band-stop IIR filter (p-value < 0.0001). On a grand scale, the power spectra of the denoised signals were highly correlated. This is an indication that there was a relatively similar response at the powerline noise frequency for neural signals.
Electrocardiogram signals without powerline noise were corrupted with artificial 60 Hz powerline noise. Thereafter, they were reconstructed via band-stop Butterworth filtering and the proposed method. The results obtained were almost identical in the time series (Fig. 6A), however there was a substantial difference circa 60 Hz. Similar to the results obtained with neural signals, the power of the extract using the proposed method at 60 Hz was lower than that of the Butterworth filter. One of the many properties of ECG signals sought after is the waveform. In accordance with the knowledge that the ECG signal used was a concatenation of three data sets with two electrode recordings each, it is expected that at least six different clusters ECG waveforms can be detected. The waveforms were detected via simple thresholding and were partitioned via the k-means algorithm. Naturally, clustering such high dimensional data requires dimensionality reduction via principal component analysis or laplacian eigenmaps. This clustering process was done without dimensionality reduction in order to accurately compare the effect of the proposed approach and that of the Butterworth filter on the reconstructed signal. As shown in Fig. 6B, although the labels of the classes were different, their elements were not. This further proves that the relative temporal morphology of each reconstructed signal is preserved.
Natural powerline noise removal
To evaluate the utility of the proposed method in a natural setting, an EEG data set with 50 Hz powerline noise was used. As shown in Fig. 7, the proposed framework and the band-stop 4th order Butterworth filter both removed the powerline noise. This further proves the potency of the proposed method.
Comparison with other approaches
The performance of the proposed approach was compared with that of a 4th order band-stop Butterworth filter, EEMD, ICA and a combination of EEMD and ICA (EEMD-ICA). The log-mean squared error between the original signal and the signal recovered after eliminating powerline noise was used as a parameter to compare the performance of the approaches. With Fig. 8, it is shown that the proposed approach performed better than the approaches mentioned previously. It can be noted that EEMD, ICA and EEMD-ICA had comparable results from a broad perspective; however ICA performed slightly better under high SNRs. From this, it is plausible that the reason why EEMD-ICA had lower log-mean squared errors at higher SNRs than EEMD alone was due to the effect of ICA. The mentioned infinite impulse response filter—which is a state-of-the-art technique—proved to be better at eliminating powerline noise than all the mentioned approaches except the proposed framework. In contrast with the other techniques, the change in performance (from a log-mean squared error point of view) plateaued rather quickly at SNRs greater than 0dB; however, it still maintained the highest performance across all SNRs evaluated. This is a strong indication that the proposed approach works relatively better at low SNRs.
A framework for the elimination of powerline noise in biomedical signals has been introduced. This adaptive method does not make assumptions on linearity, stationarity nor time-invariance and is virtually void of the need for self-correction mechanisms. Pattern recognition of the extracted features by wavelet analysis provides an enhancement to this procedure by making it fully automatable and completely unsupervised. It is worth noting that this approach is better suited for offline analysis due to its sensitivity to the time window.