Next Article in Journal
Heterogeneous Iris One-to-One Certification with Universal Sensors Based On Quality Fuzzy Inference and Multi-Feature Fusion Lightweight Neural Network
Next Article in Special Issue
DataSpoon: Validation of an Instrumented Spoon for Assessment of Self-Feeding
Previous Article in Journal
High Sensitive Biosensors Based on the Coupling Between Surface Plasmon Polaritons on Titanium Nitride and a Planar Waveguide Mode
Previous Article in Special Issue
A Systematic Approach to the Design and Characterization of a Smart Insole for Detecting Vertical Ground Reaction Force (vGRF) in Gait Analysis
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Photoplethysmographic Time-Domain Heart Rate Measurement Algorithm for Resource-Constrained Wearable Devices and its Implementation

by
Marek Wójcikowski
* and
Bogdan Pankiewicz
Faculty of Electronics, Telecommunications and Informatics, Gdańsk University of Technology, 80-233 Gdańsk, Poland
*
Author to whom correspondence should be addressed.
Sensors 2020, 20(6), 1783; https://doi.org/10.3390/s20061783
Submission received: 17 February 2020 / Revised: 20 March 2020 / Accepted: 21 March 2020 / Published: 23 March 2020
(This article belongs to the Special Issue Low-Cost Sensors and Biological Signals)

Abstract

:
This paper presents an algorithm for the measurement of the human heart rate, using photoplethysmography (PPG), i.e., the detection of the light at the skin surface. The signal from the PPG sensor is processed in time-domain; the peaks in the preprocessed and conditioned PPG waveform are detected by using a peak detection algorithm to find the heart rate in real time. Apart from the PPG sensor, the accelerometer is also used to detect body movement and to indicate the moments in time, for which the PPG waveform can be unreliable. This paper describes in detail the signal conditioning path and the modified algorithm, and it also gives an example of implementation in a resource-constrained wrist-wearable device. The algorithm was evaluated by using the publicly available PPG-DaLia dataset containing samples collected during real-life activities with a PPG sensor and accelerometer and with an ECG signal as ground truth. The quality of the results is comparable to the other algorithms from the literature, while the required hardware resources are lower, which can be significant for wearable applications.

1. Introduction

Advancements in modern technologies enabled field monitoring of some parameters of human health [1]; for example, heart monitoring is used for off-hospital monitoring and in fitness and professional sport activities. The most commonly measured value is the heart rate (HR), although advanced applications also use other values, e.g., pulse irregularity, as well as biometric identification or analysis of accurate electrical signals that cause heart contraction, i.e., electrocardiography (ECG) [2]. Accurate ECG requires connecting electrodes to the patient’s body in several different places, which is inconvenient for the patient, and it can be used only in certain situations. A much more convenient method is measuring the pulse on the wrist by using photoelectric methods. The skin of the wrist is irradiated with single or multicolor light, and then the reflected light is measured. The intensity of the reflected light depends on the absorption of the skin, which depends on the blood volume supplied to the tissues. In this way, the received signal contains information about the current blood supply to the vessels near the measuring device. This method, introduced by Hertzman [3], is known as photoplethysmography (PPG). Unfortunately, PPG signals obtained from a moving person’s wrist are weak, distorted, and contain noise. The noise level is often higher than a usable PPG signal. Correct analysis of a low-quality PPG signal is a very challenging task and can consume significant processing time, energy, and resources.
Heart rate estimation from wrist PPG is now a popular area of investigation, and many algorithms for such a task have been proposed [4]. Some of them require significant computing power and memory usage, blocking their application in small portable low-power devices. There are two main approaches used: time-domain sophisticated filtering and frequency-domain processing. Both are often accompanied by a movement sensor (accelerometer) for movement-based artefacts/spectrum removal.
A straightforward approach toward heart rate estimation is peak detection in a periodic signal. One of the simplest possibilities is to use threshold or auto-threshold values in the signal time window [5,6]. Another way is to use transforms such as continuous wavelet [7] or Hilbert [8]. In [9], a nonlinear filter bank was used, with variable cutoff frequencies. In [10], the detection algorithm was based on a time-varying autoregressive model, with a Kalman filter used for autoregressive parameter estimations. Hidden Markov models were used in [11] for combining structural and statistical knowledge of the signal in a single parametric model. A neural network adaptive whitening filter to model the lower frequencies of the signal is presented in [12]. Much recent work is concentrated on methods using time-frequency spectra [13,14] and deep-learning approaches [15,16,17]. A review of the current state-of-the-art signal-processing techniques for HR estimation from a wrist PPG signal can be found, for example, in [4,18].
The proposed solution of time-domain heart rate measurement algorithm (TDHR) consists of three main blocks: signal conditioning, peak detection, and heart-rate-measuring blocks. In this work, the modified Automatic Multiscale-based Peak Detection (AMPD) algorithm from [19], together with a bandpass filter/limiter, was used for finding the HR from a wrist-based PPG signal. All of the signal processing was done while taking into account the need for low power, low resources, and computing power utilization, necessary for a self-sufficient mobile sensor. The main contributions of this paper are as follows:
  • A two-stage input-signal-conditioning digital nonlinear filtering block with limiter;
  • Application of the AMPD algorithm for HR peak detection;
  • Modification of the AMPD algorithm toward efficient implementation in low-power resource-constrained hardware;
  • The proposition of a time-domain heart-rate-measuring algorithm with an accelerometer-based false measurement removal.
This paper is organized as follows. In Section 2, the signal conditioning block is described, which is followed by the peak detection algorithm described in Section 3 and the heart rate calculation block presented in Section 4. The proposed algorithm was evaluated by using a multi-hour dataset; the results of the evaluation and the comparison to the other solutions from the literature are presented in Section 5. Section 6 contains the details of the implementation in a low-power wearable device. Finally, conclusions are provided in Section 7.

2. Signal Conditioning

PPG-signal measuring is usually done with an infrared (IR), red, or green light-emitting diode (LED) as the light source and a photodetector (PD) receiving the reflected or transmitted light. Often, instead of one, a few LEDs and/or PDs are used, providing the possibility to choose the best observed signal or to use additional preprocessing, such as, for example, averaging. Two measure modes can be used: reflectance and transmission. In the first case, the LED and PD are placed on the same skin surface, close to each other; the typical spacing between the PD and LED is in the range of 5–15 mm. The PD measures the light reflected from the tissue. In transmission mode, the LED and PD are placed on the opposite sides of a human body part, and the light reaching the PD goes through the whole body part. The commonly used places on a human body for such PPG devices are the fingertip, wrist, earlobe, forehead, torso, ankle, and nose [18,20]. The wrist is the most convenient place to mount the monitoring device for everyday life, but this is not optimal regarding the strength of the signal. Due to the relatively large thickness and the presence of bones, only the reflectance method of measurement is practical on the wrist. In Figure 1, the waveforms of the signal received from the sensor in reflectance mode on the index finger and on the outer wrist, where a watch is usually worn, are presented. A useful signal carrying information about the pulse is present in the form of the peaks with the period of about 30–40 samples superimposed on the curve. As shown in Figure 1, the amplitude of these peaks for the signal measured at the finger is about 200, while the amplitude measured at the wrist is much lower, about 30–50, while the average signal level (baseline) is about 11,100.
As can be seen from the waveforms from Figure 1, the signal from the wrist is much weaker; moreover, each signal has a variable offset, and it contains noise and interference. To extract the heart rate in a time-domain, the signal must be preprocessed, so that the peak detection algorithm can find the peaks. The simplest approach is to use a band-pass filter that can help to reduce the noise and eliminate the constant component of the signal. Such a solution is in common use [4]. However, the signal at the input of the filter can have significant changes of the constant component, as well as signal fluctuations caused by wrist movement, resulting in large peaks and oscillations at the output of the band-pass filter, as can be seen in Figure 2a,b.
The proposed approach for PPG signal acquisition and processing is presented in Figure 3. The raw signal obtained from the PD is fed, in the first step, to the band-pass biquad section, with an internal limiter. The output value of the standard biquad section using direct version I is given by Equation (1):
y i = 1 a 0 ( b 0 x i + b 1 x i 1 + b 2 x i 2 a 1 y i 1 a 2 y i 2 ) ,
where yi is i-th output sample; xi is i-th input sample; a0, a1, and a2 and b0, b1, and b2 are, respectively, the denominator and nominator coefficients of the biquad transfer function. The biquad section with the internal limiter works in two steps: First, it calculates candidate yC,i for the output value, according to Equation (1), and then it uses Equation (2) to update the actual output value:
y i = { y C , i   if   L L y C , i L H L L   if   y C , i < L L L H   if   y C , i > L H ,
where LL and LH are the limiter’s parameters. These parameters should be selected so that the PPG signal from the heart contractions will remain intact, while distortions resulting from hand movements should be cut off. Figure 2a shows the PPG signal containing distortions taken from a wrist. For this case, the limiter parameters should therefore cut off these high distortions, while the small sawtooth waveform should not be affected. The signal after the bandpass filter does not contain a constant component, as is presented in Figure 2b. In the case of the signal from Figure 2b, good choices for the LL and LH parameters’ values could be −150 and 150, respectively. Figure 2c shows the results of preprocessing for such a limiter, together with the band-pass biquad section. In general, these values for a given PPG system should be selected as about 100–150% of the negative and positive amplitude of the useful signal, respectively.
In the second stage of the proposed preprocessing block, a typical fourth-order band-pass filter built from two biquads was employed. The band-pass of both stages was set to 0.5–2.5 Hz.
The introduction of the limiting section in the form of a single biquad stage significantly reduces the rapid changes of the signal at the input of the band-pass filter and improves its recovery after a rapid change of the constant component in the input signal. An example of distorted signal, together with the accelerometer readouts, is presented later in this paper.

3. Peak Detection

The conditioned signal from the detector was used as an input to the block responsible for finding the peak values of the signal; the heart rate can then be easily calculated from the detected peaks. The detection of peaks is based on the AMPD algorithm [19]. The AMPD algorithm has the capability to work with noisy periodic and quasi-periodic signals. It needs the input signal to be linearly detrended, but the use of the input filter of band-pass characteristic with the limiter described in the previous section satisfies this requirement. The AMPD algorithm performs well for the filtered PPG signal, but it is computationally expensive, which can be unacceptable for wearable devices. The need to calculate a large matrix with real-valued elements, where moving windows are used, can be avoided due to the modifications of the algorithm proposed in further parts of this paper. This section starts with the detailed description of the AMPD algorithm, and then the proposed modifications are introduced. In this way, the authors wanted to make it easier for the reader to track changes that were applied to the original algorithm, without having to refer to the reference.
The main part of the AMPD algorithm is the Local Maxima Scalogram (LMS) matrix M of elements, mk,i, which is calculated for the discrete uniformly sampled signal x = [x1, x2, …, xN] in the analyzed window, where N is the constant number of samples, and scale k defines the moving window of varying length, wk, according to Equation (3):
m k , i = { 0 x i 1 > x i k 1 x i 1 > x i + k 1 r + 1 otherwise ,
where wk = 2k | k = 1,2, …, L, k is k-th scale of the signal, L = ceil(N/2) − 1, and r is a uniformly distributed random number of values [0,1]. The values of mk,i are calculated for every scale k and i = k + 2, k + 3, …, Nk + 1.
Having calculated the LMS matrix, the next step in the AMPD algorithm is to calculate the scale-dependent distribution of zeros in the LMS, by calculating vector γ = [γ1, γ2, …, γk]:
γ k = i = 1 N m k , i   for   k   { 1 ,   2 , ,   L }
and the global minimum of γ, λ = arg min(γk). The value of γ is used to obtain the matrix Mr, which is the matrix M with deleted bottom rows for k > λ. The peaks are then found for indexes, i, for which the column-wise standard deviation, σi, is equal to zero:
σ i = 1 λ 1 k = 1 λ [ ( m k , i 1 λ k = 1 λ m k , i ) 2 ] 1 2 .
Parameter λ enables the determination of how many rows of matrix M should be used for calculating the matrix Mr. The noisier the signal from the PPG sensor, the larger the value of λ.
In this paper, the authors propose modifying the AMPD algorithm toward a more efficient implementation. From practical observation, it has been inferred that the signal is too noisy, and it is of no use for peak detection and heart rate calculation, when we have the following:
λ > λmax
The value of λmax = 17 was found empirically to be a good choice. Therefore, it is practical to assume a priori, that the signals for which the condition (6) is true are of no use and the full matrix M is not used. To save the memory, only the matrix Mr, instead of M, can be calculated and stored together with the vector γ.
Processing and storing floating point values, mk,i, requires storage space and processing resources. To simplify the processing, the authors propose replacing real-valued elements of matrix Mr with the matrix Mr’ containing 1-bit binary values, mk,i, as follows:
m k , i = { 0 x i 1 > x i k 1 x i 1 > x i + k 1 1 otherwise .
This saves a lot of the device’s memory, requires only integer operations, and results in another simplification: Instead of calculating the column-wise standard deviation σi from Equation (5), calculating the column-wise summation presented in Equation (8) can be used:
s i = k = 1 λ m k , i .
As the mk,i are 1-bit binary values, si can be calculated with fast integer summation. The indices of peaks pi can be located by finding all indices, i, for which si = 0. The values of σi and si, together with the values of γk were shown in Figure 4b,c, for an exemplar input signal (Figure 4a).
The PPG signal from a wrist is weak, as can be seen in Figure 2c, and when sampled, the situation of a “flat” peak with two equal values for samples ti-1 and ti, as shown in an example in Figure 5, can sometimes occur. Such a peak would neither be detected by the AMPD algorithm nor its modified version proposed in this paper. To improve the peak detection in such a case, an additional rule was introduced: The peak is also detected at time ti.
t i = t i 1 + t i 2
for which the following simple condition (10) holds for si from Equation (8):
si-2 > 1 ∧ si-1 = 1 ∧ si=1∧ si+1 > 1.
The proposed simplified peak detection algorithm was compared to the original AMPD algorithm presented in [19]. For this purpose, the PPG signals from the PPG-DaLia database [15] were used, and the detected peaks were compared. As the original AMPD peak detection algorithm uses a random variable to fill the LMS matrix, its results are slightly different from run to run, depending on the seed value of the random number generator. The results are very similar, while the modified version can be implemented more efficiently. An example is presented in Figure 6: the arrows indicate the missing peaks, which were exclusively detected by the other algorithm. It can be seen that only a few peaks are differently detected.
The presented example uses an unfiltered PPG signal to present the robustness of the algorithm. The signal after filtration would be more regular; thus, the differences in the results between the two versions of the algorithm would be even smaller.

4. Heart Rate Calculation

In an ideal situation, where the patient is not moving, the detected peaks from the peak detection algorithm can be directly used to calculate the heart rate. However, the PPG signal can be seriously distorted by the movement of the patient. Each movement of the patient can cause a change in the saturation of the tissues with blood, and a sensor displacement on the wrist, which results in artefacts in the signal received by the optical detector. To mitigate this problem, the accelerometer is used for detecting the movement of the patient’s hand. The movement data are aligned with the detected peaks, and the periods between the peaks which are affected by movements are discarded from the calculation of the pulse period. The example of the elimination of the patient’s movement is presented in Figure 7.
To find the time periods that will be excluded from period calculations, first the acceleration values from the accelerometer are differenced to obtain the movement indicator, di, according to the following equation:
d i = 1 G ( g x , i g x , i 1 ) 2 + ( g y , i g y , i 1 ) 2 + ( g z , i g z , i 1 ) 2 ,
where i is the index of the sample; gx,i, gy,i, and gz,i are the acceleration values read from the accelerometer for axes X, Y, and Z, respectively; and G is a constant value used for normalization to obtain di ∈ [0,1] for all i and for the acceleration values in gx,i, gy,i, and gz,i in the accelerometer’s full measuring range. The movement values, di, are then compared to the constant, H, to obtain a digital binary signal Vi:
V i = { 1 d i > H 0 otherwise
Signal Vi is then filtered in time, to eliminate spikes longer than Ts. The values of H and Ts were experimentally set to H = 0.0025 and Ts = 500 ms for the accelerometer range ±2 g.
For the final result of the heart rate, the inverse of the median of the periods between the peaks not affected by the movement is calculated. For a valid result, a minimal number of peaks, P, is required, and it is calculated as shown in Equation (13):
P = f l o o r ( N B P M min 60 f s ) ,
where fs is the sampling frequency, and BPMmin is the minimal heart rate required to be measured by the system.

5. Evaluation on Dataset

For the purpose of evaluating the algorithm, the PPG-DaLia database [15], containing more than 35 h of data recorded from 15 persons, was used. The database contains signals collected from a PPG sensor, accelerometer, and ECG, where the ECG is used as ground truth. As described in [15], ground truth heart rate values were obtained from an ECG signal processed by R-peak detection algorithm [21]. Then, the detection results were manually inspected and corrected, mainly for a few cases where significant motion was observed. The ECG signal was segmented with a shifted window; ground truth heart rate was finally calculated as the mean heart rate within each window. The signals were collected during eight different types of typical daily life activities, under controlled but close-to-real-life conditions. This dataset is the longest one available to the authors. The accuracy of the algorithm was evaluated by using the mean absolute error (MAE) metric of beats per minute (bpm), calculated by using the sliding window approach of length 8 s, with a 2 s shift, according to the following equation:
M A E = 1 W j = 1 W ( B P M e s t ( j ) B P M r e f ( j ) ) ,
where W is the total number of windows, BPMest(w) is the heart rate in bpm for window j, and BPMref(j) is the reference heart rate obtained from the ECG ground truth signal for the same window j. This evaluation method is commonly used in related work [14,15,22,23].
The quality of the signal and the possibility of the measurement are checked in the TDHR algorithm. When the measurement is not possible due to not satisfying Equation (6), i.e., λ > λmax or due to movement detected by the accelerometer affecting all of the periods between the detected peaks, the measurement result is indicated as invalid. The evaluation of the algorithm is performed in two ways: (i) When the result is not available, the last valid result is used; (ii) only valid results are used to calculate the performance metric; and then the percentage of the valid samples is also given. The results of the evaluation, together with the performance of the SpaMa [14], SpaMaPlus [15], Schaeck2017 [23], CNN average, and CNN ensemble [15] algorithms are presented in Table 1. The TDHR algorithm was evaluated for several lengths of the sliding window: N = 1024 (32 s), N = 512 (16 s), N = 256 (8 s), and N = 128 (4 s).
Table 2 presents the comparison of the computational cost of several algorithms. For the TDHR algorithm, the number of operations per second was estimated from a manual analysis of the C code as the number of arithmetic operations needed to obtain a single heart-rate result. The TDHR algorithm requires only a few parameters, as opposed to the CNN algorithms, but it needs storage memory for calculating the LMS array; the size of this memory depends on the window length, N.
As can be seen from Table 1 and Table 2, the performance of the proposed TDHR algorithm is similar to SpaMa, SpaMaPlus, and Schaeck2017, while it is worse than CNN average, CNN ensemble, and CNN constrained. The computational costs of SpaMa, SpaMaPlus, and Schaeck2017 are not known, but they can be high, as those algorithms require the calculation of the power spectral density and the analysis of the PPG spectrum. The CNN-based algorithms, even CNN constrained, require larger computational costs than most versions of the TDHR algorithm. The proposed algorithm uses the mechanism of removing the time periods with body motion registered by the accelerometer, so there can be gaps between valid measurements in the case of long-lasting motions. The presented results in Table 1 and Table 2 can help to select N to achieve a compromise between the accuracy and computational cost.

6. Implementation

The proposed algorithm was implemented in low-power wearable hardware and tested. The hardware consists of the main board, power supply PCB, with external coil for wireless inductive charging, and a LiPo battery of capacity 110 mAh. A block diagram of the hardware is presented in Figure 8. As the processing unit, a PSoC6 microcontroller from Cypress was used. This is an ultra-low-power microcontroller with dual processor architecture: Arm Cortex M4 and M0+ cores. A BH1790GLC optical sensor from Rohm detects the PPG signal and also drives the four green 527 nm LEDs, equally placed in a circle of diameter of 6 mm, with the optical sensor placed in the center. The detector detects the light emitted by the LEDS and reflected from the patient’s skin. An LSM6DSL accelerometer from STMicroelectronics is mounted on the same board as the optical sensor. The optical sensor and the accelerometer are connected to the microcontroller via an I2C bus, used for configuration and data transfer.
The software for heart rate measurements runs on the Cortex M4 processing core of the PSoC6 microcontroller. The algorithm presented in this paper was written in C as two tasks running on the FreeRTOS operating system: The first task constantly reads data from the sensors, realizes filtering, and buffers data. The second task starts periodically and executes the peak detection algorithm, using data from the buffer. The calculated heart rate is sent wirelessly, using the Bluetooth Low Energy protocol. For this purpose, a built-in BLE transceiver in PSoC6 device was used. According to the BLE nomenclature, in the proposed solution, the BLE transceiver block was configured to perform a peripheral (a device constrained in resources such as energy and computing power) and server role (a device working as a data source and sending that data to the remote master device). As a data format, standard BLE Heart Rate Profile was used. The values of the measured heart rate are transmitted periodically; the user can receive the transmitted values by connecting any Bluetooth receiver compatible with Heart Rate Profile.
The PPG signal is sampled with 14-bit resolution and with frequency fs = 32 Hz, using a buffer of N = 128… 1024 samples, which enables data from 4 to 32 s to be analyzed, depending on the selected value of N, compromising accuracy versus computational cost and memory usage.
The prototype device was built and installed inside a custom-made 3D-printed case with a rubber strap, as shown in Figure 9. The size of the printed case (without the rubber strap) is 11.7 × 26 × 46 mm. The cost of the components, including the battery, was about $60, at retail prices, for 15 pieces. The proposed implementation was capable of continuous measuring of the heart rate for more than 24 h.

7. Discussion and Conclusions

Nowadays, there are many smart watches on the market that are capable of measuring the heart rate. The details regarding the optical part of the measurement, such as the number of sensors or the wavelengths of the LEDs used during measurement, are often revealed. However, the details regarding the algorithms used are not available. Most of the publications focus on measuring the accuracy of the popular devices. In [24], the authors measured the performance of Apple’s iWatch Apple Watch Sport 42 mm (first generation), during cardiopulmonary exercise test (CPET). They observed MAE from 6.34 to 7.55. It is difficult to exactly compare this to our results, as we use different and longer test conditions. In fact, the details of the PPG algorithms can only be found in the scientific publications, where the authors try to increase the accuracy of the measurement by using novel ideas and powerful techniques.
In this paper, a time-domain algorithm for the real-time detection of the heart rate was presented. The algorithm is aimed at wearable, resource-constrained devices, where battery capacity, processor speed and memory size are constrained. The algorithm processes raw data in a time-domain and requires only a few parameters. The approach is simple, the proposed algorithm has a reasonable accuracy, and it can be implemented in a typical (not DSP) microcontroller.
The algorithm consists of a two-stage input-signal-conditioning block with a limiter, the peak detection block, and period calculation block. The two-stage input conditioning block is built out of two digital bandpass filters, where the first filter has been modified to process data nonlinearly, to provide fast recovery after large signal transients. The use of band-pass filters is very simple to implement; it appeared sufficient and very effective for conditioning the signal; therefore, other methods such as wavelet-based baseline removal would not need to be considered.
The peak detection block enables operation at significantly lower processing and implementation costs, compared to the original AMPD algorithm. The period calculation block uses the median to calculate the heart rate based on the time differences between the peaks, with the use of an accelerometer to exclude the time periods affected by body movement.
In this study, most of the parameters were set experimentally, as we targeted on a simple and economical hardware implementation. It would be interesting to provide a method of automatic and dynamic adjustment of the parameters to further reduce the computational cost, basing on the input signal quality and the movement readouts. This will be a topic of our further research.
The proposed algorithm was evaluated in several variants, for different sliding window lengths, N, providing the possibility to compromise the accuracy versus lower operational costs. The authors mainly used N = 1024, which seems to be a good compromise between the calculation cost and the accuracy. The proposed algorithm was compared to the other algorithms from the literature. The achieved accuracy is comparable to the other algorithms at smaller computational costs. The proposed solution was also implemented in a low-power wrist-wearable device.

Author Contributions

Conceptualization, B.P. and M.W.; methodology, M.W.; software, M.W.; validation, M.W.; formal analysis, M.W.; investigation, B.P. and M.W.; resources, M.W.; data curation, M.W.; writing—original draft preparation, M.W. and B.P.; writing—review and editing, M.W. and B.P.; visualization, M.W.; supervision, M.W.; project administration, M.W.; funding acquisition, M.W. All authors have read and agreed to the published version of the manuscript.

Funding

This research was partially supported by the DS Programs of Faculty of Electronics, Telecommunications and Informatics, as well as National Centre for Research and Development, Poland, project “E-Pionier—using the potential of universities to improve the innovation of ICT solutions in the public sector”, No. 17/02/2018/UD.

Acknowledgments

The authors would like to thank Małgorzata Szczerska for her initiative in organizing research and raising funds for research.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Rault, T.; Bouabdallah, A.; Challal, Y.; Marin, F. A survey of energy-efficient context recognition systems using wearable sensors for healthcare applications. Pervasive Mob. Comput. 2017, 37, 23–44. [Google Scholar] [CrossRef] [Green Version]
  2. Szczepański, A.; Saeed, K. A Mobile Device System for Early Warning of ECG Anomalies. Sensors 2014, 14, 11031–11044. [Google Scholar] [CrossRef]
  3. Hertzman, A.B. Observations on the finger volume pulse recorded photoelectrically. Am. J. Physiol. 1937, 119, 334–335. [Google Scholar]
  4. Biswas, D.; Simões-Capela, N.; Van Hoof, C.; Van Helleputte, N.; Simues-Capela, N. Heart Rate Estimation from Wrist-Worn Photoplethysmography: A Review. IEEE Sensors J. 2019, 19, 6560–6570. [Google Scholar] [CrossRef]
  5. Jacobson, M. Auto-threshold peak detection in physiological signals. In Proceedings of the 2001 Conference Proceedings of the 23rd Annual International Conference of the IEEE Engineering in Medicine and Biology Society, Istanbul, Turkey, 25–28 October 2001. [Google Scholar] [CrossRef]
  6. Pan, J.; Tompkins, W.J. A Real-Time QRS Detection Algorithm. IEEE Trans. Biomed. Eng. 1985, 32, 230–236. [Google Scholar] [CrossRef] [PubMed]
  7. Du, P.; Kibbe, W.A.; Lin, S.M. Improved peak detection in mass spectrum by incorporating continuous wavelet transform-based pattern matching. Bioinform 2006, 22, 2059–2065. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  8. Benitez, D.; Gaydecki, P.; Zaidi, A.; Fitzpatrick, A.P. The use of the Hilbert transform in ECG signal analysis. Comput. Boil. Med. 2001, 31, 399–406. [Google Scholar] [CrossRef]
  9. Aboy, M.; McNames, J.; Thong, T.; Tsunami, D.; Ellenby, M.; Goldstein, B. An Automatic Beat Detection Algorithm for Pressure Signals. IEEE Trans. Biomed. Eng. 2005, 52, 1662–1670. [Google Scholar] [CrossRef] [PubMed]
  10. Tzallas, A.T.; Oikonomou, V.P.; I Fotiadis, D. Epileptic Spike Detection Using a Kalman Filter Based Approach. In Proceedings of the 2006 International Conference of the IEEE Engineering in Medicine and Biology Society, New York, NY, USA, 30 August–3 September 2006; Volume 1, pp. 501–504. [Google Scholar] [CrossRef]
  11. Coast, D.A.; Stern, R.; Cano, G.G.; A Briller, S. An approach to cardiac arrhythmia analysis using hidden Markov models. IEEE Trans. Biomed. Eng. 1990, 37, 826–836. [Google Scholar] [CrossRef] [PubMed]
  12. Xue, Q.; Hu, Y.H.; Tompkins, W.J. Neural-network-based adaptive matched filtering for QRS detection. IEEE Trans. Biomed. Eng. 1992, 39, 317–329. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  13. Mashhadi, M.B.; Asadi, E.; Eskandari, M.; Kiani, S.; Marvasti, F.; Boloursaz, M. Heart Rate Tracking using Wrist-Type Photoplethysmographic (PPG) Signals during Physical Exercise with Simultaneous Accelerometry. IEEE Signal Process. Lett. 2015, 23, 227–231. [Google Scholar] [CrossRef] [Green Version]
  14. Salehizadeh, S.M.A.; Dao, D.; Bolkhovsky, J.; Cho, C.; Mendelson, Y.; Chon, K.H. A Novel Time-Varying Spectral Filtering Algorithm for Reconstruction of Motion Artifact Corrupted Heart Rate Signals During Intense Physical Activities Using a Wearable Photoplethysmogram Sensor. Sensors 2016, 16, 10. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  15. Reiss, A.; Indlekofer, I.; Schmidt, P.; Van Laerhoven, K. Deep PPG: Large-Scale Heart Rate Estimation with Convolutional Neural Networks. Sensors 2019, 19, 3079. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  16. Sumukha, B.N.; Kumar, R.C.; Bharadwaj, S.S.; George, K. Online peak detection in photoplethysmogram signals using sequential learning algorithm. In Proceedings of the 2017 International Joint Conference on Neural Networks (IJCNN), Anchorage, AK, USA, 14–19 May 2017; pp. 1313–1320. [Google Scholar] [CrossRef]
  17. Vijaya, G.; Kumar, V.; Verma, H.K. ANN-based QRS-complex analysis of ECG. J. Med Eng. Technol. 1998, 22, 160–167. [Google Scholar] [CrossRef] [PubMed]
  18. Ghamari, M.; Castaneda, D.; Esparza, A.; Soltanpur, C.; Nazeran, H. A review on wearable photoplethysmography sensors and their potential future applications in health care. Int. J. Biosens. Bioelectron. 2018, 4, 195–202. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  19. Scholkmann, F.; Boss, J.; Wolf, M. An Efficient Algorithm for Automatic Peak Detection in Noisy Periodic and Quasi-Periodic Signals. Algorithms 2012, 5, 588–603. [Google Scholar] [CrossRef] [Green Version]
  20. Nasal ALAR SpO2 Pulse Oximetry Sensor. Available online: https://www.pentlandmedical.co.uk/critical-care/nasal-alar-spo2-pulse-oximetry-sensor/ (accessed on 4 February 2020).
  21. Hamilton, P. Open source ECG analysis. Computers in Cardiology 2003, 101–104. [Google Scholar] [CrossRef]
  22. Zhang, Z. Photoplethysmography-Based Heart Rate Monitoring in Physical Activities via Joint Sparse Spectrum Reconstruction. IEEE Trans. Biomed. Eng. 2015, 62, 1902–1910. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  23. Schack, T.; Muma, M.; Zoubir, A.M. Computationally efficient heart rate estimation during physical exercise using photoplethysmographic signals. In Proceedings of the 2017 25th European Signal Processing Conference (EUSIPCO), Kos, Greece, 28 August–2 September 2017; pp. 2478–2481. [Google Scholar]
  24. Falter, M.; Budts, W.; Goetschalckx, K.; Cornelissen, V.; Buys, R.; Shcherbina, A.; Goessler, K.; Boudreaux, B.; Goris, J. Accuracy of Apple Watch Measurements for Heart Rate and Energy Expenditure in Patients With Cardiovascular Disease: Cross-Sectional Study. JMIR mHealth uHealth 2019, 7, e11889. [Google Scholar] [CrossRef] [PubMed]
Figure 1. Photodetector(PD) waveforms: raw signals obtained from photoplethysmography (PPG) sensor in reflective configuration placed on the index tip (dotted line) and on the wrist (continuous line).
Figure 1. Photodetector(PD) waveforms: raw signals obtained from photoplethysmography (PPG) sensor in reflective configuration placed on the index tip (dotted line) and on the wrist (continuous line).
Sensors 20 01783 g001
Figure 2. Waveforms of the signals from the PPG detector: (a) the raw signal; (b) the raw signal from (a) filtered only by a band-pass filter; and (c) the raw signal from (a) filtered first with the limiting section described in the paper and then by the same biquadratic filter as used in (b).
Figure 2. Waveforms of the signals from the PPG detector: (a) the raw signal; (b) the raw signal from (a) filtered only by a band-pass filter; and (c) the raw signal from (a) filtered first with the limiting section described in the paper and then by the same biquadratic filter as used in (b).
Sensors 20 01783 g002aSensors 20 01783 g002b
Figure 3. Block diagram of PPG signal acquisition, preprocessing, and final HR estimation by TDHR algorithm.
Figure 3. Block diagram of PPG signal acquisition, preprocessing, and final HR estimation by TDHR algorithm.
Sensors 20 01783 g003
Figure 4. Sample of (a) filtered input PPG signal, (b) calculated and normalized values of σi from Equation (5) and si from Equation (8), and (c) calculated values of γk from Equation (4).
Figure 4. Sample of (a) filtered input PPG signal, (b) calculated and normalized values of σi from Equation (5) and si from Equation (8), and (c) calculated values of γk from Equation (4).
Sensors 20 01783 g004
Figure 5. An example of special case of a “flat” peak consisting of two equal values at ti-1 and ti.
Figure 5. An example of special case of a “flat” peak consisting of two equal values at ti-1 and ti.
Sensors 20 01783 g005
Figure 6. Comparison of the results of the peak detection algorithm proposed in this paper (a) with the AMPD peak detection algorithm (b) from [19] on the same waveform. The arrows indicate the peaks not detected by one of the algorithms but detected by the other one. All other peaks were detected at the same positions, by both methods.
Figure 6. Comparison of the results of the peak detection algorithm proposed in this paper (a) with the AMPD peak detection algorithm (b) from [19] on the same waveform. The arrows indicate the peaks not detected by one of the algorithms but detected by the other one. All other peaks were detected at the same positions, by both methods.
Sensors 20 01783 g006
Figure 7. Elimination of the patient’s movement from pulse-signal detection: (a) the PPG signal distorted with the patient’s movement; the areas in gray are excluded from pulse period calculation due to the detected movement; (b) acceleration values read from the accelerometer placed together with the heart rate detector on the wrist, (c) and the values of movement di from Equation (11); (d) signal di after threshold, as in Equation (12), and with the eliminated short-term peaks (the eliminated peak is shown with a dotted line).
Figure 7. Elimination of the patient’s movement from pulse-signal detection: (a) the PPG signal distorted with the patient’s movement; the areas in gray are excluded from pulse period calculation due to the detected movement; (b) acceleration values read from the accelerometer placed together with the heart rate detector on the wrist, (c) and the values of movement di from Equation (11); (d) signal di after threshold, as in Equation (12), and with the eliminated short-term peaks (the eliminated peak is shown with a dotted line).
Sensors 20 01783 g007
Figure 8. Block diagram of the device implementing the algorithm for measuring the heart rate.
Figure 8. Block diagram of the device implementing the algorithm for measuring the heart rate.
Sensors 20 01783 g008
Figure 9. The prototype implementation of a wrist sensor for heart rate measurements with the implemented TDHR algorithm: (a) the picture of the internal modules and parts; (b) the picture of the final device with strap. The main board has the dimensions of 20 × 30 mm.
Figure 9. The prototype implementation of a wrist sensor for heart rate measurements with the implemented TDHR algorithm: (a) the picture of the internal modules and parts; (b) the picture of the final device with strap. The main board has the dimensions of 20 × 30 mm.
Sensors 20 01783 g009
Table 1. Comparison of the accuracy of the presented algorithm to the algorithms from the literature on the large PPG-DaLia dataset as MAE (bpm). The proposed TDHR algorithm was tested in several versions, with various numbers of analyzed samples N (for all other algorithms N = 256), sampled with a frequency of 32 Hz, with a BPMmin = 40 bpm. The measurements were done every 2 s. For the TDHR algorithm, the accuracy was also calculated for valid measurements, and then the percentage of valid measurements was also given.
Table 1. Comparison of the accuracy of the presented algorithm to the algorithms from the literature on the large PPG-DaLia dataset as MAE (bpm). The proposed TDHR algorithm was tested in several versions, with various numbers of analyzed samples N (for all other algorithms N = 256), sampled with a frequency of 32 Hz, with a BPMmin = 40 bpm. The measurements were done every 2 s. For the TDHR algorithm, the accuracy was also calculated for valid measurements, and then the percentage of valid measurements was also given.
S1S2S3S4S5S6S7S8S9S10S11S12S13S14S15All
SpaMa [14]11.8614.759.5317.239.2816.7815.8815.217.199.0821.6312.639.510.7312.2315.56 ± 7.5
SpaMaPlus [15]8.869.676.414.124.0611.346.3111.2516.046.1715.1512.038.57.768.2911.06 ± 4.8
Schaeck 2017 [23]33.0527.8118.4928.8212.648.7220.6521.7522.2512.621.0522.7427.7112.0516.420.45 ± 7.1
CNN average [15]8.457.925.967.8618.9713.555.1611.4910.656.079.879.955.255.855.258.82 ± 3.8
CNN ensemble [15]7.736.744.035.918.5112.883.9110.878.794.039.229.354.294.374.177.65 ± 4.2
TDHR
N = 1024
8.107.9810.5111.8220.6012.117.6211.7114.794.9220.058.767.889.158.4510.96 ± 4.49
TDHR
N = 1024
only valid measurements (%)
5.32
74%
4.99
86%
6.47
71%
6.57
72%
16.26
43%
8.29
55%
5.40
80%
10.28
80%
6.93
64%
3.01
81%
9.17
66%
5.22
61%
4.70
83%
4.84
79%
4.12
75%
6.77 ± 3.26
72%
TDHR
N = 512
8.617.4911.5511.9320.7914.247.9611.9114.955.8319.608.868.229.108.6611.31 ± 4.41
TDHR
N = 512
only valid measurements (%)
6.00
74%
5.60
86%
7.72
73%
7.55
71%
16.05
58%
9.62
65%
6.01
81%
10.50
79%
8.79
65%
3.71
80%
12.15
71%
5.65
62%
5.27
84%
5.34
78%
4.82
75%
7.65 ± 3.30
74%
TDHR
N = 256
11.0810.8411.6314.0621.6715.638.8613.3015.127.2721.0810.499.5410.249.0412.66 ± 4.25
TDHR
N = 256
only valid measurements (%)
8.40
74%
8.21
86%
9.40
75%
10.10
73%
17.02
69%
12.85
72%
7.57
83%
12.04
78%
11.75
68%
5.14
79%
13.86
74%
7.73
62%
7.33
84%
7.45
81%
6.63
74%
9.70 ± 3.21
76%
TDHR
N = 128
13.1313.5713.1115.5925.5518.0510.2516.1919.2310.3219.3912.6411.6312.3312.8014.92 ± 4.16
TDHR
N = 128
only valid measurements (%)
12.71
80%
12.25
89%
11.93
82%
14.05
80%
20.84
82%
16.44
79%
9.53
86%
15.65
83%
17.53
78%
8.75
85%
17.13
82%
12.02
72%
10.29
86%
10.66
85%
10.73
80%
13.37 ± 3.46
82%
Table 2. Comparison of the performance versus computational cost of the TDHR algorithm and the algorithms from the literature. The computational cost of the TDHR algorithm was calculated as the number of arithmetical operations and the number of memory bytes needed for algorithm realization in a microcontroller.
Table 2. Comparison of the performance versus computational cost of the TDHR algorithm and the algorithms from the literature. The computational cost of the TDHR algorithm was calculated as the number of arithmetical operations and the number of memory bytes needed for algorithm realization in a microcontroller.
AlgorithmPerformance
Mean MAE ± STD
(% = Only Valid Measurements)
Computational Cost
Number of Parameters/Memory BytesOperations Per Second
CNN average8.82 ± 3.88.5 M34.5 M
CNN ensemble7.65 ± 4.260 M240 M
CNN constrained9.99 ± 5.926 K190 K
TDHR N = 102410.96 ± 4.4966 k2.4 M
TDHR N = 10246.77 ± 3.26 (72%)
TDHR N = 51211.31 ± 4.4116 k598 k
TDHR N = 5127.65 ± 3.30 (74%)
TDHR N = 25612.66 ± 4.254 k152 k
TDHR N = 2569.70 ± 3.21 (76%)
TDHR N = 12814.92 ± 4.161 k40 k
TDHR N = 12813.37 ± 3.46(82%)

Share and Cite

MDPI and ACS Style

Wójcikowski, M.; Pankiewicz, B. Photoplethysmographic Time-Domain Heart Rate Measurement Algorithm for Resource-Constrained Wearable Devices and its Implementation. Sensors 2020, 20, 1783. https://doi.org/10.3390/s20061783

AMA Style

Wójcikowski M, Pankiewicz B. Photoplethysmographic Time-Domain Heart Rate Measurement Algorithm for Resource-Constrained Wearable Devices and its Implementation. Sensors. 2020; 20(6):1783. https://doi.org/10.3390/s20061783

Chicago/Turabian Style

Wójcikowski, Marek, and Bogdan Pankiewicz. 2020. "Photoplethysmographic Time-Domain Heart Rate Measurement Algorithm for Resource-Constrained Wearable Devices and its Implementation" Sensors 20, no. 6: 1783. https://doi.org/10.3390/s20061783

APA Style

Wójcikowski, M., & Pankiewicz, B. (2020). Photoplethysmographic Time-Domain Heart Rate Measurement Algorithm for Resource-Constrained Wearable Devices and its Implementation. Sensors, 20(6), 1783. https://doi.org/10.3390/s20061783

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop