Extraction for fetal ECG using single channel blind source separation algorithm based on multi-algorithm fusion

: Nowadays, detecting fetal ECG using abdominal signal is a commonly used method, but fetal ECG signal will be affected by maternal ECG. Current FECG extraction algorithms are mainly aiming at multiple channels signal. They often assume there is only one fetus and did not consider multiple births. This paper proposed a single channel blind source separation (SCBSS) algorithm based on source number estimation using multi-algorithm fusion to process single abdominal signal. The method decomposed collected single channel signal into multiple intrinsic mode function (IMF) utilizing Empirical Mode Decomposition (EMD), mapping single channel into multiple channels. Four multiple channel source number estimation (MCSNE) methods (Bootstrap, Hough, AIC and PCA) were weighting fused to estimate accurate source number and the particle swarm optimization algorithm (PSO) was employed to determine weighted coefficient. According to source number and IMF, nonnegative matrix was constructed and nonnegative matrix factorization (NMF) was employed to separate mixed signals. Experiments used single channel signal mixed by four man-made signals and single channel ECG mixed by two to verify the proposed algorithm. Results showed that the proposed algorithm could determine number of independent signal in single acquired signal. FECG could be extracted from single channel observed signal and the algorithm can be used to solve separation of MECG and FECG.


INTRODUCTION
Fetal ECG (FECG) has an important significance for fetal monitoring. The commonly used method to detect FECG is using abdominal signal. However, the obtained FECG is easy to be affected by maternal ECG (MECG) and noise and pure FECG cannot be acquired [1] [2].
In order to eliminate the interference of MECG and noise and extract accurate FECG, many scholars have proposed methods to solve this problem, such as, autocorrelation techniques, adaptive filtering and wavelet transform. Recently, blind source separation (BSS) method has been introduced into FECG extraction domain. Especially, single channel blind source separation (SCBSS) has attracted attention because of its less require for observed signal. Singular value decomposition (SVD) and analysis based on the singular value ratio (SVR) spectrum was applied, followed by an iterated application of independent component analysis (ICA) on the principle components [3]. But SVD method requires that decomposition path must be orthogonal. And the selection for separated signals needs prior knowledge. Method combined wavelet decomposi-tion with BSS algorithm independent component analysis (ICA) to extract FECG based on single-channel recordings was developed [4] [5]. But the selection of mother wavelet in wavelet analysis is needed, and mother wavelet determination is difficult with lack of prior knowledge. Reference [6] proposed to adopt ensemble empirical mode decomposition (EEMD) to decompose single-channel abdominal signal, and applied FastICA to obtain recovered FECG. But FastICA in this method is sensitive to initial value, and it could not converge if selection of initial value is not proper. Meanwhile, all above algorithms only consider one fetus, without considering multiple births, having poor practicability. This paper proposed a single channel blind source separation algorithm to separate maternal and fetal ECG. The method employs the principle of multiple channels mapping, uses Empirical Mode Decomposition (EMD) to map single channel signal into multiple intrinsic mode function (IMF) adaptively; four different source number estimation method for multiple channels are applied to estimate independent component number, and four obtained number are fused utilizing Particle Swarm Optimi-

MATHEMATIC MODEL OF SCBSS
Lathauwer had proved separation of FECG is blind source separation problem based on linear instantaneous mixed model [7]. Under such model, received single channel mixed signal x(t) is composed of N independent sources. The mathematical model can be expressed as:

SOURCE NUMBER ESTIMATION FOR SINGLE CHANNEL SIGNAL
It is commonly known that a signal is a function of time, and is thus considered dynamic. The acquisition system is a dynamic system when multiple sources are sampled by one sensor. Single channel dynamic signal can be mapped into multiple channels signal based on dynamical systems theory [8]. In this way, methods for multiple channels source number estimation can be used to effectively determine the number of single-channel mixed signals.

Multiple channel mapping based on empirical mode decomposition
EMD is a self-adaptive signal decomposition method based on local features [9]. EMD assumes that signals are composed of different and simple non-sine signals. EMD acts to decompose a signal into a series of linear, steady intrinsic mode functions (IMF) adaptively. IMF must meet the following two conditions: (1) the mean value of the envelope defined by local maxima and the envelope defined by the local minima is zero at any point; and (2) the number of extrema and number of zero crossings must either equal or differ only by one in the entire data set.
3. Extract the detail h(t)=x(t)-xm(t). 4. Set x(t)=h(t) and iterate the process to meet IMF requirements.
Employing EMD on a single channel signal may produce multiple IMF with different frequencies. Utilizing one IMF as a signal, single channel can then be mapped into multiple channels.

Source number estimation methods for multiple channels signal
Bootstrap, AIC, Hough, and PCA are MCSNE methods, considering signal characteristics in the time domain and frequency domain.

Bootstrap Method Based on Hypothesis
Testing Bootstrap method takes resample on actual data to acquire samples. It can estimate the confidence intervals of the statistics under the condition where statistic distribution is unknown, doing judgment on hypothesis testing. Bootstrap method was applied on multi-hypothesis test and a method for estimating the source number was proposed [12]. Multi-hypothesis test The process is described in detail as follows: (1) Take resample on data for B times in time domain.
(2) Calculate eigenvalues ( ) l b of covariance matrix * ( ) R b in frequency domain. Estimate eigenvalue source (useful signal or noise signal) according to the value of ( ) ( ) ( ) , 1, , (3) Set D as the confidence level of ki,j, and (4) If the significant of Hk satisfies k P k D t , assumption Hk is accepted, otherwise, refused. Then continue to test Hk-1 until estimating the source number.
Source number estimation based on Bootstrap always is combined with hypothesis testing, which will be affected by subjective factors, resulting in unstable results. There is great difference between estimated number and real source number with poor practicability.

Hough Transform Method
The next method is based on Hough transform. Sorted eigenvalues of covariance matrix for observed signal can be transformed into Hough parameter space using equation (4). Source number is determined according to detecting accumulated peaks in parameter space [13].
(2) Discretize the plane ( , ) U T and determine the Divide the Hough plane into several small squares (an equally divided part of T and U is N1 and N2 respectively,) and set AN1 N2 be 0.
(3) Use Eq. (4) to calculate ik U , which corresponds to the midpoint of small squares for each sinusoid. If (4) For k=1,…,p, repeat the step (3). The original number of source signals is obtained according to 1 2 max{ ( , )}, [1, ], The method based on Hough transform performs better than AIC with low SNR. But maximum number estimation is restricted and it cannot be more than P-3 (P is the number of array elements).

Akaike Information Criteria
Recently, the information theory approach used for source number estimation was proposed [14]. A related estimation model commonly applied is AIC expressed as follows: Where, N is sampling points, m is the number of array elements, and k is the source number to be estimated. 2k(2m-k) is the penalty function that ensures unbiased estimation. The main steps are as follows: (1) Calculate the covariance matrix S=(sij), p p of sampled data, where (2) Calculate the eigenvalues 1 2 , ,..., p (3) Calculate the value of Eq. (6) according to loglikelihood function. The optimal number k is obtained when the value of Eq. (6) reaches its minimum.
AIC compensates for weaknesses in hypothesis testing, and performs well in engineering applications. Estimation error exists in results, though, when SNR (Signal to Noise Ratio) is low.

Principal Component Analysis
PCA is used to reduce the dimensions of data [15]. Because the multidimensional vector is composed of sample characteristics that some elements is similar, elements with large variance must be identified and elements with little change must be eliminated. Thus, all remaining characteristics are useful. The main steps of the process are as follows: 1) Calculate the covariance matrix ( ) ¦ ¦ is big enough (generally larger than 90%). 4) Calculate the score of n samples on principal component r using equation (8): 1 1 ( 1, , ) PCA mainly depends on the cumulative contribution rate of the characteristics to determined source number, which means that human error will quite easily affect its accuracy.
In above four methods, Bootstrap and Hough transform methods estimate source number based on covariance in frequency domain, while AIC and PCA are based on time domain covariance. In order to make more efficient use of the characteristics in time and frequency domains, this study uses the weighted fusion method described in the next section. The advantages and disadvantages of above estimation methods have been discussed. This work proposes a new method based on multi-algorithm weighted fusion to estimate the number of sources.

Mathematical Model
After studied above four methods, we can obtain the optimal solution by optimizing the estimation error for the optimization problem The objective function is:

Particle swarm optimization algorithm
In this study, the weighted coefficient in the mathematical model above can be optimized by PSO, which regards a particle as an available solution [16], then guides the particle to its optimal position using itself and its neighboring particles' information.

A. Definition of Position and Speed
In PSO, position of every particle represents an available scheme for optimizing weights. xij is the weight that algorithm Sj have in fusion algorithm, where the sum of all weights must be equal to 1. Speed is used to calculate the probability of particle positions transformation. The speed of the particle i is V={vij}, where max max

B. Update of Position and Speed
In PSO, the position and speed are updated as follows: x v (12) where c1, c2 are learning factors, w is the weighted parameter, rand() are random numbers, and

C. Specific procedure of PSO
(1) Initialize every particle, particle number, maximum iteration number and maximum speed, randomly generate initial position and initial speed; (2) Generating a new position of each particle according to current position and speed; While (Iterations < maximum iterations) do (3) Compute fitness function value for every particle according to f0, current fitness value of particle is compared with the most optimized value and individual optimal value pb is updated; (4) According to the pb of each particle, update the group optimal value gb; (5) According to formal (11), update new speed, this is limited within maximum speed; (6) According to formal (11), update new position; End As discussed above, the proposed source number estimation algorithm for single channel based on algorithm fusion can be described as Algorithm 1.

Algorithm 1 Estimating the Number of Single Channel Source (ENSCS)
1. Input: single channel signal x(t).

Nonnegative multiple channels signal reconstruction
Following independent components estimation in a single channel signal, combined with IMF, multiple channel signals reconstruction is required for transforming underdetermined blind source separation into a well-posed problem.
The signal based on IMF and the number of independent components r is constructed, as follows: x t c c c c c r (13) As negative value will appear in EMD results, the constructed matrix is not a nonnegative matrix and cannot meet the requirements for input of NMF, thus a positive matrix must be added to allow the reconstructed matrix to become a nonnegative matrix x 2 (t), that is, 2 1 ( ) ( ) .* ( , ) x t x t positive value ones r l (14) Where r is the number of independent components, and l is the length of data.
This method only increases the amplitude of signals while not influencing signal information and results of separation.

Nonnegative matrix factorization
NMF was proposed by Lee (17) Thus, optimization problem can be transformed into nonlinear programming problem with constraints. If standard NMF algorithm is used to solve blind source separation problem without any constraint on W and H, it is difficult obtaining desired separation results. This paper employs determinant constraint [18] on basis matrix W, and applies minimum correlation constraint on coefficient matrix H [19].

A. Determinant constraint
Firstly, let us explain a definition. Let P(W) is a space which is spanned with w1,w2,…,wn. If W is square matrix, volume of P(W) can be described as: Then determinant constraint can be defined as: if the space P(W) is spanned with vector w1,w2,…,wn, and the volume vol(P) of this space is minimum, then vectors w1,w2,…,wn have uniqueness.

B. Minimum correlation constraint
When NMF is applied to solve BSS problem, coefficient matrix H represents sources in BSS problem. If source signal is not related to each other, absolute value of correlation coefficients between mixed signals is larger than source signals. That is, if mixed signals are not separated completely, every row of coefficient matrix will contain components of other source signals. Comparing completely separated signal with not completely separated signal, the absolute value of correlation coefficients for the former is less than the latter. So minimum correlation constraint can be applied to coefficient matrix H, making correlation coefficient of separated signals the least. Expression of correlation coefficient is: Determinant constraint on W and minimum correlation constraint on H are introduced simultaneously into Eq. (17), a new cost function can be obtained:    3. Estimate source number according to algorithm based on multi-algorithms fusion proposed above, the result is r; 4. Construct nonnegative matrix based on r and IMFs. Reconstruct the multiple channels signal according to equation (14), then construct nonnegative matrix x2(t) according to equation (15) Four channels signal were reconstructed according to nonnegative matrix construction principle described above. NMF was employed to separate reconstructed signals and the results are shown in Fig.4. The four channel man-made signals experiment revealed the algorithm proposed in this paper could determine source number accurately and achieve single channel signal separation with more than two channel sources.

B. ECG experiment
This section selected two channels ECG signal from MIT-BIH to verify the feasibility of proposed algorithm. The two ECG signals are added together as a single observed signal. The waveform of original and mixed signal is shown in Figure.   Estimation of number for independent component: firstly, EMD was employed to single channel signal to derive a set of IMF whose results are depicted in Fig. 7. Figure.7 a set of IMFs obtained using empirical mode decomposition. They are arranged from high frequency to low frequency, from left to right, and from top to bottom Then AIC, PCA, Bootstrap and Hough were utilized to estimate source number in single channel signal, deriving the result 2, 2, 1, 3 respectively.
At last, the four results were used to do weighted fusion to determine the accurate source number is 2 which is no difference with the experimental setup. During the weighted fusion, weighting coeffi-cient was determined by PSO.
Number estimations of AIC, PCA, Bootstrap and Hough are compared with fusion algorithm result in Ta Two channels signal were reconstructed according to nonnegative matrix construction principle described above. NMF was employed to separate reconstructed signals and the results are shown in Fig.8. The experiment revealed the algorithm proposed in this paper can be used to separate ECG signals and it can determine source number in single ECG signal, having practical application value.
Typical man-made signals and ECG signals are chosen as experimental objects to verify the effectiveness and practicability of proposed SCBSS algorithm based on source number estimation using multi-algorithm fusion. Results show that proposed algorithm could estimate source number accurately in single channel mixed signal and could solve blind source separation with two or more sources, which can be applied in engineering practices. The deficiency exists in this algorithm is construction of nonnegative matrix will affect amplitude of signals. Corresponding disposal, however, under specific amplitude requirement circumstances, may act to eliminate these defects.

CONCLUSION
A single channel blind source separation algorithm based on source number estimation using multi-algorithm fusion was proposed in this paper, which is used to solve signal mixture problem for maternal and fetal ECG exists in clinical fetal monitoring. The algorithm employed EMD to decompose single channel signal into multiple IMF, mapping single channel into multiple channels. Then AIC, PCA, Bootstrap and Hough transform were utilized to estimate source number in single channel signal, and PSO was employed to optimize weight during weighted fusion, deriving accurate number of independent components. Multiple channel signals were reconstructed and NMF was used to achieve separation of signals. Experiments utilized single channel speech whose source is two and single channel man-made signal whose source is four to verify the feasibility of proposed algorithm. Results indicate that proposed algorithm could determine source number in single channel signal accurately, and it could separate fetal ECG from single channel collected abdominal signal.
Meanwhile, the proposed algorithm is suitable for multiple births and provides a total solution for separation of maternal and fetal ECG.