Email updates

Keep up to date with the latest news and content from BMC Medical Informatics and Decision Making and BioMed Central.

Open Access Technical advance

An adaptive Kalman filter approach for cardiorespiratory signal extraction and fusion of non-contacting sensors

Jerome Foussier1*, Daniel Teichmann1, Jing Jia2, Berno Misgeld1 and Steffen Leonhardt1

Author Affiliations

1 Philips Chair for Medical Information Technology, Aachen University, Pauwelsstraße 20, 52074 Aachen, Germany

2 Philips Medizin Systeme Böblingen GmbH, Hewlett-Packard-Straße 2, 71034 Böblingen, Germany

For all author emails, please log on.

BMC Medical Informatics and Decision Making 2014, 14:37  doi:10.1186/1472-6947-14-37

The electronic version of this article is the complete one and can be found online at: http://www.biomedcentral.com/1472-6947/14/37


Received:5 November 2013
Accepted:29 April 2014
Published:9 May 2014

© 2014 Foussier et al.; licensee BioMed Central Ltd.

This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly credited. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Abstract

Background

Extracting cardiorespiratory signals from non-invasive and non-contacting sensor arrangements, i.e. magnetic induction sensors, is a challenging task. The respiratory and cardiac signals are mixed on top of a large and time-varying offset and are likely to be disturbed by measurement noise. Basic filtering techniques fail to extract relevant information for monitoring purposes.

Methods

We present a real-time filtering system based on an adaptive Kalman filter approach that separates signal offsets, respiratory and heart signals from three different sensor channels. It continuously estimates respiration and heart rates, which are fed back into the system model to enhance performance. Sensor and system noise covariance matrices are automatically adapted to the aimed application, thus improving the signal separation capabilities. We apply the filtering to two different subjects with different heart rates and sensor properties and compare the results to the non-adaptive version of the same Kalman filter. Also, the performance, depending on the initialization of the filters, is analyzed using three different configurations ranging from best to worst case.

Results

Extracted data are compared with reference heart rates derived from a standard pulse-photoplethysmographic sensor and respiration rates from a flowmeter. In the worst case for one of the subjects the adaptive filter obtains mean errors (standard deviations) of -0.2 min −1 (0.3 min −1) and -0.7 bpm (1.7 bpm) (compared to -0.2 min −1 (0.4 min −1) and 42.0 bpm (6.1 bpm) for the non-adaptive filter) for respiration and heart rate, respectively. In bad conditions the heart rate is only correctly measurable when the Kalman matrices are adapted to the target sensor signals. Also, the reduced mean error between the extracted offset and the raw sensor signal shows that adapting the Kalman filter continuously improves the ability to separate the desired signals from the raw sensor data. The average total computational time needed for the Kalman filters is under 25% of the total signal length rendering it possible to perform the filtering in real-time.

Conclusions

It is possible to measure in real-time heart and breathing rates using an adaptive Kalman filter approach. Adapting the Kalman filter matrices improves the estimation results and makes the filter universally deployable when measuring cardiorespiratory signals.

Keywords:
Adaptive Kalman filter; Sensor fusion; Heart rate; Breathing rate; Signal processing

Background

The increase in both life quality expectancy and quality of life, together with improvements in medical support, has led to an increase in the mean age of the population in developed lands. The American Administration on Aging (AoA) predicts that, compared with the year 2000, the absolute number of people aged >65 years will double by 2030 and represent 19% of the total U.S. population [1]. The increase in the elderly population will lead to additional strains on the health-care system. Therefore, there is an explicit need to transfer clinical measurement devices and technology to the patient’s home. The goal of personal health care is to relieve clinicians and the clinical infrastructure by means of technical improvements, but without reduction in diagnostic and rehabilitation performance, e.g. with telemonitoring at home [2]. This implies that new technology needs to be integrated into daily activity which, compared with a well-defined clinical environment, poses considerable challenges in terms of signal acquisition and processing. More complex processing algorithms are needed to overcome noise, artifacts and multi-sensor problems. One algorithmic approach is the use of the Kalman filtering technique.

Since its introduction in the 1960s, the Kalman filter has become a well accepted and state-of-the-art approach for many applications, especially in the technical domain. Initially, computational capacity was limited and costly, making it difficult to use Kalman filters in real-time applications. However, over the years computational power has increased, also outside the personal computer domain. Today, high performance devices are available, including microcontrollers, digital signal processors or specialized computational units such as Field Programmable Gate Arrays (FPGA). Therefore, it is now possible to move complex mathematical computations into such devices working in a self-sufficient way. Many examples of technical applications using Kalman filters in real-time have been described [3,4]. In bio-medical engineering, Kalman filters are often applied to smooth or extract physiological signals, such as respiration and cardiac activity [5,6]. Also in the domain of Electrical Impedance Tomography (EIT), the Kalman filter is able to track fast changes in impedance [7]. The Kalman filter is adequate for the present work, as it is possible to integrate prior knowledge (e.g. about sensor or system noise or state transitions). In addition, the Kalman filter is capable to denoise, separate signals or fuse sensor data, all in one architecture. Compared to other filtering or signal separation methods the Kalman filter is able to perform in real-time with very little systemic delays.

To record evaluation data a non-contact measurement technique, called magnetic induction monitoring, was chosen. This technique is a relevant method since it comprises a variety of problems typical for non-contact monitoring of vital signs. The complete employed measurement setup is described in Section ‘Sensors and measurement setup’. Afterwards, the detailed working principle and the implementation of the Kalman filter is explained (Section ‘Kalman filter’). Then, heart and respiration rate extraction results of two different subjects are shown and discussed in Section ‘Results and discussion’. We oppose the non-adaptive, described in earlier work, to the developed adaptive Kalman filter and show the performance increase, especially when no a priori knowledge about the sensor system is present. Finally, the conclusions summarize the results and gives an outlook to future work (see Section ‘Conclusions’).

Methods

Sensors and measurement setup

The signal were acquired in voluntary self-experiments of the first two authors (JF and DT). Both signed an informed consent to participate in this study. The local ethics committee decided that this type of study with this device was not in the scope of their responsibility (internal reference number EK 013/14). In addition it was stated that the self-measured data can be used for publication. Both signed an informed consent to participate in this study. Consider, that this work is not conceived to be an extensive study but rather a proof of concept. Magnetic induction monitoring is a non-contact technique for recording thoracic activity. The technique is based on magnetic coupling between the thorax and a nearby sensor-coil, as developed by Teichmann et al.[8,9]. As the coil has no conductive contact with the skin it is called a non-contact technique. The sensor-coil is driven by an alternating current and sends out an alternating magnetic field. This field induces eddy currents within conductive objects in the vicinity of the coil. In turn, these eddy currents reinduce a secondary alternating magnetic field, which affects the primary one and thereby changes the reflective impedance of the coil. If the coil is placed near the thorax, cardiorespiratory activity modulates the impedance of the coil due to motions of inner organs and the thoracic wall and/or because of changes in conductivity, e.g. caused by more or less air in the lungs or blood shifts in the heart. In this way, respiration and pulse can be easily obtained by measuring the impedance changes of the sensor-coil.

In this work, processed data were collected by three magnetic induction sensors attached to the thorax, as described in [10,11]. The sensors were located on left breast, right breast and stomach level (Figure 1). With this setting it is possible to detect motion patterns. In addition to those three sensor signals, two reference channels acquire data in parallel: one for a pulse photoplethysmographic (PPG) finger clip sensor (ChipOx, Corscience GmbH & Co. KG) the other for a thermal mass flowmeter (Model 4040, TSI). ECG and respiratory induction plethysmography (RIP) as references have been discarded as they might influence the magnetic measurement system. All measurements have been performed in resting state, so that no artifacts due to running or other movements are induced since we want to show the performance of being able to extract heart and respiration rates. The task of detecting artifacts is not part of this work, since motion patterns have been analyzed in previous work [11]. Magnetic induction monitoring is completely contact free, i.e. no mechanical, conductive or optical contact is necessary. Although this offers a considerable advantage over other monitoring methods, there are problems related to this and other non-contact monitoring techniques. A major drawback of magnetic induction monitoring is its sensitivity to other thoracic motion not related to respiration or pulse; therefore, motion artifacts are a common problem. Also, the signal content related to respiration is much higher than the cardiac-related signal content, since the thoracic conductivity distribution is more strongly affected by respiratory motion. Often, the higher harmonics of the respiratory signal may overlap the much smaller cardiac signal; in this case signal separation by common frequency filtering is not possible. To overcome these issues (at least to some extent), an adaptive Kalman filter was developed, based on previous work [5,12]. The Kalman filter system consists of three individual sensor signals (S1 to S3), that are direct inputs of the Kalman filter (Figure 1). In addition to the filter system described earlier [12], the respiratory rate and heart rate are estimated continuously and fed back to the Kalman filter where internal states and matrices are updated. Also some of the state matrices containing information about the measurement system are updated recursively and thus, automatically adapting itself to any new configuration. This procedure is described in detail in the following sections.

thumbnailFigure 1. Kalman filter sensor fusion. Kalman filter sensor fusion (sensors: S1 to S3 with flow and PPG references) and vital signs extraction with frequency adaptation.

The raw sensor data were acquired at a sampling rate of 95 Hz, which is also the processing rate of the Kalman filter. No pre-processing of the data is performed. Figure 2 presents a 30-s representative example of the employed sensor signal with the PPG and flow reference acquired in parallel.

thumbnailFigure 2. Raw signal example with parallel Flow and PPG reference recordings. Raw signal of sensor 1-3 example (top with arbitrary units [a.u.]) with parallel flow and PPG reference recordings (bottom with arbitrary units [a.u.]) (30 s window).

In summary, the Kalman filter has to deal with the following signal properties:

● A very high offset compared to the signal amplitude

● Respiration is visible in the time domain

● Heart activity is only slightly visible in the time domain (very low ratio of heart to respiration signal level)

● Good signal-to-noise ratio (SNR) for the respiratory signal, since noise is not noticeable in the time domain for this sensor

● Higher harmonics of the respiratory signal may overlap the pulse signal in the frequency domain

Kalman filter

In 1960, R.E. Kalman reported a new method for linear filtering and solving problems related to prediction [13]. Generally, the so-called “Kalman filter” consists of mathematical equations that represent an efficient way to predict a future and/or unknown state of a system, based only on the use of the preceding step. The calculations are very efficient so that they can be performed and implemented in today’s standard default personal computers, in digital signal processors (DSP) and in microcontrollers, to work in real-time applications. When process and sensor noise have time-dependent characteristics, the filter is also called a “time-varying Kalman filter”. The filter finds a prediction value <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M1','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M1">View MathML</a> with minimum variance for a disturbed state vector x[14]. The transitions from one state to another are represented in the state transition matrix A when no disturbance is present.

System and measurement noise, which are assumed to be white (zero mean) and statistically independent from each other, are represented by w and v, with their co-variance matrices Q and R (both hermitian symmetric and positive semidefinite), respectively. External disturbances (e.g. systematic errors) are fed into the system with the control vector u and the matrix B, that describe the dynamics of the disturbance. Finally, the measurement matrix H describes the integration of a real measurement into the filter procedure.

The matrices A, Q, R, B and H are generally time invariant, but may be adaptable over time. Especially when the acquired quality of the sensor data changes over time, the values of Q, R and H should take this into consideration. When the assumed system model is exposed to changes, A and B also need to be updated. All time variant matrices and vectors are denoted with a small index k.

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M2','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M2">View MathML</a>

(1)

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M3','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M3">View MathML</a>

(2)

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M4','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M4">View MathML</a>

(3)

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M5','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M5">View MathML</a>

(4)

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M6','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M6">View MathML</a>

(5)

The prediction step, also called “time update”, is described with (1) and (2). The aim of this step is to minimize the co-variance of the estimation error which represents a degree of uncertainty of the estimation [14]. Note that both equations are only valid for k>0. The initial values of <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M7','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M7">View MathML</a> and P0 have to be determined before the first iteration. The control input uk is not used in this model and is set to the zero vector, thus the matrix Bk in (1) can be ignored.

The Kalman gain Kk defined in (3) weights the innovation with respect to the measurement error co-variance matrix Rk and the estimation error co-variance matrix <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M8','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M8">View MathML</a> directly related to Qk. Higher values in Rk give more importance to the real measurement, whereas higher values for Qk put more trust in the estimation of a state. The choice of Rk and Qk is crucial for optimal filter performance. The correction step, also called “measurement update”, integrates the innovation of a new measurement zk to the estimated measurement in (4). Finally, the co-variance matrix of the estimation error Pk is updatedin (5).

The advantage of a Kalman filter is that it does not act as a pure filter, but also signal separation and fusion are realizable in a single implementation. Signal separation is done by defining several system states that are desired as separate outputs of the filter. Fusion is realized by expanding the measurement matrix H or Hk to multiple sensors. With this, the Kalman technique is applicable to a wide field of applications, especially where no exact system model is known [14].

Signal extraction and separation

The employed system was developed to work with three sensors, each measuring heart and breathing activity separately on top of a large offset (Figure 2). Therefore, the Kalman filter was split into seven states, two for heart activity (named Xf,k and Vf,k), two for respiration (named Xs,k and Vs,k) and three sensor offset estimates (named C1,k, C2,k and C3,k), merged into one estimated state vector <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M9','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M9">View MathML</a>:

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M10','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M10">View MathML</a>

(6)

In previous work, the state transition matrix A was time invariant [12]. It was designed to extract two sinusoidal-shaped signals (with different but fixed frequencies) and three offsets from three measurement sensor signals. It has been shown that the filter model for periodic motion is able to compensate (to some extent) for incorrect assumptions of the angular frequencies ωf and ωs[5]. In conditions where either the real breathing frequency or the heart rate is different from the model assumption, the estimated states tend to measure wrong frequencies. Therefore, the model needs to be updated over time by feeding back the estimated frequency into an adaptive state transition matrix which has been derived from the work of Spincemaille et al.[5]:

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M11','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M11">View MathML</a>

(7)

with ωs,k=2π·fs,k and ωf,k=2π·ff,k being the time changing angular frequencies of the breathing frequency fs,k and the heart beat rate ff,k at time step k, respectively.

The frequency measurement is based on a simple peak detector which extracts the location of maxima and minima in the signal and then computes the mean interval length between maximal and minimal points. The main advantage of a peak detector is the reduced computational complexity compared to other frequency measurement methods, e.g. based on the spectrum analysis. To accurately measure respiration and heart rate, the algorithm needs to buffer the output signal of the Kalman filter, where the denoted “Breathing” and “Heart” feedbacks in Figure 1 are represented by the internal Kalman states Xs,k and Xf,k, respectively. Note that it is not possible to estimate the heart rate directly from the raw sensor signal, due to the small heart signal amplitude and a possible frequency overlap of respiratory harmonics. A buffer of 20 s (approx. 1900 samples) for respiration and of 10 s (approx. 950 samples) for the heart signal were chosen to guarantee at least two consecutive respiratory cycles and several heart beats within the buffer. However, having such a large buffer causes systematic delays in the frequency estimation of half of the buffer length, i.e. 10 s and 5 s, respectively. Generally, because heart and respiration rates change slowly, we do not need to adapt the system model at 95 Hz. In fact, only at every tenth sample a frequency measurement is performed, resulting in an effective estimation rate of 9.5 Hz reducing computational effort by 90%. The remaining nine samples are filled with the last valid estimation to keep the same sampling rate of 95 Hz as the rest of the Kalman filter. This high sampling rate is needed to keep systematic errors at a minimum. At 95 Hz, the maximum absolute time resolution to measure the time between two signal peaks is <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M12','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M12">View MathML</a>, inducing a measurement error increasing with the heart rate with 0.0175% per beat per minute (bpm) (e.g. at 80 bpm the relative error is up to 1.4% equivalent to an absolute error of 1.12 bpm). The frequency measurement procedure has been validated with simulated sinusoidal signals covering heart rates from 60 bpm to 120 bpm and respiration rates from 12 min−1 to 15 min−1.

For stability reasons, abrupt changes in the frequency measurement need to be avoided, e.g. false estimations of the frequencies or artifacts, directly fed into the state transition matrix Ak of the adaptive Kalman filter. Therefore, the estimated frequencies representing the heart and the respiration rates are each lowpass filtered with first order infinite impulse response (IIR) Butterworth filters, designed with the Matlab R2011a (MathWorks) “Signal Processing Toolbox” with cutoff frequencies of 0.1 Hz and 0.05 Hz, respectively.

In our previous work [12], the measurement noise co-variance matrix R was also time-invariant, fixed to the standard deviation of each sensor signal determined before the real measurement. Therefore, no continuous adaptation to changes in the sensor signal quality was possible. Since such adaptation is important for optimal performance and, generally, this type of a priori calibration is inconvenient for our work, the measurement noise was also estimated leading to a time variant matrix Rk:

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M13','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M13">View MathML</a>

(8)

with <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M14','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M14">View MathML</a> as ideal measurement variances of sensor 1-3 at time step k, respectively. Generally, the ideal value is unknown and different for each sensor. So, at each discrete time step k, the noise co-variance matrix Rk has to be updated by estimating the standard deviation (std) of a numerically differentiated (diff) short signal segment. In our case the last 0.5 s (approx. 48 samples) si,Δ of each sensor i have been taken into account:

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M15','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M15">View MathML</a>

(9)

Note that all estimations that are described in this work have been computed every tenth sample and lowpass filtered in the same way as the heart rate estimation value as described previously (first order Butterworth filter with 0.1 Hz cutoff frequency). Also, the development of all the described estimations first has been applied to simulated signals, where noise and standard deviations of the heart and respiratory signals were known.

Signal filtering

Spincemaille et al. have shown [5] that same ratios of Rk and Qk also give similar filtering performances. So, changing the Rk matrix also means adapting the system noise co-variance matrix Qk appropriately:

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M16','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M16">View MathML</a>

(10)

with variances of the long-term trend <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M17','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M17">View MathML</a> of sensors 1-3.

The co-variance matrix Q0 is initialized with the heart and breathing frequencies ωf,k and ωs,k and the trend variations <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M18','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M18">View MathML</a> initialized by the user at time step k=0, respectively. We will show, that those values do not need to be very accurate, since they are automatically estimated over time. However, as the exact trend variance is unknown it is also estimated using the raw sensor signal:

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M19','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M19">View MathML</a>

(11)

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M20','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M20">View MathML</a>

(12)

with <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M21','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M21">View MathML</a> as standard deviation estimation of the respiration using the long-term segment li,Δ over the last 20 s (approx. 1900 samples). We defined that the trend variation corresponds to 1% of the respiratory variation. Lower entries in the matrix (10) mean that more trust is placed on the predicted model states than on the raw measurement. Higher frequencies also cause higher amplitudes in the differentiated Kalman state vector signals Vf,k and Vs,k and are proportional to <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M22','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M22">View MathML</a> and <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M23','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M23">View MathML</a>, respectively. A high long-term variance in the sensor signal is assumed to be caused by respiration and heart beat. The values for the states Xs,k and Xf,k defined in eq. (6) normally are smaller than for Vs,k and Vf,k. With this, more trust is given to the system model than to the measurement resulting in smoothed sinusoidal-shaped signals with angular frequencies of around ωs,k and ωf,k, respectively. This simplifies the detection of peaks in the signal and improves frequency estimation. In fact, the states Vs,k and Vf,k contain more actual sensor information, thus more noise and possible artifacts are present; these states would need a processing step before being used for frequency estimation.

Sensor fusion

In general, we assumed that the sensor configuration is able to measure vital signs. The first four columns of the measurement matrix Hk, representing the heart and respiration states of the estimated state vector <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M24','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M24">View MathML</a>, defined in (6), are updated by using the noise estimation previously described in (9), (12) and (15):

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M25','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M25">View MathML</a>

(13)

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M26','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M26">View MathML</a>

(14)

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M27','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M27">View MathML</a>

(15)

with <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M28','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M28">View MathML</a> and <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M29','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M29">View MathML</a> as sensor specific scaling factors for measuring heart and respiratory activity. This is the only a priori information that has to be defined by the user. Generally, the general performance of the filter is not very sensitive to this scaling as the Kalman filter itself adapts to the average of all three sensors. Nevertheless, it might enhance the filtering performance. In our case we know that sensor 1 mainly contains cardiac information since it is attached right above the heart (see Figure 1) and all other sensors are located far away from the heart. Those are scaled to 10% as well as the respiratory data content of sensor 1:

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M30','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M30">View MathML</a>

(16)

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M31','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M31">View MathML</a>

(17)

If no a priori information is given, e.g. when the sensor location is unknown, the scalf and scals settings should be chosen as (1 1 1)T to give the same scaling to all sensors.

The last three columns in Hk indicate which offset in the estimated state vector <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M32','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M32">View MathML</a> is linked to a specific sensor, finally leading to:

<a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M33','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M33">View MathML</a>

(18)

However, the offset links cannot be changed as they act totally independently. In (3) and (4) the sensor fusion is performed, as all measurement channels are joined into the estimated state vector <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M34','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M34">View MathML</a> with the Kalman Gain Kk.

Results and discussion

In this section, the performance of the non-adaptive and the adaptive Kalman filter is evaluated by comparing the estimated respiration rates and heart rates with the reference flow and PPG signals. Instead of the peak detection and filtering methods described in Section ‘Signal extraction and separation’ we applied a spectrum based frequency estimator to the reference flow and PPG signals to obtain much higher frequency resolution. This estimator is not adequate for real-time application since it is computationally much more complex. Note that in the first 22.5 s (i.e. the minimal length of the buffers used for the respiration frequency and signal variance estimations), frequency estimation neither for the respiration rate nor for the heart rate is computed. Especially for the adaptive Kalman filter, it is important to regulate the frequencies in the model only when both frequencies are really measurable thus only when the Kalman filter is settled. This avoids the filter reaching an irreversible incorrect state. Therefore, during further analysis, the first 22.5 s are ignored as generally no valid and comparable data are available.

For the performance analysis between the adaptive and non-adaptive Kalman filter, we defined three different settings:

● SE: When the Sensor Estimation settings is used, the whole raw sensor signal is first analyzed, giving the optimal estimations of variances. Note that this setting only works in a post-processing way, but it is used to create the best-case setting in this work.

● DS: The Default Settings should adequately work in most cases.

● BS: Bad Settings that give a very inappropriate definition of the system - the worst-case scenario when the user misconfigures the Kalman filter.

In Table 1 the initialization parameters for each setting and subject are shown. Note that for BS, especially the noise variances have been set to very high and the respiration/heart variances to relatively low values. This setting reflects a possible misconfiguration of the Kalman filter, e.g. when wrong a priori information are set by the user. With BS we want to demonstrate that the user does not have to care about the correct settings when using the adaptive Kalman filter contrarily to the non-adaptive filter. The SE setting, so the optimal settings determined by preprocessing the raw acquisition data, shows that the sensor amplitudes between both subjects are very different making it very difficult to estimate correct parameters from case to case. The ratios between the estimated heart and respiration deviation <a onClick="popup('http://www.biomedcentral.com/1472-6947/14/37/mathml/M35','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1472-6947/14/37/mathml/M35">View MathML</a> confirm the fact that sensor 1 contains most of the heart signal information with respect to the respiratory signal compared to all other sensors.

Table 1. Initialization parameters using different settings for subject 1 and 2

The available signal segment lengths for subject 1 are 144.8 s and 163.4 s for subject 2. The results of the frequency estimations of heart and respiration over those segments are shown in Figure 3, Figure 4 and Figure 5 after having applied the procedures described in Section ‘Kalman filter’ with the adaptive and non-adaptive Kalman filter for the SE, DS and BS settings, respectively. The configurations SE and DS show very good frequency estimation performance for both Kalman filters, the adaptive and non-adaptive version, and for both subjects. The estimated rates follow the reference lines very well. To support the visual analysis, the mean difference and standard deviation between the reference and the estimated frequency have been computed which are listed for each setting and subject in Table 2. Settings SE and DS generate errors that are all in the same range, showing almost equal performances between both Kalman filters. The standard deviations are slightly smaller for the adaptive filter. The mean deviations are within the systematic error range resulting from the 95 Hz sampling rate as described in Section ‘Signal extraction and separation’. However, the BS setting gives a totally different result. The non-adaptive Kalman filter is not able to extract the heart rate any more but since the respiration signal is much higher in amplitude than the heart signal, even this bad configuration is usable for respiration rate estimation. The adaptive Kalman filter compensates the bad settings by adapting to the sensor recursively and allowing it to extract heart rate estimations much more accurately. Table 2 also clearly confirms this finding where the heart signal mean error for the adaptive Kalman filter is about sixty times and the standard deviation about four times smaller than for the non-adaptive Kalman filter.

thumbnailFigure 3. Heart and respiration rates determined by adaptive and non-adaptive Kalman filter with PPG and flow reference using the SE setting. Heart (top) and respiration (bottom) rates determined by adaptive (red) and non-adaptive (blue) Kalman filter with PPG and flow reference (gray) using the SE setting for subject 1 (left column) and 2 (right column). The dotted gray lines indicate the starting frequencies of the heart and respiration rates defined in the setting at startup.

thumbnailFigure 4. Heart and respiration rates determined by adaptive and non-adaptive Kalman filter with PPG and flow reference using the DS setting. Heart (top) and respiration (bottom) rates determined by adaptive (red) and non-adaptive (blue) Kalman filter with PPG and flow reference (gray) using the DS setting for subject 1 (left column) and 2 (right column). The dotted gray lines indicate the starting frequencies of the heart and respiration rates defined in the setting at startup.

thumbnailFigure 5. Heart and respiration rates determined by adaptive and non-adaptive Kalman filter with PPG and flow reference using the BS setting. Heart (top) and respiration (bottom) rates determined by adaptive (red) and non-adaptive (blue) Kalman filter with PPG and flow reference (gray) using the BS setting for subject 1 (left column) and 2 (right column). The dotted gray lines indicate the starting frequencies of the heart and respiration rates defined in the setting at startup.

Table 2. Mean Error (standard deviation) between flow/PPG reference and estimated respiration/heart rates for the non-adaptive and adaptive (Ad.) Kalman filter

Another performance indicator is the ability to separate the sensor signal into an offset, heart and respiration signals adequately. To evaluate the offset estimation performance we compute the mean error between the offset and the raw sensor signal shown in Table 3. High mean errors imply a high deviation from the original signal and that offsets are not correctly estimated. Note that the error can only be compared within one sensor and subject, since the sensors are not calibrated to the same amplitude range. In all settings most of the mean error and standard deviation values of the adaptive Kalman filter are under or at least very near to the non-adaptive Kalman filter, indicating that it outperforms the non-adaptive filter for offset estimation. To give an insight in the heart and respiration signal separation capabilities we show 20 s of the Xf and Xs values in the middle of the acquired signal in Figure 6, Figure 7 and Figure 8 for the three settings. It is expected that Xf and Xs follow a sinusoidal shaped signal with angular frequencies of ωf and ωs, respectively. For Xs, the signals between the adaptive and non-adaptive Kalman filter resemble through all settings as expected. On the contrary, we see a large amount of respiration on top of the heart rate signal in the non-adaptive Xf which is not present in the Xf of the adaptive filter, regardless of the used configuration. This interaction of Xf with the respiratory signal degrades the performance of detecting the peaks that are needed for heart rate estimations. For BS, the Xf of the non-adaptive filter does not contain any visible heart signal any more, making it impossible to estimate heart rates with the employed peak detection procedure.

thumbnailFigure 6. Internal statesXf andXs determined by adaptive and non-adaptive Kalman filter with PPG and Flow Reference using the SE setting. Internal states Xf (top) and Xs (bottom) determined by adaptive (red) and non-adaptive (blue) Kalman filter with PPG and Flow Reference (gray) using the SE setting for subject 1 (left column) and 2 (right column).

thumbnailFigure 7. Internal statesXf andXs determined by adaptive and non-adaptive Kalman filter with PPG and Flow Reference using the DS setting. Internal states Xf and Xs determined by adaptive (red) and non-adaptive (blue) Kalman filter with PPG and Flow Reference (gray) using the DS setting for subject 1 (left column) and 2 (right column).

Table 3. Mean Error (standard deviation) between estimated offsetC1,k- C3,kand raw sensor signals S1-S3 for the non-adaptive and adaptive (Ad.) Kalman filter

thumbnailFigure 8. Internal statesXf andXs determined by adaptive and non-adaptive Kalman filter with PPG and Flow Reference using the BS setting. Internal states Xf and Xs determined by adaptive (red) and non-adaptive (blue) Kalman filter with PPG and Flow Reference (gray) using the BS setting for subject 1 (left column) and 2 (right column).

To assess whether the described procedure is able to work in real-time, computational time measurements were conducted. The filter procedure was implemented in Matlab R2011a (MathWorks) on a personal computer with a 3.10 GHz dual core processor and 4 GB RAM. In total the adaptive Kalman filter (including all initialization processes and filtering, but excluding reference frequency measurements) in average 35.9 s (24.79% of the total 144.8 s) and 40.5 s (24.80% of the total 163.4 s) processing times are required, respectively for subject 1 and 2. Without adaptation, the procedures take in average 24.2 s (16.73% of the total 144.8 s) and 24.7 s (15.08% of the total 163.4 s). Hence, performing one regular Kalman filter step requires about 1.7 ms out of 10.5 ms at 95 Hz sampling rate. The effective adaptation rate is reduced to 9.5 Hz (every tenth sample equal to 105 ms sampling time) and in average about 50-60% more computational time is needed, i.e. approximately 2.6 ms per sample. Note, that the Kalman filter computations still are at 95 Hz sampling rate. The surplus of 9 ms within the 105 ms window (ten samples) compared to the non-adaptive filter only occurs during the first adaptation sample which takes 1.7 ms+9 ms=10.7 ms in total. This is slightly more than the sampling time of 10.5 ms resulting in a maximum lag of one sample every tenth sample. Therefore we conclude that the intrinsic computations of both filters are real-time compliant for separating respiratory and cardiac signal activity. However, as described in Section ‘Signal extraction and separation’, buffer sizes of 20 s and 10 s are still necessary to estimate breathing and heart rates, respectively. Note that, until now, no special code optimization has been performed, which would further reduce the required computational time and make it possible to integrate the procedures into microcontrollers or digital signal processing units.

Although the signals were acquired on healthy young men in this proof of concept, a target application might be the telemonitoring of elderly at home. The main advantage of this technique is the contactless and unobtrusive way the measurement are performed. The sensors may be arranged in a shirt, a bed or in a chair just to give some examples. Nevertheless, it has to be examined how daily activities (e.g. cooking, vacuuming, mowing the lawn or walking around) affect the performance of detecting breathing and heart rates. The main drawback of the magnetic induction measurement method is that movements relative to the coils cause large artifacts in the signal on top of the desired signal. These may be detected or in a good case even be compensated. We do not expect that the performance is age dependent. Rather the position of the sensors is more crucial for optimal performance. During resting activities (e.g. watching TV, sitting in a chair reading or sleeping) where the number of body movements is small, the developed filter should be able to perform as well as described since it adapts automatically to the measurement conditions.

Conclusions

The above analysis shows that estimations of respiration and heart rate based on the implemented adaptive Kalman filter perform very well in a real-time acquisition scenario employing contactless sensors measuring cardiorespiratory signals. The direct evaluation of the respiration and heart rate only based on the Kalman filter states and the direct feed-back into the time variant state transition matrix Ak improved the results compared with the non-adaptive procedure. It is also important that the filter does not require a specific signal shape as long as it contains periodic content for respiration and heart activity; therefore, the filter is suitable for all biomedical signals containing two periodic motions on top of large signal offsets. The delays caused by the buffer windows of 20 s for measurement of the respiration rate and of 10 s for the heart rate are necessary to correctly estimate the frequencies with the described peak detection method. The frequency measurement resolution of 0.0175% per bpm inherently generates errors, especially with increased heart rates. For shorter delays and a better frequency resolution, the frequency estimation method needs to be improved, e.g. by up-sampling the signal in the buffer or employing more complex detection algorithms but without neglecting the real-time ability.

In general, this implementation of the adaptive Kalman filter is kept as simple as possible; this allows to port the code to embedded processing units, e.g. microcontrollers or digital signal processors. It is also possible to expand the system to more than three sensors by adding one row and one column in the state transition matrix Ak, the system noise co-variance matrix Qk and the measurement noise co-variance matrix Rk and only one row in the measurement matrix H.

Since this analysis is a proof of concept, a larger study with additional subjects needs to be performed to validate all of the above methods.

Competing interests

The authors declare that they do not have any competing interests.

Authors’ contributions

DT and JJ developed the employed hardware measurement system and carried out the data acquisition. JF designed the signal analysis and processing with the support of BM and SL, who contributed with their long lasting knowledge in feedback control and automation. All authors carefully revised the manuscript and contributed to the development of the described system and software. All authors read and approved the final manuscript.

Authors’ information

Jérôme Foussier was born in Cologne, Germany, in 1984. He holds a Dipl.-Ing. degree in Electrical Engineering from RWTH Aachen University, Germany. Currently, he is pursuing the Dr.-Ing. (Ph.D.) degree at the Chair of Medical Information Technology, RWTH Aachen University, where he is also working as a Research Assistant. His research interests include signal processing and classification as well as physiological measurement techniques.

Daniel Teichmann was born in Essen, Germany, in 1982. He holds a Dipl.-Ing. degree in Electrical Engineering from RWTH Aachen University, Germany. Currently, he is pursuing the Dr.-Ing. (Ph.D.) degree at the Chair of Medical Information Technology, RWTH Aachen University, where he is also working as a Research Assistant. His research interests include non-contact monitoring techniques and signal processing.

Jing Jia was born in Shanghai, China, in 1985. She holds a Dipl.-Ing. degree in Electrical Engineering as well as in Economic Science from RWTH Aachen University, Germany. Currently, she is working at ’Philips Medizin Systeme Böblingen GmbH’ in Böblingen, Germany.

Berno J.E. Misgeld received the Dipl.-Ing. (FH) in Electrical and Automation Engineering from University of Applied Sciences, Aachen, Germany, and the M.Sc. degree from Coventry University, Coventry, U.K., in 2003, respectively. In 2007 he received the Dr.-Ing. degree in Biomedical and Control Engineering from Ruhr-University Bochum, Bochum, Germany. From 2006 to 2011 he was a research and development engineer for guidance and flight control systems at Diehl-BGT-Defence, Ueberlingen, Germany. Since 2011 he is a Senior Scientific Engineer in Biomechatronical Systems and Rehabilitation Robotics at the Chair of Medical Information Technology at RWTH Aachen University, Aachen, Germany. His research interests include feedback control and filtering with application to biomedical systems, robotics and medicine.

Steffen Leonhardt was born in Frankfurt, Germany, in 1961. He holds an M.S. in Computer Engineering from SUNY at Buffalo, NY, USA, a Dipl.-Ing. in Electrical Engineering and a Dr.-Ing. degree in Control Engineering from the Technische Universität Darmstadt, Germany, and a M.D. in Medicine from J. W. Goethe University, Frankfurt, Germany. He has 5 years of R&D management experience in medical engineering industry and was appointed Full Professor and Head of the Philips endowed Chair of Medical Information Technology at RWTH Aachen University, Germany, in 2003. His research interests include physiological measurement techniques, personal health care systems and feedback control systems in medicine.

Acknowledgements

We would like to thank Laraine Visser-Isles for her insightful comments and the language review of the manuscript.

References

  1. Administration on Aging: Aging statistics - department of health & human services.

    2012.

    [http://www.aoa.gov/aoaroot/aging\_statistics/index.aspx webcite]

  2. Meystre S: The current state of telemonitoring: a comment on the literature.

    Telemed J e-Health 2005, 11:63-69. Publisher Full Text OpenURL

  3. Moreno V, Pigazo A (eds.): Kalman Filter: Recent Advances and Applications. Rijeka (Croatia): InTech; 2009. OpenURL

  4. Chui CK, Chen G: Kalman Filtering with Real-Time Applications. Springer Ser. Info. Sci., Vol. 17. Berlin-Heidelberg (Germany): Springer; 2009. OpenURL

  5. Spincemaille P, Nguyen TD, Prince MR, Wang Y: Kalman filtering for real-time navigator processing.

    Magn Reson Med 2008, 60:158-168. PubMed Abstract | Publisher Full Text OpenURL

  6. Sayadi O, Shamsollahi MB: ECG denoising and compression using a modified extended Kalman filter structure.

    IEEE Trans on Biomed Eng 2008, 55(9):2240-8. OpenURL

  7. Vauhkonen M, Kaipio JP, Karjalainen Pa: A Kalman filter approach to track fast impedance changes in electrical impedance tomography.

    IEEE Trans on Biomed Eng 1998, 45(4):486-493. Publisher Full Text OpenURL

  8. Teichmann D, Foussier J, Leonhardt S: Respiration monitoring based on magnetic induction using a single coil. In Biomedical Circuits and Systems Conference (BioCAS), 2 . Paphos, Cyprus: IEEE; 2010:37-40. OpenURL

  9. Teichmann D, Foussier J, Jia J, Leonhardt S, Walter M: Noncontact monitoring of cardiorespiratory activity by electromagnetic coupling.

    IEEE Trans on Biomed Eng 2013, 60(8):2142-2152. OpenURL

  10. Teichmann D, Foussier J, Buscher M, Walter M, Leonhardt S: Textile integration of a magnetic induction sensor for monitoring of cardiorespiratory activity. In IFMBE Proceedings of World Congress on Medical Physics and Biomedical Engineering, May 26-31, 2012, Beijing, China, Volume 39. Berlin-Heidelberg (Germany): Springer; 2013:1350-1353. OpenURL

  11. Teichmann D, Kuhn A, Leonhardt S, Walter M: Human motion classification based on a textile integrated and wearable sensor array.

    Physiol Meas 2013, 34(9):963-975. PubMed Abstract | Publisher Full Text OpenURL

  12. Foussier J, Teichmann D, Jia J, Leonhardt S: Fusion of non-contacting sensors and vital parameter extraction using Kalman filtering. In Proceedings of the World Congress on Engineering 2011, WCE 2011, July 6 - 8 2011, London, U.K., Volume II. Hong-Kong (China): IAENG; 2011:1592-1596. OpenURL

  13. Kalman R: A new approach to linear filtering and prediction problems.

    Trans ASME-J Basic Eng 1960, 82:35-45. OpenURL

  14. Welch G, Bishop G: An introduction to the Kalman filter. Tech. rep1995

Pre-publication history

The pre-publication history for this paper can be accessed here:

http://www.biomedcentral.com/1472-6947/14/37/prepub