=Paper= {{Paper |id=Vol-2843/paper34 |storemode=property |title=Heart rate intellectual analysis by structural decomposition methods of photoplethysmography signals (paper) |pdfUrl=https://ceur-ws.org/Vol-2843/paper034.pdf |volume=Vol-2843 |authors=Andrey Borzov,Alexey Kasikin,Leonid Labunets,Mariya Ryakhina }} ==Heart rate intellectual analysis by structural decomposition methods of photoplethysmography signals (paper)== https://ceur-ws.org/Vol-2843/paper034.pdf
     Heart rate intellectual analysis by structural
decomposition methods of photoplethysmography signals*

       Andrey Borzov1, Alexey Kasikin1, Leonid Labunets1,2 and Mariya Ryakhina1
      1 Bauman Moscow State Technical University, 5, 2-nd Baumanskaya, Moscow, 105005,

                                            Russian Federation
    2 Autonomous non-profit organization of higher education Russian new university, 22, Radio

                              street, Moscow, 105005, Russian Federation
                                       labunets@bmstu.ru



          Abstract. In this paper, various algorithms for the analysis of
          photoplethysmography signals are considered. The sequence of performing the
          intelligent analysis of the measured signals is presented. The study of heart rate
          variability was carried out, based on the assessment of the instantaneous
          frequency and the instantaneous period of the pulse wave. The estimates of heart
          rate and indicators of variation pulsometry were obtained.

          Keywords: Photoplethysmography,               Multiresolutional      analysis,    Signal
          processing, Heart rate variability.


1         Introduction

An important subject area of data mining (ADM) is the development of medical
decision support systems (DSS) for recognizing a patient's diagnosis based on a
presented set of symptoms [1-2]. One of these promising areas is
photoplethysmography is a remote PhotoPlethysmography - rPPG. This method allows
you to quickly identify abnormalities in the physical state of a person and is also an
effective methodology for remote express diagnostics in the telemedicine format. Of
practical interest is the analysis of biological wave processes in the human body - pulse,
respiratory, myogenic, etc.
   In solving the problem of identifying biological wave processes, an important role
is played by the formation of the space of informative features based on the measured
initial data. The basis for solving this problem is the use of artificial intelligence
methodology in the form of a natural symbiosis of methods for structural and spectral
analysis of the dynamics of non-stationary time series (TS).
   Remote PhotoPlethysmography monitors changes in the intensity of light reflected
from a subject's skin. An RGB video camera is used as a sensor, and the area of interest
is most often the face, less often the palm or wrist.

*
    Copyright © 2021 for this paper by its authors. Use permitted under Creative Commons License Attribution
4.0 International (CC BY 4.0).
   The aim of this work is the intellectual analysis of heart rate variability (HRV) by
various algorithms for processing photoplethysmography signals. It is important to note
that HRV analysis refers to the pre-medical stage, i.e. does not make a diagnosis, but
forms the basis for its statement.


2      Materials and methods

2.1    Dataset description

As the initial data for the study, we took a database recorded by volunteers of the
Eindhoven University of Technology (Netherlands), and designed to study the methods
of photoplethysmography [3]. The database (DB) contains a set of video recordings, as
well as reference HRV indicators, taken simultaneously with a photoplethysmogram.
There are videotapes of the facial area of three subjects with different skin colors. For
each participant in the experiments, a recording was carried out after physical activity,
a set of recordings in a calm state, but under different lighting intensities, as well as
recording with head turns (movement). Thus, based on this database, it is possible to
study the operation of algorithms under various conditions, that is, their adaptability.
   In this work, we analyze relatively short recordings of rPPG data - recording duration
no more than 15 minutes at a frequency of 30 frames / sec. Such a measurement mode
is typical for mass prophylactic examinations or during preliminary outpatient and
clinical examinations, which allows registering the general state of health of the subject
and finding out some significant indicators.
   An area of the face is selected from the current frame - the so-called. the area of
interest of certain areas of the skin of the face. Based on the results of measuring the
degree of light scattering, the initial signals are formed in the form of TS. Each pixel in
the image represents the red (R), green (G), and blue (B) color components. Thus, the
face area can be represented as an array, the rows of which are frames, and the columns
are color components. To obtain a pulse signal, which, in turn, does not depend on the
spectral characteristics of a stationary illumination source and its brightness level, it is
necessary to normalize each color channel by dividing the TS realizations by their
average values (trends).
   The paper considers promising methods for analyzing the dynamics of indicators of
non-stationary TS - stabilization of the pulse and respiration after physical exercises.

2.2    Methods

To determine a person's health status, an rPPG signal is required, which can be
compared with a cardiogram or rhythmogram. This is due to the need to isolate the
sequence of RR intervals for HRV analysis. ADM models and algorithms are a reliable
methodological basis for the formation of an informative component, i.e. an aggregated
signal based on a set of RGB signals. The general algorithm for analyzing a
photoplethysmogram consists of the following steps:
─ For each frame, a face area is selected according to the Viola-Jones algorithm [4]
  and skin pixels according to the skin tone. The intensity value in each color channel
  is calculated as the average of the pixel intensity of the region of interest.
─ Structural decomposition of RGB non-stationary TS is performed using
  multiresolution analysis (MRA) in the discrete wavelet transform basis [5]. Removes
  color channel trend estimates.
─ An informative pulse wave frequency range from 0.667 to 4 Hz is identified
  according to the Analysis Mode Decomposition (AMD) method [6].
─ An informative pulse wave signal is generated in the form of a combination of color
  signals according to the algorithms.
─ Estimate the instantaneous phase and frequency of the pulse wave based on the
  Hilbert-Huang transform.
─ Assess the indicators of heart rate variability.

   Allocation of the frequency range of the pulse wave using the AMD method. AMD
extracts informative harmonic components from a wideband signal as follows:
𝑥𝑠 (𝑡) = 𝑠𝑖𝑛(2𝜋𝑓𝑏 𝑡) ∙ 𝐻[𝑐𝑜𝑠(2𝜋𝑓𝑏 𝑡)𝑥0 (𝑡)] − 𝑐𝑜𝑠(2𝜋𝑓𝑏 𝑡) ∙ 𝐻[𝑠𝑖𝑛(2𝜋𝑓𝑏 𝑡)𝑥0 (𝑡)] (1)

                              𝑥𝑓 (𝑡) = 𝑥0 (𝑡) − 𝑥𝑠 (𝑡)                               (2)
   Where x0 (t) is the original signal; fb - boundary frequency; xf (t) - fast changing
signal; xs (t) is a slowly changing signal.
   Extraction of photoplethysmogram from TS color channels. CHROM algorithm -
provides invariance of informative signals to the spectral composition of lighting as a
result of aggregation of components based on color difference, i.e. chrominance,
signals. The idea of the method consists in orthogonal projection of TS samples of color
channels onto the plane specified by the orts X = (3, -2, 0)/√13 and Y = (-1.5, -1,
1.5)/√5.5 in RGB-space, informative orthonormal basis (X, Y, Z) [7], where Z = (1, 1,
1)/√3 is the unit of skin tone [8]:
                                              𝜎(𝑋𝑓 )
                            𝐶𝐻𝑅𝑂𝑀 = 𝑋𝑓 +               𝑌𝑓 .                          (3)
                                              𝜎(𝑌𝑓 )
   Here Xf = (3R - 2G)/√13 and Yf = (-1.5R - G+ 1.5B)/√5.5 are two orthogonal chroma
signals at the output of the FIR filter with a bandwidth equal to the frequency range of
the pulse wave; σ is the operator of a robust estimate of the standart deviation (MSD).
   In the POS algorithm, the unit vectors of the informative plane orthogonal to the skin
tone unit Z in RGB space are chosen as follows: X = (0, 1, -1)/√2 and Y = (-2, 1, 1) /
√6. Accordingly, the projections of the readings of the TS-centered color channels X =
(G - B) /√2 and Y = (-2R + G + B)/√6 onto this plane form an estimate of the pulse wave
[8]:
                                            𝜎(𝑋)
                               𝑃𝑂𝑆 = 𝑋 +         𝑌,                                  (4)
                                            𝜎(𝑌)
  Where σ(X) and σ(Y) are robust standard deviations of nonstationary TS X and Y.
Estimates of the instantaneous phase and frequency of the pulse signal based on
the Hilbert - Huang transform. The Hilbert transform provides adaptive estimates of
the instantaneous envelope, phase, and frequency of a signal. The total phase of the
analytical signal is:
                                                         𝑠̂ (𝑡)
                       𝜑𝑠 (𝑡) = arg(𝑧(𝑡)) = 𝑎𝑟𝑐𝑡𝑔                                   (5)
                                                         𝑠(𝑡)
    The instantaneous frequency is defined as the time derivative of the total phase:
                          d        𝑠̂ (𝑡) 𝑠̂ ′(𝑡)𝑠(𝑡) − 𝑠̂ (𝑡)𝑠′(𝑡)
               𝜔𝑠 (𝑡) =      𝑎𝑟𝑐𝑡𝑔       =                          ,               (6)
                          dt       𝑠(𝑡)         𝑠̂ 2 (𝑡) + 𝑠 2 (𝑡)
    Where 𝑠(𝑡) is a pulse signal; 𝑠̂ (𝑡) - the signal obtained as a result of the Hilbert
transformation of the pulse signal.
    The energy of the Hilbert - Huang spectrum is the square of the pulse signal
envelope 𝐸(𝑡) = 𝑠 2 (𝑡) + 𝑠̂ 2 (𝑡).
Assessment of indicators of heart rate variability. The RR time sequence of the pulse
wave intervals is very informative. Based on this characteristic of the pulse wave,
various HRV indicators are calculated, one of which is heart rate. It is defined as the
number of samples for the total recording time.
                                                     n
                              HR = 60 1000 n
                                                                                    (7)
                                              NN ( мс)
                                              i =1
                                                     i



    Where NN is the normal RR interval.


3       Results

3.1     Allocation of the informative component of color channels

The TS extracted from the RGB video are subjected to multiresolution analysis (MRA).
In our computational experiments, the 40th order Debeschy wavelets were used as a
basis for the expansion. Ten levels of decomposition made it possible to form an expert
model of the structural components of non-stationary TS color channels of rPPG,
namely:

─ Trends (Figure 1) as the sum of the eighth, ninth, tenth detailing and tenth
  approximating components of the MRA;
─ Respiratory, myogenic and other wave processes with frequencies not exceeding
  0.667 Hz in the form of the sum of the slow part of the fifth, sixth and seventh
  detailing components of the MRA;
─ Centered pulse waves of each color channel (Figure 2) in the frequency range from
  0.667…4 Hz as the sum of the slow part [6] of the second, third, fourth and fast parts
  of the fifth detailing components of the MRA normalized to the trend value;
─ Noise components with frequencies of at least 4 Hz as the sum of the first and fast
  parts of the second detailing components.
               а)                              b)                              c)
          Fig. 1. Time series of RGB channels and their trends: R - a), G - b), B - c).




               а)                              b)                              c)
   Fig. 2. Fragments (10 c) of centered pulse waves of color channels: R - a), G - b), B - c).


   Figure 3 demonstrates the topology of a 3D scatter diagram of the TS samples of
centered pulse waves (Figure 2) in RGB space in projections onto the coordinate planes
R-G, R-B and G-B, respectively.
   It is important to note that the absence of anomalous values in the dynamics of the
TS of the centered pulse waves of the color channels makes it possible to obtain an
effective assessment of the aggregated pulse wave and its correlation - spectral
characteristics according to the corresponding algorithm. In particular, the
implementation of the non-stationary TS of the pulse wave pw and its averaged
periodogram estimates Psd of the Schuster and Thomson power spectral density (PSD)
obtained using the POS algorithm are shown in Figure 4.
   The presence of dominant harmonics at frequencies of 2.50 is characteristic; 2.08;
1.71; Hz, (Figure 4, b), demonstrating the dynamics of stabilization of instantaneous
frequencies and, accordingly, heart rate after physical exercises of the subject.
   For a more detailed analysis of this dynamics, the time interval of measurements
with a duration of 300 s, containing 9000 TS counts, was divided into 200 overlapping
segments of 1024 counts each. Figure 5 illustrates the dependence of the sample
estimates of autocorrelation functions (ACF) and Thomson periodograms of the pulse
wave pw on the segment number ns.
Fig. 3. Projections of the 3D scattering diagram of the TS counts of the centered pulse
                  waves onto the coordinate planes of the RGB space.




                   а)                                              b)
Fig. 4. Fragment (10 s) of the pulse wave - a) and its PSD - b), averaged over the entire
                                 implementation of TS.




                   а)                                              b)
    Fig. 5. ACF of overlapping segments of the pulse wave - a) and their Thomson
                                  periodogram - b).
   It is clearly seen that as ns increases, the intervals between the zeros of the ACF
increase, and the largest values of the PSD of the segments shift to the region of lower
frequencies.

3.2      Assessment of indicators of heart rate variability

A simplified version of the heart rate assessment for a subject in a static state does not
require the restoration of the pulse waveform by the rPPG method. The frequency
range, including the fundamental heart rate harmonic, is 0.8 ... 1.6 Hz [11]. Also,
considering this simplification, the values of the instantaneous periods are taken as R-
R intervals, which are the basis for the HRV analysis. Average heart rate values for
each algorithm are shown in Table 1. Reference values are presented in Table 2.

                           Table 1. Heart rate values for each algorithm.

        Algorithm                              POS                                 CHROM
         HR, bpm                               133                                  131

                                  Table 2. HR ground truth values.

                                             Standart
                           Mean                               Minimal HR           Maximum HR
                                             deviation
      Value, bpm            107                7.69                  90               134

   The shift of the heart rate to the region of high frequencies is associated with the
peculiarity of the TS example chosen for the study.
   The data in Table 3 demonstrate the dynamics of heart rate for 5 non-overlapping
video segments.

                    Table 3. HRV indices for video with division into intervals.

Algorithm                           No.                  HR, bpm                   R-R, s
                                     1                     141                      0.42
                                     2                     142                      0.42
CHROM                                3                     147                      0.40
                                     4                     138                      0.43
                                     5                     150                      0.40
                                     1                     133                      0.45
                                     2                     131                      0.45
POS                                  3                     138                      0.43
                                     4                     138                      0.43
                                     5                     146                      0.41
   There is a significant correlation between the heart rate values, the instantaneous
frequency values and their correspondence to the visual assessment of the subject's
condition on video.
   One of the significant disadvantages of estimating (1) using the Hilbert transform
hilbert(pw) is, as is known, the presence of negative values of instantaneous
frequencies. The physically non-interpretable result is due to errors in the calculation
of the accumulated instantaneous phase (unwrap phase). In our opinion, this drawback
can be easily eliminated by dividing the accumulated instantaneous phase into two
components. The first is a continuous, differentiable, and monotonically non-
decreasing function of time. The second component contains discontinuities of the first
kind, i.e. is a non-decreasing "step function", thus taking into account the errors in the
estimates of the accumulated instantaneous phase.
   Within the framework of this approach, the continuous component of the
accumulated instantaneous phase is described using shape-preserving piecewise cubic
interpolation by Hermite polynomials [12]. The numerical estimate of the instantaneous
frequency is the first derivative of this model (Figure 6, a).




                     а)                                             b)
        Fig. 6. Dynamics of instantaneous frequency - a) and energy - b) pulse wave.


   The solid line in Figure 6.a represents the nuclear robust locally weighted quadratic
regression of Cleveland [13] and demonstrates the stabilization of the average
instantaneous frequency of the pulse wave in time after the end of physical exercise of
the subject.


4      Discussion

Correct determination of HRV and heart rate indices requires a chronological sequence
of RR intervals, similar to ECG data [14, 15]. Otherwise, there is a significant deviation
of the indicators from the norm. To solve this problem, it is necessary to restore the
shape of the pulse wave using the rPPG method.
5      Conclusion

The paper describes algorithms for extracting an informative component - a pulse signal
from noisy data. For the pulse wave, its main characteristics, required for HRV analysis,
were identified, and an example of analysis for one video sequence was given. In the
example, both a complete record and its splitting were considered to identify the
dynamics of indicators.
   Based on the results obtained, it can be concluded that a structural decomposition is
necessary in the analysis of rPPG signals, since interference at extraneous frequencies
has a high-power relative to the required spectral range. This will allow restoring the
shape of the pulse wave, obtaining a cardiogram, excluding contact with the patient.
   It should be noted that the problem of reducing the degree of mixing of modes
because of signal decomposition is essential for reconstructing the pulse waveform.


References
 1. Atkov O.Yu., Kudryashov Yu.Yu., Prohorov A.A., Kasimov O.V.: Medical decision support
    systems. Medical decision support systems, 6, 65-75 (2012).
 2. Zaripova G.R., Bogdanova Yu.A., Galimov O.V., Kataev V.A., Bikkinina G.M.: Modern
    models of medical descision support systems in surgery. Bashkortostan medical bulletin, 11,
    6(66), 96-101.
 3. Addressing Reproducability Issues in Remote Photoplethysmography (rPPG) Research: An
    Investigation of Current Challenges and Release of a Public Algorithm Benchmarking
    Dataset. https://osf.io/xjf7u/, last accesed 2020/09/20.
 4. Viola P., Jones M.J.: Rapid object detection using a boosted cascade of simple features.
    IEEE Conf. on Computer Vision and Pattern Recognition. Kauai, Hawaii, USA, 1, 511–518
    (2001).
 5. Labunets L.V. Automatic Intellectual Signal Processing in Subsurface Radar Systems.
    Journal of Communications Technology and Electronics, 60, 4, 362 – 374 (2015).
 6. Zhongzhe Ch., Baqiao L., Xiaogang Y., Hongquan Y.: An Improved Signal Processing
    Approach Based on Analysis Mode Decomposition and Empirical Mode Decomposition.
    Energies, 12 (2019).
 7. deHaan, G., Jeanne, V.: Robust pulse-rate from chrominance-based rPPG. IEEE
    Transactions on Biomedical Engineering, 60, 10, 2878-2886 (2013).
 8. Wang, W., den Brinker, A.C., Stuijk, S., de Haan, G.: Algorithmic principles of remote-
    PPG. IEEE Transactions on Biomedical Engineering, 64, 7, 1479–1491 (2017).
 9. Baevskiy, R.M., Ivanov, G.G.: Heart rate variability: theoretical aspects and clinical
    potential. Ultrasound and functional diagnostics, 3, 106-127 (2001).
10. Hertzman, B.: Photoelectric plethysmography of the fingers and toes in man. Exp.Biol.
    Med., 37(3), 529-534 (1937).
11. Krupatkin, A.I., Sidorov, V.V.: Laser Doppler Flowmetry of Blood Microcirculation: A
    Guide for Physicians. Medicine, Moscow (2005).
12. Fritsch F.N., Carlson R.E.: Monotone Piecewise Cubic Interpolation. SIAM Journal on
    Numerical Analysis, 17, 238-246 (1980).
13. Cleveland W.S.: Robust locally weighted regression and smoothing scatterplots. Journal of
    the American Statistical Association, 74, 368, 829-836 (1979).
14. Sapunova, D. A.: Vascular rigidity and cardiovascular diseases. Medical Bulletin, 4-5, 581-
    582 (2012).
15. Mathematical methods of heart rate variability analysis. In: Larin, V.V., Baevskiy R. M.
    (eds.) 1st All-Soviet Union Symposium, Science, Moscow (1968).