1. Introduction

2. Related Work

3. Method

3.1. Preprocessing

Preprocessing encompasses the set of techniques to obtain the plethysmographic signal from the frames that make up the video sequence. This signal serves as input for models responsible for extracting heart rate and blood pressure.

The first step involves extracting the frames from the videos in the dataset. Next, face landmarks are obtained for each frame. We used the Face Landmark Detection model from MediaPipe[8]. This model uses Convolutional Neural Networks (CNNs) to locate key facial landmarks, allowing the extraction of information from strategic regions. In this case, the right cheek is selected as the region of interest due to its lower detection error.

Once the signal is extracted, it is decomposed into RGB components, as the process for obtaining heart rate uses only the green channel, while blood pressure estimation uses all components separately.

The next step is to filter the signal, since the relevant information for obtaining heart rate and blood pressure is within a specific range of the frequency spectrum. A healthy human pulse ranges between 40 and 240 beats per minute (bpm) [9], which, converted to seconds, corresponds to 0.67 and 4 beats per second (bps) or Hz. To restrict the information to this range, a Chebyshev Type II filter is used, as it is a low-pass filter characterized by an attenuation band with a ripple zone and a flat passband with a rapid transition. After filtering, a clean and precise signal is obtained, which is essential for the subsequent analysis and extraction of physiological parameters by the estimation models.

3.1.1. Dataset

The dataset consists of two data sets: one composed of employees from our company and the other from the UBFC-rPPG dataset [10]. In both datasets, most participants reported being normotensive and having a heart rate within normal parameters, that is, between 60 and 100 beats per minute at rest. At the start of the study, participants completed the profile on the specially created iOS app for data collection. The following information was recorded for each participant:

• Room temperature (measured with a thermometer in ºC);
• Age (determined by date of birth);
• Gender (female/male);
• Weight (in kg);
• Height (in cm);
• Disposition (whether the participant is active or relaxed during the measurement);
• Body Mass Index (auto-completed with weight and height).

After completing these steps, participants were instructed to take a seat with their feet on the floor and position themselves in front of the tablet, ensuring that their head was centered. The next step involved placing the index finger of the right hand in a pulse oximeter to measure heart rate. A cuff was also applied for blood pressure measurement. Once both devices were properly placed, the video recording began simultaneously with the heart rate measurement. The duration of the video is approximately one minute.

Figure 1. Recording Process Image.

Figure 1. Recording Process Image.

3.1.2. Motion Artifacts

During the video recording, the participant is given several instructions, one of which is to remain still in front of the camera. However, small movements or changes in illumination are inevitable. These movements and changes in illumination are known as ’motion artifacts’ and can interfere with the signal of interest in the frequency spectrum. This issue is significant because heart rate is estimated specifically from the peak of highest amplitude. In the case of blood pressure, these movements can lead to errors in the selection of components during the signal segmentation process.

The solution developed to mitigate this disadvantage is to track the face to obtain a movement signal along the X and Y axes [11]. This strategy involves recording the coordinates of two facial points for each frame and transforming the movement into a signal similar to the plethysmographic signal.

To extract this signal, the movement of each eye is calculated across each video frame, that is, the landmark of each eye in every frame along the X and Y axes - Figure 2. Once these signals are obtained, one per eye, normalization is performed. The normalization of the signal, prior to the FFT, ensures that the amplitudes of the signals are comparable. A Min-Max normalization has been selected, which allows scaling the signal amplitude data between 0 and 1. The following formula has been applied for this purpose:

$$x_{norm}={\displaystyle \frac{x-x_{min}}{x_{max}-x_{min}}}$$

(1)

Let $x_{norm}$ be the value of the normalized signal, $x_{min}$ the minimum amplitude value of the signal, and $x_{max}$ the maximum amplitude value of the signal.

Each movement in the X or Y axis results in two signals, one for each eye. To combine these signals into a single signal per axis, the average of the values from both eyes is calculated for the X axis and the Y axis for each frame. Once the movement signals are obtained, their impact on the plethysmographic signal is evaluated. This signal lies within a specific frequency range, between 0.67 and 4 Hz, corresponding to the human pulse interval. To determine if the previous movement signals interfere within this frequency range, they are converted to the frequency domain using a Fast Fourier Transform (FFT).

To eliminate the effect of movement, the final step will be to subtract this noise signal from the frequency spectrum of the plethysmographic signal of interest. This approach will help to attenuate the peaks generated by the noise and enhance those that contain relevant information for the final estimation.

$$Resulting Signal = Plethysmographic Signal - (factor*noise)$$

(2)

However, the reduction of noise will depend on each specific signal, requiring more or less attenuation as needed. To determine the amount of reduction required, an Artificial Intelligence (AI) model has been developed to predict the attenuation factor based on the plethysmographic signal. Since the input of the model is the plethysmographic signal in the frequency spectrum, an architecture capable of relating a vector of values to infer patterns and features that help predict the required reduction factor was needed. The chosen architecture was a Convolutional Neural Network (CNN), as it meets all the established requirements.

To extract the factor values used for the network’s learning, a Grid Search strategy was followed. This technique involves setting a range of parameter factor values and evaluating it based on the resulting errors of the method. The factor value with the lowest error for each of the input signals was set as the label or ground truth. Figure 4 shows the input signal values, the target factor value used for training the model, and the error associated with that configuration.

The model training is performed using the Leave-One-Out Cross-Validation (LOOCV) technique [5]. This strategy maximizes data usage and avoids partition bias, as each data point is used in both the training and test sets. The network has the following hyperparameters:

• Optimizer: Adam;
• Loss function: Mean Squared Error (MSE);
• Epochs: 20;
• Batch size: 16;

Figure 5 represents the comparison between the actual and predicted factor values by the CNN model for motion artifact removal:

3.1.3. Peaks Average

The preprocessing methods analyzed in the literature generally have the same goal: to process the signal so that, in the frequency spectrum, there is a peak with significantly higher amplitude than the others. This peak, when multiplied by 60, corresponds to the heart rate. The multiplication is necessary because the plethysmographic signal being used is filtered between 0.67 and 4 Hz, which corresponds to the normal human pulse frequency in beats per second. To convert it to the international unit, it needs to be transformed into beats per minute (bpm).

$$Heart Rate = PeakFreq*60\left[bpm\right]$$

(3)

In summary, the frequency with the highest amplitude is chosen, converted to beats per minute (bpm), and selected as the heart rate. This choice is based on the assumption that this frequency represents the subject’s pulse, as it has a stationary component that pulses with considerable consistency. When analyzing the signal in the frequency domain, this consistency manifests as a prominent peak compared to other values [12].

This ideal scenario is rarely encountered in practice, as various sources of noise, such as lighting, movement, and camera quality, can affect the signal. When noise is present, the signal no longer shows isolated and well-defined pulses, and other smaller peaks emerge around the main peak. In such cases, one can either average the frequencies of these peaks to approximate the true result or use the most prominent peak directly. The choice of method depends on the characteristics of the signal, but both approaches can keep the error within acceptable limits. However, the problem arises when, after filtering the signal within the range of interest, multiple high-amplitude peaks appear at disparate frequencies within this band, caused by noise. To address this challenge, an iterative process has been developed to progressively identify the peaks in the signal using a threshold that is iteratively adjusted. Initially, the threshold is set to zero and is progressively adjusted as peaks are detected within a defined range of peak mean ± standard deviation.

Figure 6. Peak Detection Process

Figure 6. Peak Detection Process

The process stops once a predetermined number of peaks have been detected [11]. The higher this value, the more peaks will be considered for averaging, which will affect the final result. To optimize the stopping threshold value, an artificial intelligence (AI) model has been developed to predict the appropriate value, as the number of peaks to be considered will depend on the specific characteristics of the signal. The previous figure shows a real example in which the model has estimated that the most optimal value for peak selection is 7. This averaging allows the prediction to properly adjust to the actual value, resulting in a heart rate of 77.3 bpm, which is very close to the actual value of 79.4 bpm.

Similar to the previous motion noise reduction model, the input to the model consists of the discrete amplitude values of the plethysmographic signal in the frequency domain. Since the input data have different lengths, it is necessary to normalize them to a uniform size. For this task, a Convolutional Neural Network (CNN) is used, capable of identifying complex patterns and accurately estimating the result.

Like the AI model described in the previous section, this one also trains using the Grid Search technique. Figure 4 shows the input data: on one hand, the plethysmographic signal in the range of interest within the frequency domain, and on the other, the target value, $N_{peaks}$ associated with each signal. This $N_{peaks}$ value represents the optimal number of peaks to consider to obtain the minimum error when estimating the heart rate, and it is therefore the value used as the label or ground truth when training the CNN.

Likewise, the training is carried out using the Leave-One-Out Cross-Validation (LOO-CV) strategy, due to the limitation in the amount of available data. The network has the following hyperparameters:

• K-fold: 72;
• Optimizer
• Loss function: Mean Squared Error (MSE);
• Epochs: 50;
• Batch size: 16;

Subsequently, in Figure 7, some comparisons between the actual values and those predicted by the model are presented.

3.3. Blood Pressure

To reduce the effect of motion artifacts as well as other factors that can negatively impact the acquisition of the plethysmographic signal (rPPG), from which blood pressure estimation is derived, this section outlines the set of preprocessing techniques implemented in the AIVA project.

Figure 11 shows a series of processes discussed in section () Preprocessing, as it is common to both the heart rate estimation method and the blood pressure estimation method. However, once the RGB signals are obtained separately, the next step differs from the heart rate model; in this case, a ‘detrending and moving average filtering’ is performed.

3.3.1. Detrending & Moving Average Filtering

According to what is mentioned in [12], the signals resulting from the decomposition of RGB channels are filtered using a 5-point moving average filter to smooth them, i.e., reduce noise and make the underlying trends more visible. After that, the stationary component is removed using the mean scaling and centering algorithm. The detrending and moving average filtering process is shown in the following equation [16]:

$$x\left(t\right)={\displaystyle \frac{x_i\left(t\right)-m(t,l)}{m(t,l)}}$$

(4)

Where m(t,l) = ${\displaystyle \frac{1}{L}}$${\sum }_{k=0}^{L-1}{x}_i(t-k)$ , with L representing the moving average of L points from the color vectors for channel i (in our case, $L=5$).

The primary objective of applying this equation is to eliminate the DC component of the signal in the resultant signals from the RGB channel decomposition, as well as to suppress high-frequency noise, thereby preparing the signal for subsequent analysis.

3.3.2. SPA Detrending

After performing the detrending of the signal as described in the previous section, a more flexible detrending process known as the Smoothness Prior Detrending Approach (SPA) is then applied. This technique, unlike the moving average method which uses a fixed set of points, offers improved accuracy by explicitly modeling the trend, even for complex and nonlinear trends. The procedure is as follows:

• Definition of the Smoothness Prior:
  a. In our case, the smoothness of the signal is penalized using the squared difference of the signal components, that is, ${\sum }_t{(x\left(t\right)-x\left(t\right))}^2$
• Definition of the Detrending Function
  a. The detrending function combines the fidelity to the original signal with the smoothness of the trend, and is defined as follows:

  $$x_{SPA\left(t\right)}=\sum _t{(x\left(t\right)-\widehat{x}\left(t\right))}^2-\alpha \sum _t{\left({\displaystyle \frac{d^{2}x\left(t\right)}{d{t}^2}}\right)}^2$$

  (5)

  i. The first term aims to remove the linear trend or slope of $x\left(t\right)$ by subtracting the original signal values $x\left(t\right)$ from the estimated values $\widehat{x}\left(t\right)$, where $\widehat{x}\left(t\right)$ is the result of modeling $x\left(t\right)$ as a linear regression.
  ii. The second term represents the square of the second derivative of the function and penalizes rapid variations in the signal.
  iii. $\alpha$ is the smoothing parameter that controls the balance between the fidelity to the signal and the smoothness of the fit, effectively preventing excessive smoothing of the function.

Figure 12. SPA Detrending Process

Figure 12. SPA Detrending Process

To determine the parameter $\alpha$ that provided the best performance, a grid search was conducted. Specifically, for a range of $\alpha$ values between 0.01 and 0.0001, the model’s performance was evaluated. The results demonstrated that $\alpha =0.001$ offers the best trade-off between signal fidelity and smoothness of the fit.

Table 1. Grid Search Results for the Performance Parameter $\alpha$

Table 1. Grid Search Results for the Performance Parameter $\alpha$

$\mathit{\alpha }$ value	0.001	0.08	0.001	0.0015	0.0001
Systolic Error	0.029019	0.012015	0.059270	0.056503	0.039644
Diastolic Error	0.063691	0.077609	0.032232	0.033700	0.051969

3.3.3. Fast ICA

After removing the trend from the signal using the proposed methods, the next step is to apply Blind Source Separation (BSS) techniques. These techniques enable the extraction of original source signals from a set of mixed signals, where the sources are combined. In other words, from mixed signals that result from a set of source signals, the goal is to recover the original sources from these mixtures, with minimal or no information about the sources or the mixing process.

There are various methods for blind source separation. In this paper, we focused on Independent Component Analysis (ICA) as it provided good results and is the most frequently used method in the literature. ICA is a computational method used to separate a mixed signal into additive subcomponents that can potentially reconstruct the original signal. The method starts with the mixed signal, assuming it has a Gaussian-like behavior, and applies techniques to reduce the Gaussianity of this mixed signal, under the assumption that the subcomponents are non-Gaussian and statistically independent. The following figure shows the distribution functions of the amplitudes of the R, G, and B channel signals before applying FastICA, and it is observed that they exhibit some Gaussian behavior.

Figure 13. Gaussianity of Signals Before Applying Fast ICA

Figure 13. Gaussianity of Signals Before Applying Fast ICA

Several versions of ICA exist; however, we have focused on the Fast ICA algorithm [17]. This algorithm is an optimized version that makes the ICA process faster and more efficient. In our case, Fast ICA is applied to the extracted signal RGB channels (3.1), resulting in 3 components. To determine which component represents the PPG signal and which are sources of noise, each of these components is transformed into the frequency domain using FFT and filtered within the range of interest using a Chebyshev Type II filter, as it is a low-pass filter characterized by an attenuation band with a ripple zone and a flat passband with a rapid transition. The component with the most prominent peak in the frequency band of interest is selected as the one containing the photoplethysmography information.

Figure 14. Independent components

Figure 14. Independent components

3.3.4. Peak and Valley Detection

After applying FastICA and selecting the component with the maximum peak within the frequency spectrum, this stage involves counting the peaks and valleys of the signal.

According to the Radial Resonance Theory [18], the arterial system can be considered similar to a blood wave transmission system. This approximation allows for the development of a model capable of predicting blood pressure based on equations and properties of the arteries. Using the captured frames, the arterial pulse signal can be filtered, and the peaks and valleys, which are considered systolic and diastolic pressures, can be continuously estimated.

Along with the periodic cardiac contraction activities, the light incident on the face is absorbed and attenuated by the skin, muscles, and blood, so that the intensity of the reflected light detected by the optical receiver is weakened. During heart contraction, the volume of peripheral blood is greater, leading to significant light absorption, and the detected light intensity is at its minimum; whereas when the heart is in a diastolic state, in contrast, the detected light intensity is at its maximum.

Figure 15. Blood Pressure signal

Figure 15. Blood Pressure signal

Using the BMI measure [19] as a correction coefficient and fitting an equation described in the Radial Resonance Theory, where the heart is described as a cyclic pressure pump, the following equation is obtained:

$$f(x,bmi)=a_0+a_1*x+a_2*bmi+a_3*x*bmi$$

(6)

Where $f(·)$ is the function to be fitted; $bmi$ is the BMI correction coefficient; and x is the average of all peaks or valleys in the signal. The averages of peaks and valleys are extracted using the following equations:

$$E_{peak}=\sum _{n_1}H{D}_{n1},E_{valley}=\sum _{n_2}H{L}_{n2}$$

(7)

Where $n_1$ and $n_2$ represent the number of peaks and valleys, respectively; $HD$ and $HL$ are the amplitude values of the peaks and valleys. When counting the number of peaks and valleys, those whose separation is less than the following equation have been omitted:

$$minPeakDistance=bps*fps$$

(8)

Where $bps$ represents the beats per second of the subject, derived from the heart rate, and $fps$ represents the frames per second at which the video is captured.

According to the experiments defined in [20], the parameters $a_i(i=0,1,2,3)$ of the function $f(·)$ have been determined using a least squares fit of the blood pressure measurements present in our dataset, resulting in the following equations for systolic blood pressure (SBP) and diastolic blood pressure (DBP):

$$SBP=23.7889+95.4335*E_{peak}+4.5958*bmi-5.109*E_{peak}*bmi$$

(9)

$$DBP=-17.3772-115.1747*E_{valley}+4.0251*bmi-5.2825*E_{valley}*bmi$$

(10)

4. Results and Discuss

Figure 16. Comparison of Results for Heart Rate

Figure 16. Comparison of Results for Heart Rate

Figure 21. Histogram of Systolic Blood Pressure Errors

Figure 21. Histogram of Systolic Blood Pressure Errors

Figure 22. Histogram of Diastolic Blood Pressure Errors

Figure 22. Histogram of Diastolic Blood Pressure Errors

Figure 23. Bland-Altman Blood Pressure

Figure 23. Bland-Altman Blood Pressure

5. Conclusions

References

1. Poh, M. Z. , McDuff, D. J., & Picard, R. W. Non-contact, automated cardiac pulse measurements using video imaging and blind source separation. Optics express, 2010, 18, 10762–10774. [Google Scholar]
2. Kado, S. , Monno, Y., Moriwaki, K., Yoshizaki, K., Tanaka, M., & Okutomi, M. Remote heart rate measurement from RGB-NIR video based on spatial and spectral face patch selection. In 2018 40th annual international conference of the IEEE engineering in medicine and biology society (EMBC), 2018, 5676–5680.
3. Yang, F.; He, S.; Sadanand, S.; Yusuf, A.; Bolic, M. Contactless measurement of vital signs using thermal and RGB cameras: A study of COVID 19-related health monitoring. Sensors 2022, 22, 627. [Google Scholar] [CrossRef] [PubMed]
4. Secerbegovic, A. , Bergsland, J., Halvorsen, P. S., Suljanovic, N., Mujcic, A., & Balasingham, I. Blood pressure estimation using video plethysmography. In 2016 IEEE 13th international symposium on biomedical imaging (ISBI), 2016, 461–464.
5. Iwashita, Y. , Nagumo, K. , Oiwa, K., & Nozawa, A. Estimation of resting blood pressure using facial thermal images by separating acute stress variations. Artificial Life and Robotics, 2021, 26, 473–480. [Google Scholar]
6. Luo, H. Yang, D., Barszczyk, A., Vempala, N., Wei, J., Wu, S. J., ... & Feng, Z. P. Smartphone-based blood pressure measurement using transdermal optical imaging technology. Circulation: Cardiovascular Imaging,2019, 212(8).
7. Wu, J. , Liang, H., Ding, C., Huang, X., Huang, J., & Peng, Q. Improving the accuracy in classification of blood pressure from photoplethysmography using continuous wavelet transform and deep learning. International journal of hypertension, 2021,2021.
8. Google (2010). Mediapipe [Computer Software]. Google. Available online: https://ai.google.dev /edge/mediapipe/solutions/guide?hl=es-419.
9. Bousefsaf, F. , Maaoui, C., & Pruski, A. Continuous wavelet filtering on webcam photoplethysmographic signals to remotely assess the instantaneous heart rate. Biomedical Signal Processing and Control, 2013, 8(6), 568–574.
10. S. Bobbia, R. S. Bobbia, R. Macwan, Y. Benezeth, A. Mansouri, J. Dubois. Unsupervised skin tissue segmentation for remote photoplethysmography. Pattern Recognition Letters, 2017.
11. Wu, B. F. , Huang, P. W., Lin, C. H., Chung, M. L., Tsou, T. Y., & Wu, Y. L. Motion resistant image-photoplethysmography based on spectral peak tracking algorithm, 2018, 6, 21621–21634. [Google Scholar]
12. A. Bhattacharjee and M. S. U. Yusuf. A Facial Video based Framework to Estimate Physiological Parameters using Remote Photoplethysmography, 2021, 1–7.
13. S. Hara et al. Optimización de parámetros del sensor de frecuencia cardíaca basado en PPG con cancelación de artefactos de movimiento mediante validación cruzada, 2017, 6, 73–76.
14. Sánchez, G. L. Sánchez, L. L., & Suárez, A. D. Composición corporal y variabilidad de la frecuencia cardiaca: relaciones con edad, sexo, obesidad y actividad física. SPORT TK-Revista EuroAmericana de Ciencias del Deporte, 2015, 4(2), 33–40.
15. Pedregosa, F. , Varoquaux, Ga"el, Gramfort, A., Michel, V., Thirion, B., Grisel, O., … others. (2011). Scikit-learn: Machine learning in Python. Journal of Machine Learning Research, 12(Oct), 2825–2830.
16. Gerard de Haanand Vincent Jeanne. "Robust Pulse Rate from Chrominance-Based rPPG”, IEEE Transactions on Biomedical Engineering, 2013, 2873–2886.
17. A. Hyvärinen and E. Oja. A Fast Fixed-Point Algorithm for Independent Component Analysis. Neural Computation, 1997, 9(7), 1483-1492.
18. Wang, Y. Y. L. , Jan, M. Y., Wang, G. C., Bau, J. G., & Wang, W. K. Pressure pulse velocity is related to the longitudinal elastic properties of the artery. Physiological Measurement, 2004, 25(6), 1397.
19. M. Higuchi, K. Sorachi, and Y. Hata. Health checkup data analysis focusing on body mass index, 2017, 8, 1634–1641.
20. Yimin Zhou, Member, IEEE, Haiyang Niy, Qi Zhang and Qingtian Wu, Member, IEEE. The Noninvasive Blood Pressure Measurement based on Facial Images Processing, 2015, 14, 1–11.