This is a self-archived – parallel-published version of an original article. This version may differ from the original in pagination and typographic details. When using please cite the original. AUTHOR Elnaggar Ismail, Hurnanen Tero, Sandelin Jonas, Lahdenoja Olli, Kaisti Matti, Koivisto Tero TITLE Multichannel Bed Based Ballistocardiography Heart Rate Estimation Using Continuous Wavelet Transforms and Autocorrelation YEAR 2022 DOI http://dx.doi.org/10.22489%2FCinC.2022.364 VERSION Publisher’s pdf LICENSE Copyright in each article is held by its authors, who grant permission to copy and redistribute their work with attribution, under the terms of the Creative Commons Attribution License. Multichannel Bed Based Ballistocardiography Heart Rate Estimation Using Continuous Wavelet Transforms and Autocorrelation Ismail Elnaggar, Tero Hurnanen, Jonas Sandelin, Olli Lahdenoja, Antti Airola, Matti Kaisti, Tero Koivisto University of Turku, Department of Computing, Turku, Finland Abstract Bed based ballistocardiography (BCG) is a prime candidate for at home and nighttime monitoring especially in the growing elderly population because co-operation from the user is not required to be able to record signals. One issue with BCG is that the signal quality has intra- and inter-person variability based on factors such as age, gender, body position, and motion artifacts, making it challenging to accurately measure heart rate. A rule-based algorithm which considers all eight available BCG channels simultaneously from a given time epoch was developed using continuous wavelet transform (CWT) to extract the localized time-frequency representation of each epoch and then an averaging method was applied across the different scales of the CWT to produce a 1-dimensional array. Autocorrelation was then applied to this array to produce a heart rate estimate based on the lag between the autocorrelation maximum and the first side peak. This method does not require identification of individual heart beats to estimate heart rate and does not require annotated training data. This model produces an average mean absolute error (MAE) of 1.09 bpm across 40 subjects when compared to heart rate derived from ECG. This method produces competitive results without the need for annotated training data, which can be challenging to collect. 1. Introduction In the United States non-invasive wearables are a part of a $700 billion direct to consumer health care industry [1]. Despite the evidence pointing to the benefits of wearable health devices, only about 3.3% of users of wearable devices in the US are 65 years old or older [2]. This could be due to a number of reasons, such as lack of knowledge of wearable health technology, costs associated with such technology, or the perceived or actual difficulty in using wearables depending on an individual’s physical and mental abilities [3]. A possible alternative to wearable health technology for heart rate monitoring is non-invasive ballistocardiography (BCG). Ballistocardiography is defined as a measure of the ballistic forces generated by the heart. This is measured by the movement of the body produced by the recoil force produced from the aorta by blood passing through it. An advantage of bed based BCG is that it requires no real user input from a subject outside the initial set up of the device. This is advantageous for patient groups that suffer from diseases that affect the memory or dexterity which also happen to become more common as individuals age. This also allows for truly long term multiyear remote measurements, which medical grade wearable patches are unable to provide because of device limitations. Other options such as implantable loop recorders can be invasively implanted under the skin but are only able to take recordings of a few minutes in duration at a time [4]. Different methods of heart rate detection from BCG signals have been proposed using time domain or time- frequency domain approaches such as template matching[5], continuous wavelet transform (CWT) based peak detection [5], as well as machine learning based approaches [6]. Some of the drawbacks of these methods such as requiring annotated data to train a model [6] can be challenging when considering long-term measurements, where it may be impractical to collect enough time synced ECG data which would allow for correct beat labeling. Template matching may also prove to be challenging when we consider the fact that the morphology of a BCG waveform can vary from subject to subject as well as within one subject’s recording because of changes in body position during the course of one recording [7]. In this study, we provide a simple method that requires no training data and does not require the identification of individual heart beats from the BCG signal to accurately estimate heart rate. We believe that the proposed method could be used across different sensing configurations because it is not dependent on any dataset specific properties. This method can be used as a baseline heart rate estimate for BCG signals when developing methods that are more sophisticated when there is a lack of annotated heart rate data or external ECG data to compare heart rate estimates to. The proposed method is tested on a 40 subjects taken from an open dataset [8] against heart rate estimates produced by a single lead ECG that is time synced with the BCG signals. Computing in Cardiology 2022; Vol 49 Page 1 ISSN: 2325-887X DOI: 10.22489/CinC.2022.364 2. Methods and Materials Figure 1 gives an overview of the overall data processing pipeline used in this study to evaluate the custom bed system compared to signal lead ECG. Figure 1: Overall pipeline of proposed method 2.1. Dataset The dataset used in this study is an open access dataset from Carson et al. [8] describing a custom-made bed based BCG system. The dataset contains ECG, PPG, blood pressure waveforms as well as eight channels of BCG waveform recordings from two different types of BCG sensors from 40 subjects (17 males, mean age 34 +/- 15 years). In this study only the BCG and ECG waveforms were considered. Four electomechanical films placed under the mattress and four load cells placed under the bedposts recorded the eight channels of BCG waveforms. The participants were instructed to lay in the same supine position. More information about the bed system, such as sensor locations and data collection procedures can be found in [8]. 2.2. Signal processing The ECG and BCG signals were originally captured using a 1 kHz sampling rate. Both signals were first low- pass filtered with a 4th order Butterworth filter with a cut- off frequency value of 100 Hz before down sampling to 200 Hz. For each channel of the newly down sampled BCG signals, a 4th order Butterworth bandpass filter was applied with cut of frequencies of 5 and 35 Hz. The signals were subtracted by their mean value then normalized by dividing each channel by its standard deviation. The ECG signal was high pass filtered using a 5th order Butterworth filter with a cutoff frequency of 0.5 Hz. This was performed using the default parameters from the ecg_clean function from the NeuroKit2 python library. All the signal processing was performed in python using the open source SciPy and NeuroKit2 libraries. No other preprocessing methods such as motion artifact removal was applied to the ECG and BCG signals in order to determine the algorithms capability during noisy segments present in the recordings. 2.3. CWT Multi-Scale Averaging and Heart Rate Estimation using Autocorrelation In order to reliably extract heart rate estimates from the BCG signal, the signals were first segmented into 4-second segments to be processed individually. Figure 2 below shows a typical 4-second segment from one subject featuring BCG signal from one of the sensors used in the bed monitoring system. Figure 2: Example 4-second segment of ECG and BCG signals from dataset used in this study. After the BCG signals were filtered and segmented, a CWT using a morlet wavelet is applied individually to each channel of BCG signal. The following steps outline the algorithmic steps of producing a 1-dimensional time series array from a given CWT from one BCG signal channel. 1. Before applying the CWT, we allow the length of the given BCG signal segment to be extended one second on both sides of the segment. After application of the CWT to the now 6-second segment, the first and last second will later be discarded to reduce any boundary effects produced by the CWT. 2. The scales of the CWT are selected to represent frequencies between 5 and 35 Hz with a scale spacing of 0.1. 3. After the CWT is produced, an averaging function is applied across all scales for each time sample in the CWT. 4. A 100ms time window is centered around each time sample t and then the mean value is computed for each s-by-w matrix, where s represents all scales between 5 and 35 Hz and w represents time samples centered around t. 5. The output of this function is a 1-D array, CWT1-D, where each sample represents the mean of each s-by-w matrix considered. 6. CWT1-D is then normalized between 1 and 0 and then a gradient function is applied creating CWT1-D Gradient. 7. After which, a rolling median filter is applied to CWT1-D Gradient with a window of 100ms to further smooth the signal. The signal is re-normalized between 1 and 0 creating CWT1-D Clean Page 2 8. The first and last second of CWT1-D Clean is then discarded to produce the final 4-second segment. After step eight, the autocorrelation of CWT1-D Clean is computed resulting in an autocorrelation curve from which a heart rate estimate is found using the difference in time between the zero-point and the first side peak in the autocorrelation curve. This measure was inspired by the work done in [9] proposed earlier to extract heart rate from seismocardiographs using an autocorrelation curve. Other autocorrelation based methods have also been proposed for different BCG bed sensor systems such as [10]. Figure 3 below shows plots of the algorithm at different stages using the 4-second BCG segment in Figure 2. Figure 3: Visualization of algorithm steps to produce a heart rate estimate. From Figure 3 we see the steps taken from a filtered BCG signal segment to the resulting CWT, from which the CWT1-D Clean is extracted. The autocorrelation curve produced from CWT1-D Clean produces a heart rate estimate of 89.0 bpm from the distance shown between the red markers. For this example, the bpm estimate produced from the ECG R-peak locations in figure 2 was 90.56 bpm. This process is repeated for each of the eight BCG channels recorded by the bed sensor system to produce eight different heart rate estimates for each channel for each 4- second segment. With eight heart rate estimates per time epoch, a method is needed to select the optimal heart rate estimate for each time epoch. 2.3. Periodicity Check In order to accurately estimate the heart rate, a rule- based method is considered. First periodicity of each autocorrelation curve produced from the eight channels of BCG is checked. This is done by considering the maximum difference in the differences between peak locations detected in the autocorrelation curves represented in equation 1. 𝑀𝑎𝑥(|𝐷𝑖𝑓𝑓(𝐷𝑖𝑓𝑓( PeakLocations ))|) (1) For a given autocorrelation curve, if the result of equation 1 is greater than 20 samples at sampling frequency of 200 Hz or 100ms, the resulting autocorrelation curve is consider non-periodic, meaning that the heart rate produced by the channels autocorrelation curve cannot be considered as a reliable heart rate estimate. Figure 4 below shows an example of two auto correlation curves. Figure 4: Example of Periodic and Non-Periodic autocorrelation curves. From Figure 4, the distance between peaks varies greatly in a non-periodic segment compared to the periodic segment in green. This non-periodic segment could be caused by either poor sensor contact to the body or motion artifacts in the signal. For this given time epoch, the true heart rate estimate produced by the ECG R-peaks is 58.82. The periodic window produced a heart rate estimate of 58.25 bpm while the non-periodic window produced an estimate of 117.64 bpm. The value of 100ms was selected as a cut off because it can be seen as a representation of heart rate variability (HRV). As the target population for this type of monitoring system would be elder subjects, a 100 ms threshold was seen as reasonable because HRV tends to decrease with age and is on average less than 100 ms during sinus rhythm in elderly populations [11]. 2.4. BCG Channel Selection Figure 5 below shows the logic used to select the optimal heart rate estimate for each 4-second segment. A reference heart rate value was created by taking the median HR found within all periodic segments from the BCG channel that contained the maximal amount of periodic segments defined as BCGMaxPeri in Figure 5. This decision was made because of the relative short length of recordings across the 40 subjects, which were on average 7 minutes in length. For recordings that are multiple hours long, it is suggested that a reference HR is created by windowing the total recording into smaller segments of 5 minutes for example. This is because heart rate fluctuations can occur for example during different stages of Non-REM and REM Page 3 sleep [12]. By applying this method to shorter local windows of time, better time resolution for the estimation of heart rate can be produced when using a local reference heart rate. Figure 5: Logic used during BCG channel selection to produce accurate heart rate estimates. 3. Model Evaluation and Results The algorithm presented in section 2 was applied to all 40 subjects individually. The first second and last four seconds of each recording was omitted from analysis because of the windowing required in the proposed method. When producing the heart rate estimates from the single lead ECG recordings, an R-peak detection algorithm was selected from the NeuroKit2 library. The method selected is from Kalidas et al. [13]. This method was selected over the typical Pan-Tompkin R-peak detection method because it provided more accurate R-peak detection results when visually inspecting the ECG signals for R-peak detection errors. A heart rate estimate for each 4-seconds segment was produced by the proposed algorithm and the R-peak locations found within the given 4-second segment. The mean absolute error (MAE) was then calculated across the two heart rate estimation arrays and a result was then produced. This method produced an average MAE of 1.09 bpm across the 40 subjects with the minimum MAE being 0.33 bmp and the maximum MAE being 3.35 bpm. 4. Conclusion The proposed method produced on average a MAE of 1.09 bpm when compared to heart rate estimates produced form ECG R-peak detection. This method is significant because it does not require any annotated training data nor identification of individual heart beats to produce a reliable heart rate estimate. This method can prove to be a baseline heart rate estimator in the presence of data captured with no true reference signal such as single lead ECG. Some limitations of this study are a limited number of subjects (n=40) and a relatively short average record length of 7 minutes. References [1] Cohen et al., “Direct-to-Consumer Digital Health”, Lancet Digit. Health, vol. 2, no. 4, pp. 163-165, April 2020. [2] Wurmser, Y., “Wearables 2019 – eMarketer Trends, Forcasts, and Statistics https://www.emarketer.com/content/wearables-2019 , Retrieved August 7, 2022. [3] Kekade et al., “The Usefulness and Actual Use of Wearable Devices Among the Elderly Population”, CMPB, vol. 153, pp. 137-159, Jan. 2018. [4] Rashkovska et al., “Medical-Grade ECG Sensor for Long- Term Monitoring”, Sensors, vol. 20, no. 6, March 2020. [5] I. Sadek and B. Abdulrazak, “A Comparision of three Heart Rate Detection Algorithms over Ballistocardiogram Signals”, BSPC, vol. 70, Sept. 2021. [6] Mai et al., “Non-Contact Heartbeat Detection Based on Ballistocardiogram Using UNet and Bidirectional Long Short-Term Memory”, IEEE JBHI, vol. 26, no. 8 , pp. 3720- 3730, Aug. 2022. [7] Sadek et al., “Ballistocardiogram Signal Processing: a Review”, Health Inf Sci Syst., vol.7, May 2019. [8] Carson et al., “Bed-Based Ballistocardiography: Dataset and Ability to Track Cardiovascular Parameters”, Sensors, vol. 21, Dec. 2020. [9] Hurnanen et al., “Automated Detection of Atrial Fibrillation Based on Time–Frequency Analysis of Seismocardiograms”, IEEE JBHI, vol. 21, no. 5, pp. 1233- 1241, Sept. 2017. [10] Laurino et al., “Moving Auto-Correlation Window Approach for Heart Rate Estimation in Ballistocardiography Extracted by Mattress-Integrated Accelerometers”, Sensors, vol. 20, Sept. 2020. [11] Almeida-Santos et al., “Aging, Heart Rate Variability, and Patterns of Autonomic Regulation of the Heart”, AGG, vol. 63, April 2016. [12] Eduardo E. Benarroch, “Control of the Cardiocascular and Respiratory Systems During Sleep”, Autonomic Neuroscience, vol. 218, pp. 54-63, May 2019. [13] V. Kalidas and L. Tamil, “Real-Time QRS Detector using Stationary Wavelet Transform for Automated ECG Analysis”, BIBE, vol. 17, pp. 457-461, Oct. 2017. Address for correspondence: Ismail Elnaggar Kiinamyllynkatu 10, 20520, Turku, Finland imelna@utu.fi Page 4