Extraction of cardiac and respiration signals in electrical impedance tomography based on independent component analysis

ECG gated method is the first approach of separation of cardiac and respiration changes in electrical impedance tomography by McArdle et al [10]. It is based on two significant criteria. First, the impulse component of perfusion signal in the lung is in synchronization with cardiac changes and therefore the ECG. Second, the ventilation changes in impedance occur at a much lower rate compared to the ECG. Thus, by recording the changes of impedance during the cardiac cycle, triggered by a separated ECG module, the lung perfusion signal can be recorded. Furthermore, it has been shown in various studies that averaging the results over more than 100 cardiac cycles can reduce the interference of respiration. Averaging the measurements over 100 cardiac cycles is used in [11] and [12], where some cases averaging has been done for over 200 cardiac cycles [13]. ECG-gated acquisition demands higher frame rates from the EIT instrument to capture the changes within a cardiac cycle. New generation EIT systems have overcome these demerits. The averaging over a large number of cardiac cycles requires long duration of data acquisition and therefore may not be suitable for continuous bedside or real-time monitoring. Furthermore, as Deibele et al [14] pointed out, this averaging over a large data set could result in the fact that ‘sudden changes or irregular anomalies cannot be detected’. In addition, with this method of acquisition, the ventilation component of the signal is lost, which is certainly not the optimal solution.

These methods are based on the fact that the heart rate and the respiratory rate are very well separated in frequency. The fast Fourier transform (FFT) of the global signal can provide a good approximation of the required cut-off frequency for filtering as shown in Fig. 3 showed by Grant et al [15]. Accordingly, FIR or IIR filters can be implemented as a high-pass filter to extract the cardiac-related component from the mixed signal. Similarly, the respiratory component of the signal can be extracted using a low-pass filter, possibly with the same cut-off frequency. This filtration process can be carried out on the voltage signals or the reconstructed images [16].

Fig. 3

In spontaneously breathing subjects, the respiratory impedance changes are rarely sinusoidal, leading to frequency components across a wider range of spectrum than anticipated. Frerichs et al [17] and Grant et al [15] used frequency domain filtering to first separate components in the heart-rate domain and the respiratory rate domain, followed by linear regression fit to identify the regions of interest (ROIs). Whilst [17] used mechanically ventilated subjects to ensure no overlapping in the frequency spectrums of the ventilation and perfusion components of the acquired signal, [15] successfully separated the two components in spontaneously breathing subjects by combining the phase component of the signals with the algorithm suggested by [17].

A Principal Component Analysis (PCA) - based method for separating the respiratory and the perfusion-related signal was suggested in [14]. Their method successfully separated the two signals, allowing continuous beat-to-beat updates after a 20 s period of initiation. Time domain filtering with template functions was used to initially separate the cardiac and the ventilation signal. The results were then combined to obtain the best fit with the original signal, computed by performing least mean square (LMS) between the approximation and the actual signal.

In order to achieve proper separation of the respiratory and cardiac-related signals, the template functions for time domain filtering were identified by performing the PCA on the original data matrix. In [18] each of the aforementioned techniques compares the results achieved by each of these approaches in their experiment. Their results showed that the PCA and the frequency domain filtering methods produced highly correlated results. On the other hand, the correlation between these two methods with the ECG-gating alternative was reported considerably less. However, the correlations between each of these methods and a known standard method for perfusion imaging were not computed, rendering the results from their study inconclusive regarding as to which method is more reliable.

In electrical impedance tomography, we preprocess the acquired data to decompose them into several independent activation maps. The source for each activation map will then be localized independently [9]. To apply ICA algorithm in EIT, we have made some assumptions [9, 19, 20, 21]. At First, the sources must be independent and is valid for our multiple conductivity distribution source of localization problem. The mixture is linear. There is no delay in signal propagation from sources to sensors and finally, the total number of independent signal sources does not exceed the number of electrodes.

To apply ICA in EIT we need to do some preprocessing [9]. A useful preprocessing technique in ICA is to whiten the observed signals. Whiten refers such transformation to preliminary sphere the data where the observed vector x islinearly transformed to a new vector x˜$\tilde{x}$which is white i.e. its components are uncorrelated and the covariance matrix of x equals the identity matrix.

The second step in processing the raw EIT data x˜$\tilde{x}$is to decompose it into signal and noise subspaces by applying the PCA method in order to reduce the dimensionality of the data and remove some of its noise.

The projection of the data on the signal subspace will be referred to as V _source . The signal subspace V _source will then be decomposed into statistically independent components Vsourcej$V_{source}^{j}$. Then the EIT recordings V _source , which are the result of linear combinations of the source signals S , can therefore be expressed as:

Where M is the so-called mixing matrix. If U. V _source = S or U₌M⁻¹ , then S can be found. Here U is the weighting matrix or unmixing matrix, which

inverts the mixing process. M is unknown and the only data we have is the V _source matrix. The ICA process consists of two phases: learning and processing. During the learning phase, the ICA algorithm finds a U which minimizes the mutual information among channels. There exist several different ways to estimate the U matrix. Here the natural gradient vector is used to speed up solution convergence by avoiding the complex computation of matrix inverses [19]. Weighting matrix U can be constructed iteratively by:

where the vector p_i is defined as,

and for the nonlinear function G we used,

In the equation (3), α_i is a learning rate and I is the identity matrix. The learning rate decreases during the iteration and we stop when α_i smaller than a predefined tolerance.

The second phase of the ICA algorithm is the actual source separation. Independent components can be computed by applying the U to the signal subspace data.

Projection of independent activation maps S back onto the electrodes can be done by:

where Vsourcej$V_{source}^{j}$is the set of sensing voltage due to just the j th source. λ and R _source are the signal subspace singular values and singular vectors. For each electrode voltage just the j th source voltage map Vsourcej$V_{source}^{j}$, we can now reconstruct the different conductivity distribution in different regions using the sensitivity theorem for EIT image.

In order to test the proposed method, data was collected from a human subject. EIT data were acquired using 16 channels KHU Mark2.5 with 16 electrodes. The electrodes were placed circumferentially (equally spaced) around the thorax just below the level of the axilla. An electrical current of 1mA was injected at 50 kHz through a pair of neighboring electrodes. Differential voltages were measured between the other non-injecting pair of electrodes with neighboring pattern. Following this initial current injection, the electrical current was then injected sequentially via the next pair of electrodes and repeated until all electrodes had served for current injection. The data for one complete cycle produced a so-called “frame,” and they were saved in a raw data file for later processing. The equipment has an SNR of 80 dB.

Considering all current patterns, 256 measurements were obtained. A 2-D finite element mesh with 499 linear triangular elements was built, with shape and size similar to the cross section of the volunteer’s thorax. Using the sensitivity matrix algorithm, 1650 measurement vectors were used to generate the transition matrix. Truncated singular valued decomposition (tSVD) algorithm was applied to get the reconstruction images.

In the experimental section, some intermediate results of our algorithm stages is shown and analyzed. For our experiment we first analyzed the input signals of our EIT data. We applied some preprocessing steps like whitening, PCA for denoising and to ensure the input data to be uncorrelated. Then we applied ICA to get the independent signals from our uncorrelated raw data. In Fig. 6 we have shown the sixteen independent components (IC’s) generated from ICA algorithm. A typical input signal as can be found in the lung region is depicted. The superposition of a lower frequency ventilation signal and a much smaller higher frequency cardiac signal can be distinguished. The signals are represented as frame numbers (horizontally) along with amplitude (vertically).

Fig 6

All these IC’s are representing conductivity distribution up to some extent, which are verified by the Fast Fourier Transformation later on. Firstly, we applied FFT of all the selected IC’s generated from ICA algorithm. For our experiment, we have set the upper bound of the breathing rate which is 0.53 Hz (32 breath / min) with no modulation. The lower limit is estimated in our subjects able to breathe at rates of 0.4 Hz (24 breath / min) in most individuals. So, if we analyze the FFT of the IC’s, we will get high amplitude or peak in the region between 0.4 Hz to 0.52 Hz. From our analysis we analyze the IC’s by using FFT and observe the high amplitude of the signals. The high spectrum sustained on the above region is selected as candidate IC for our experiment. By analysis, we have selected IC_2, IC_6, IC_7 and IC_9 as our desired signal. Then we generate the reconstruction image of the IC’s and select the proper one for our desired respiration template RT. In Fig. 7 we have shown the reconstructed images of the selected IC’s. From the reconstruction images and from the rigorous survey we came to the conclusion that IC_9 represent the best respiration signal and we have chosen this independent respiration signal as RT and go further for the next step of the experiment. In Fig. 8 we have shown the magnified image of the RT and its FFT.

Fig. 7

Fig. 8

For the next step we again use the preprocessing steps and the ICA algorithm as described in Fig. 5. After the first approximation and getting the RT, now we can get the remaining signal and get the cardiac signal as our desired template. By applying ICA on the remaining signal we can get the IC’s again. From the IC’s, we choose the right one as our desired template function CT. Then we generate the reconstruction image of the IC’s and select the proper one for our desired respiration template RT. In Fig. 9 we have shown the reconstructed images of the selected IC’s (right column). Also, the remaining IC’s which are selected and used for template function as cardiac template are shown in Fig. 10.

Fig. 9

Fig. 10