Heart rate estimation method and device and electronic equipment

Abstract

The invention is suitable for the technical field of remote physiological parameter monitoring, and provides a heart rate estimation method and device and electronic equipment, and the method comprises the steps: obtaining a continuous video image sequence of the face of a to-be-detected object; determining RGB time sequence signals changing along with time according to the video image sequence; blind source separation is conducted on the RGB time sequence signals through an independent component algorithm, and independent component signals related to the heart rate are obtained; performing time-frequency analysis on the independent component signals related to the heart rate by adopting at least two time-frequency analysis methods to obtain corresponding time-frequency analysis results; and estimating the heart rate of the to-be-tested object according to each time-frequency analysis result. Through the method, the accuracy of the estimated heart rate can be improved.

Description

Technical Field

[0001] This application belongs to the field of remote physiological parameter monitoring technology, and particularly relates to heart rate estimation methods, devices, electronic devices, computer-readable storage media, and computer program products. Background Technology

[0002] Heart rate refers to the number of times the heart beats per minute, and it is a basic physiological parameter for assessing cardiovascular function. Heart rate can be detected by methods such as palpation and auscultation, or it can be detected using the Independent Component Analysis (ICA) algorithm.

[0003] When using the ICA algorithm to detect heart rate, the process typically involves the following: A continuous sequence of video images of the subject's face is acquired non-contactly using a camera device. Based on this sequence, an RGB time-series signal is determined, containing information about subtle changes in light intensity caused by blood flow variations. The ICA algorithm is then applied to perform blind source separation on the RGB time-series signal to obtain independent components related to heart rate. These independent components are then analyzed using a single time-frequency analysis method, and the corresponding heart rate (or heart rate value) is determined based on the analysis results.

[0004] However, when using the above method to determine heart rate in non-stationary scenarios, the accuracy of the determined heart rate is low. Summary of the Invention

[0005] This application provides a heart rate estimation method, apparatus, and electronic device, which can solve the problem of low accuracy in heart rate estimation by existing methods.

[0006] In a first aspect, embodiments of this application provide a heart rate estimation method, including:

[0007] Acquire a continuous video image sequence of the face of the subject under test; The RGB timing signal that varies over time is determined based on the video image sequence; The independent component algorithm is used to perform blind source separation on the RGB time-series signal to obtain independent component signals related to heart rate; At least two time-frequency analysis methods are used to perform time-frequency analysis on the independent component signal related to heart rate to obtain the corresponding time-frequency analysis results. Each of the time-frequency analysis methods is based on a different principle when performing the time-frequency analysis. The heart rate of the subject under test is estimated based on the time-frequency analysis results.

[0008] The beneficial effects of the embodiments in this application compared with the prior art are: In this embodiment, at least two time-frequency analysis methods are used to perform time-frequency analysis on the independent component signals related to heart rate, and then the heart rate of the subject is estimated based on the obtained time-frequency analysis results. Since each time-frequency analysis method is based on a different principle, the errors of the time-frequency analysis results obtained after the analysis are also independent of each other, which helps to reduce the error of the heart rate of the subject estimated based on each time-frequency analysis result, thereby improving the accuracy of the obtained heart rate.

[0009] Secondly, embodiments of this application provide a heart rate estimation device, comprising: The video image sequence acquisition module is used to acquire a continuous video image sequence of the face of the object under test; A timing signal acquisition module is used to determine the RGB timing signal that changes over time based on the video image sequence; The independent component signal determination module is used to perform blind source separation on the RGB time-series signal using an independent component algorithm to obtain independent component signals related to heart rate. The time-frequency analysis module is used to perform time-frequency analysis on the independent component signal related to heart rate using at least two time-frequency analysis methods to obtain the corresponding time-frequency analysis results, wherein each of the time-frequency analysis methods is based on a different principle when performing the time-frequency analysis; The heart rate determination module is used to estimate the heart rate of the subject under test based on the various time-frequency analysis results.

[0010] Thirdly, embodiments of this application provide an electronic device, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the method described in the first aspect.

[0011] Fourthly, embodiments of this application provide a computer-readable storage medium storing a computer program that, when executed by a processor, implements the method described in the first aspect.

[0012] Fifthly, embodiments of this application provide a computer program product that, when run on an electronic device, causes the electronic device to perform the method described in the first aspect.

[0013] It is understood that the beneficial effects of the second to fifth aspects mentioned above can be found in the relevant descriptions in the first aspect mentioned above, and will not be repeated here. Attached Figure Description

[0014] To more clearly illustrate the technical solutions in the embodiments of this application, the accompanying drawings used in the description of the embodiments or the prior art will be briefly introduced below.

[0015] Figure 1 This is a flowchart illustrating a heart rate estimation method provided in an embodiment of this application; Figure 2 This is a schematic diagram of the structure of a heart rate estimation device provided in an embodiment of this application; Figure 3 This is a schematic diagram of the structure of an electronic device provided in an embodiment of this application. Detailed Implementation

[0016] In the following description, specific details such as particular system architectures and techniques are set forth for illustrative purposes and not for limitation, in order to provide a thorough understanding of the embodiments of this application. However, those skilled in the art will understand that this application may also be implemented in other embodiments without these specific details. In other instances, detailed descriptions of well-known systems, apparatuses, circuits, and methods have been omitted so as not to obscure the description of this application with unnecessary detail.

[0017] It should be understood that, when used in this application specification and the appended claims, the term "comprising" indicates the presence of the described features, integrals, steps, operations, elements and / or components, but does not exclude the presence or addition of one or more other features, integrals, steps, operations, elements, components and / or a collection thereof.

[0018] It should also be understood that the term “and / or” as used in this application specification and the appended claims means any combination of one or more of the associated listed items and all possible combinations, and includes such combinations.

[0019] Furthermore, in the description of this application and the appended claims, the terms "first," "second," etc., are used only to distinguish descriptions and should not be construed as indicating or implying relative importance.

[0020] References to "one embodiment" or "some embodiments" in this specification mean that one or more embodiments of this application include a specific feature, structure, or characteristic described in connection with that embodiment. Therefore, the phrases "in one embodiment," "in some embodiments," "in other embodiments," "in still other embodiments," etc., appearing in different parts of this specification do not necessarily refer to the same embodiment, but rather mean "one or more, but not all, embodiments," unless otherwise specifically emphasized.

[0021] Heart rate refers to the number of times the heart beats per minute; it is also an instantaneous expression of the oscillation frequency of the cardiovascular system. Currently, heart rate is usually measured using non-contact methods to improve the convenience of testing.

[0022] Heart rate detection methods combining the ICA algorithm are commonly used non-contact methods. However, when using the ICA algorithm for heart rate detection, the accuracy of the determined heart rate is low if a single time-frequency or single frequency domain analysis is performed on the independent component signals related to heart rate. This is because single time-frequency analysis cannot effectively handle the time-varying characteristics of heart rate, and therefore its performance degrades significantly in non-stationary scenarios such as exercise, resulting in low accuracy of the determined heart rate.

[0023] To improve the accuracy of heart rate determination based on the ICA algorithm, this application provides a heart rate estimation method.

[0024] In this heart rate estimation method, at least two time-frequency analysis methods are used to perform time-frequency analysis on the independent component signals related to heart rate, and then the heart rate of the subject is estimated based on the obtained time-frequency analysis results. Since each time-frequency analysis method is based on a different principle, the errors of the time-frequency analysis results obtained after the analysis are also independent of each other, which helps to reduce the error of the heart rate of the subject estimated based on each time-frequency analysis result, thereby improving the accuracy of the obtained heart rate.

[0025] The heart rate estimation method provided in the embodiments of this application is described below with reference to the accompanying drawings.

[0026] Figure 1 A flowchart illustrating a heart rate estimation method provided in an embodiment of this application is shown. This heart rate estimation method can be applied to electronic devices, and is described in detail below: S11, acquire a continuous video image sequence of the face of the object under test.

[0027] Specifically, an RGB camera can be used to capture a user's face to obtain a continuous time-frequency image sequence (or video stream) in a color space that includes red (R), green (G), and blue (B). An RGB camera is a device that generates color images by capturing the brightness values ​​of the R, G, and B color channels.

[0028] Optionally, the color space of the aforementioned video stream may also include the Hue, Saturation, Value (HSV) color space and the Luma (Y) + Chroma-red (Cr) + Chroma-blue (Cb) (YCrCb) color space, which is not limited here.

[0029] Optionally, the acquisition frame rate of the above continuous time-frequency image sequence is not less than 30 frames per second (30 fps) to meet the requirements of heart rate measurement range and motion artifact suppression.

[0030] S12, determine the RGB timing signal that changes over time based on the above video image sequence.

[0031] Optionally, considering that the density of blood vessels in the face and occlusion conditions will affect the subsequently determined RGB time-series signal, the region of interest detection of the video image sequence can be performed based on at least one of the following conditions in the time-frequency image sequence: blood vessel density, occlusion conditions, illumination conditions, skin color conditions (which can be determined based on whether a skin color region is included), and signal quality (which can be determined based on the signal-to-noise ratio of the video image sequence).

[0032] For each region of interest in the above region of interest sequence, calculate the spatial average value of the region of interest in each channel, and determine the time-varying RGB timing signal corresponding to each channel based on the spatial average value of each region of interest in each channel.

[0033] S13, use the independent component algorithm to perform blind source separation on the above RGB time-series signal to obtain the independent component signal related to heart rate.

[0034] In this embodiment, a mixed signal matrix can be constructed based on the RGB time-series signal that changes over time. Then, based on the mixed signal matrix and a preset objective function, the separation matrix to be iterated is iterated until one of the following iteration stopping conditions is met to obtain the target separation matrix. Based on the target separation matrix, the mixed signal matrix is ​​converted into an independent component matrix, and the independent component signal related to heart rate is determined based on the independent component matrix.

[0035] Optionally, the objective function includes the negative entropy of the separated signal and a heart rate signal constraint term, which is used to ensure that the separated signal conforms to the physiological laws of heart rate. Since the objective function includes not only the negative entropy of the separated signal but also the heart rate signal constraint term, when the mixed signal matrix is ​​converted into an independent component matrix using the objective separation matrix obtained from the objective function, the separated signals of this independent component matrix not only possess independence but also satisfy the physiological laws of heart rate, thereby improving the accuracy of the heart rate calculated from this independent component matrix.

[0036] Optionally, the above-mentioned iteration stopping conditions include: the above-mentioned objective function reaches a preset convergence threshold or the number of iterations reaches a preset upper limit; Optionally, considering that among the multiple separated signals contained in the independent component matrix, usually only one is a heart rate signal, while the other separated signals are usually signals corresponding to interference components such as noise, motion artifacts, and changes in illumination, the score of the heart rate signal constraint term corresponding to each separated signal in the independent component matrix can be calculated. Based on the scores of the heart rate signal constraint terms corresponding to each separated signal, the separated signal that best matches the physiological characteristics of heart rate (such as taking the separated signal with the highest score of the heart rate signal constraint term as the heart rate signal) is selected as the independent component signal related to heart rate.

[0037] S14, At least two time-frequency analysis methods are used to perform time-frequency analysis on the independent component signals related to heart rate, and the corresponding time-frequency analysis results are obtained. Each time-frequency analysis method is based on a different principle when performing time-frequency analysis.

[0038] In this embodiment, the aforementioned principle refers to the principle of time-frequency analysis. This principle may include: Fourier transform based on a fixed window, wavelet basis function decomposition with variable scale, adaptive empirical mode decomposition and Hilbert transform, Fourier transform of the instantaneous autocorrelation function of the signal, etc. Of course, other principles may also be included, but are not limited here.

[0039] Optionally, three different time-frequency analysis methods are used to perform time-frequency analysis on the independent component signal to achieve time-frequency analysis of the fused independent component signal. The three methods are: Short-Time Fourier Transform (STFT), Continuous Wavelet Transform (CWT), and Hilbert-Huang Transform (HHT).

[0040] Among them, STFT is based on fixed-window short-time Fourier transform for time-frequency analysis of the signal, CWT is based on variable-scale wavelet basis function decomposition for time-frequency analysis of the signal, and HHT is based on adaptive empirical mode decomposition and Hilbert transform for time-frequency analysis of the signal. Since the three methods are based on different principles, their errors are independent of each other, which helps to reduce the error after fusion and improve the calculated heart rate. At the same time, when a single method fails (such as STFT failing during rapid changes), other methods can still provide effective estimates, thereby improving the robustness of subsequent heart rate techniques.

[0041] Assuming that each method yields a time-frequency analysis result after performing time-frequency analysis on the independent component signal, and these three time-frequency analysis results are respectively the first time-frequency analysis result, the second time-frequency analysis result, and the third time-frequency analysis result, then the above-mentioned time-frequency analysis results obtained by using at least two time-frequency analysis methods to perform time-frequency analysis on the heart rate-related independent component signal include: A1. The Short-Time Fourier Transform (STFT) method was used to perform time-frequency analysis on the above independent component signals related to heart rate, and the first time-frequency analysis result was obtained.

[0042] In this embodiment, a first time-frequency analysis result is obtained by performing a short-time Fourier transform on the independent component signal related to heart rate according to a preset window length. This first time-frequency analysis result can be the short-time Fourier transform time-frequency spectrum STFT(f,t), or the energy value at time t and frequency f calculated based on the STFT(f,t). The STFT described above uses a fixed window length, which can quickly generate a uniform time-frequency distribution across the entire signal. It has high computational efficiency and is therefore suitable for rough localization in the early stages of heart rate analysis, providing a priori range for subsequent fine estimation.

[0043] Optionally, considering that the first time-frequency analysis result is related to the window length, and that the first time-frequency analysis result will be used to calculate the heart rate, the aforementioned window length can be calculated based on relevant heart rate information. That is, before performing time-frequency analysis on the independent component signal related to heart rate using the short-time Fourier transform method, the following steps are also included: Based on the physiological characteristics of heart rate and the requirements for time-frequency resolution and / or statistical stability, the window length range of the above-mentioned short-time Fourier transform is determined; based on the above-mentioned window length range, the window length of the above-mentioned short-time Fourier transform is determined.

[0044] The aforementioned physiological characteristics of heart rate include: heart rate range, minimum frequency resolution to distinguish adjacent heart rates, heart rate cycle corresponding to the heart rate range, and heart rate changes slowly at rest, but may change significantly within seconds during exercise or emotional fluctuations.

[0045] The aforementioned statistical stability includes including a sufficient number of sampling points within the window length to ensure statistical stability of the power spectrum estimation and reduce the impact of noise on the frequency estimation.

[0046] In this embodiment, the heart rate range is 40~200 bpm (0.67-3.33 Hz), and the minimum frequency resolution for distinguishing adjacent heart rates is approximately 0.1~0.2 Hz. According to the Fourier transform frequency resolution formula Δf = 1 / T (where T is the window length), a window length of 2 seconds provides a resolution of 0.5 Hz, and 5 seconds provides a resolution of 0.2 Hz, which is sufficient to accurately distinguish different heart rate values. Therefore, a window length range of 2~5 seconds can meet the requirement of minimum frequency resolution for distinguishing adjacent heart rates.

[0047] A heart rate range of 40-200 bpm corresponds to a heart rate cycle of 0.3-1.5 seconds. A window length of 2 seconds can contain at least 1-6 complete heart rate cycles, and 5 seconds can contain 3-16 heart rate cycles. Therefore, a window length range of 2-5 seconds is sufficient for reliable frequency estimation.

[0048] Heart rate changes slowly at rest, but can change significantly within seconds during exercise or emotional fluctuations. A window length of 2 to 5 seconds can capture this time-varying characteristic, providing sufficient temporal resolution to detect gradual changes in heart rate.

[0049] When the video frame rate is 30fps, the window length range of 2 to 5 seconds corresponds to 60 to 150 frames, which is sufficient to provide enough sampling points, thereby making the power spectrum estimation statistically stable and reducing the impact of noise on frequency estimation.

[0050] Furthermore, considering the uncertainty of time-frequency resolution in Short Time Fourier Transform (STFT): the longer the window length, the higher the frequency resolution but the lower the time resolution; the shorter the window length, the higher the time resolution but the lower the frequency resolution. Moreover, an excessively long window length (e.g., >10 seconds) will lead to a decrease in real-time performance, failing to meet the requirement of outputting heart rate within 1 second, while an excessively short window length (e.g., <1 second) will lead to insufficient frequency resolution and increase the computational load (requiring more frequent STFT calculations). A window length of 2-5 seconds achieves the best balance between real-time performance and computational efficiency, with a single STFT calculation time of approximately 10-20ms (based on an optimized FFT algorithm), meeting the requirements for real-time processing.

[0051] In summary, the window length can be set to a range of 2 to 5 seconds.

[0052] Once the window length range is determined, any value within that range can be used as the window length for the short-time Fourier transform. Alternatively, the specific window length can be determined within that range, taking into account the current environment of the test subject. Specifically, considering that the current environment of the test subject can affect their heart rate, the window length for the short-time Fourier transform can be determined within that range based on the current environment. That is, determining the window length for the short-time Fourier transform based on the aforementioned window length range includes: The current scene of the object under test is determined, including: a resting scene, a moving scene, and a transition scene from resting to moving; within the window length range, the corresponding window length of the short-time Fourier transform is determined according to the scene.

[0053] The resting scenario described above includes situations where the subject is in an office or during sleep monitoring. Considering that the heart rate changes slowly when the subject is in a resting scenario, the requirement for time resolution is relatively low. Therefore, in order to improve frequency resolution and measurement accuracy, a longer window length can be selected, such as setting the window length to a value between 4 and 5 seconds.

[0054] The aforementioned exercise scenarios include situations where the subject is in a gym or undergoing rehabilitation training. Considering that the subject's heart rate changes rapidly during exercise, a shorter window length can be selected to improve temporal resolution and track rapidly changing heart rates. For example, a window length of 2-3 seconds can be used. While the frequency resolution decreases slightly in this case, it remains within an acceptable range.

[0055] When the subject is transitioning from rest to exercise, the window length can be dynamically adjusted. Specifically, the heart rate change rate ΔHR / Δt of the previous window is analyzed. If the heart rate change rate ΔHR / Δt exceeds the preset first change rate threshold (e.g., 5 bpm / s), the window length can be shortened to 2 seconds; if the heart rate change rate ΔHR / Δt is lower than the preset second change rate threshold (e.g., 1 bpm / s), the window length can be extended to 5 seconds.

[0056] Optionally, this window length can be used as a reference scale for subsequent wavelet transform and Hilbert-Huang transform.

[0057] A2. The above independent component signals related to heart rate were analyzed in time and frequency using the Continuous Wavelet Transform (CWT) method to obtain the second time and frequency analysis result.

[0058] In this embodiment, a scaling parameter is used to stretch or compress the wavelet basis function, and then the stretched or compressed wavelet basis function is used to perform an inner product with the independent component signal related to heart rate, thereby obtaining the local features of the independent component signal related to heart rate at different times and frequencies.

[0059] Optionally, considering that the wavelet basis function is used to process the independent component signal related to heart rate, the wavelet basis function can be generated based on heart rate information to improve the accuracy of the subsequent second time-frequency analysis result. That is, before performing time-frequency analysis on the aforementioned independent component signal related to heart rate using the continuous wavelet transform method, the following steps are also included: The wavelet basis functions of the above continuous wavelet transform are generated based on the heart rate information.

[0060] Since the wavelet basis function is generated based on the heart rate information, when performing time-frequency analysis on the independent component signal related to heart rate based on the wavelet basis function, the components related to heart rate can be enhanced and the components unrelated to heart rate can be suppressed, thereby improving the accuracy of the obtained second time-frequency analysis results.

[0061] Optionally, the wavelet basis function of the continuous wavelet transform is generated based on at least one of the following information about the heart rate: frequency, morphology, the heart rate corresponding to the heart rate, and the systolic peak of the heart rate corresponding to the heart rate.

[0062] The aforementioned frequency can be selected from typical heart rate frequencies. For example, considering that the test subject is usually in a resting environment, the frequency of the test subject's heart rate in a resting environment can be used as the center frequency for generating the wavelet basis function. When a typical frequency is selected as the center frequency of the wavelet basis function, the wavelet can be aligned with the concentrated area of ​​heartbeat energy at the frequency band level, thereby improving the probability of the wavelet basis function matching the most typical heart rate.

[0063] In this embodiment, when generating the wavelet basis function based on the heart rate morphology, the waveform corresponding to the wavelet basis function can be made highly similar to the morphology of the heart rate signal. This increases the similarity between the "wavelet basis function" and the "real instantaneous heartbeat waveform," thereby maximizing the wavelet-filter's response to the heart rate component and minimizing its response to non-heart rate components. Optionally, the wavelet waveform can be made closer to the non-sinusoidal characteristics of the real remote photoplethysmography (rPPG) signal by introducing a second harmonic component.

[0064] In this embodiment, the standard deviation of the Gaussian window of the wavelet basis function can be determined based on the heartbeat cycle corresponding to the heart rate to control the degree of temporal localization. For example, when the standard deviation of the Gaussian window is determined based on the time of a heartbeat cycle, it is beneficial for the wavelet basis function to capture the waveform characteristics of a single heartbeat.

[0065] In this embodiment of the application, when the wavelet peak value of the wavelet basis function is determined based on the systolic peak value of the heartbeat corresponding to the heart rate, it is beneficial for subsequent detection of the wavelet peak value and calculation of the RR interval (i.e., the time distance between two adjacent R wave peaks in the electrocardiogram R waveform).

[0066] In this embodiment of the application, when a wavelet basis function for continuous wavelet transform is generated based on the heart rate frequency, morphology, the heartbeat cycle corresponding to the aforementioned heart rate, and the systolic phase of the heartbeat corresponding to the aforementioned heart rate, the wavelet basis function can be in the following form: ψ_HR(t) = A × exp(2πif_c t) × exp(-[(t-t_shift) / σ_t]²) × [1 + α×cos(4πf_c t)].

[0067] In the above, ψ_HR(t) is the wavelet basis function; A is the normalization constant to ensure that the wavelet energy is 1; f_c is the center frequency, which can be the frequency of the heart rate in the resting scene, such as f_c being 1.2Hz (72bpm); t_shift is the time shift, which is used to make the wavelet peak correspond to the systolic peak of the heartbeat waveform. For example, if the time shift is 0.2 seconds, the wavelet peak corresponds to the systolic peak of the heartbeat waveform, then t_shift can be set to 0.2 seconds; σ_t is the standard deviation of the Gaussian window, which is used to control the degree of time localization. When a heartbeat cycle is 0.8 seconds, σ_t can be 0.8 seconds to capture the details of a single heartbeat without over-segmentation that would reduce the frequency resolution; α is the harmonic modulation coefficient, and (cos(4πf_c t)) is the second harmonic component, making the wavelet waveform closer to the non-sinusoidal characteristics of the real rPPG signal.

[0068] The time-frequency resolution of continuous wavelet transform adaptively varies with scale: the low-frequency band (heart rate 40~60 bpm) corresponds to a large scale with high frequency resolution, capable of distinguishing subtle differences of only 1~2 bpm; the high-frequency band (120~200 bpm) corresponds to a small scale with high time resolution, capable of tracking sudden heart rate changes with millisecond-level lag. Based on this, the wavelet basis functions constructed in this embodiment directly use the real heartbeat waveform as a template, resulting in a high degree of morphological consistency, thus reducing the transform to "template matching"—maximizing the inner product of the heart rate component and the template, significantly lowering the detection threshold. Furthermore, the template introduces harmonic structures beyond the fundamental frequency, while noise lacks this harmonic correlation, thereby suppressing the transform response and making the system additionally robust to interference such as motion artifacts and baseline drift.

[0069] After determining the wavelet basis function, time-frequency analysis of the independent component signals related to heart rate can be performed by combining the scaling parameter.

[0070] Optionally, the range of the aforementioned scale parameters can be determined based on the frequency range corresponding to the heart rate, and then time-frequency analysis can be performed on the independent component signals related to the heart rate based on the parameters within the range of the scale parameters. That is, the aforementioned time-frequency analysis of the aforementioned independent component signals related to the heart rate using the continuous wavelet transform method yields the aforementioned second time-frequency analysis result, including: A21, determine the range of scale parameters based on the frequency range corresponding to the heart rate mentioned above.

[0071] A22, for each of the above-mentioned scale parameters within the above-mentioned scale parameter range, wavelet coefficients are calculated based on the above-mentioned scale parameters, the preset time shift, the wavelet basis function of the above-mentioned continuous wavelet transform, and the above-mentioned independent component signal related to heart rate.

[0072] A23. Calculate the corresponding wavelet power spectrum based on each of the above wavelet coefficients, and use the above wavelet power spectrum as the result of the second time-frequency analysis.

[0073] In this embodiment of the application, the range of the scale parameter a is defined as follows: a∈[a_min, a_max], corresponding to the frequency range[f_max, f_min].

[0074] Based on the inverse relationship between the center frequency and the scale parameter, we determine: f = f_c / a.

[0075] Where f is the frequency corresponding to the heart rate, then the frequency range corresponding to the heart rate is 0.67-3.33Hz: a_min = f_c / f_max = 1.2 / 3.33 ≈ 0.36; a_max = f_c / f_min = 1.2 / 0.67 ≈ 1.79.

[0076] If a logarithmic scaling distribution is used, then: a_i = a_min × (a_max / a_min)^(i / N), i=0,1,...,N, and optionally, N=64. A scale of 64 provides sufficient frequency resolution.

[0077] Assuming the time shift is represented by shift b_j, then for each scale a_i and time shift b_j, the wavelet coefficients can be calculated as follows: W(a_i, b_j) = (1 / √a_i) ∫ s_HR(t) × ψ_HR [(t-b_j) / a_i] dt; in: W(a_i, b_j) represents the wavelet coefficients; 1 / √a_i represents the normalization factor, ensuring that energies at different scales are comparable; s_HR(t) represents the independent component signal related to heart rate; ψ_HR This represents the conjugate of wavelet basis functions; (t-b_j) / a_i means shifting the wavelet by b_j and scaling it by a_i.

[0078] Calculate the wavelet power spectrum based on the wavelet coefficients mentioned above: P(a_i, b_j) = |W(a_i, b_j)|², converting scale a_i to frequency f for display: P(f_i, t_j), where f_i = f_c / a_i.

[0079] In this embodiment, the scale parameter is constrained a priori by predefining the frequency band range corresponding to the heart rate, thereby focusing the energy of time-frequency analysis on the cardiopulmonary coupling frequency band, which effectively suppresses irrelevant noise and artifacts while enhancing the cardiac component.

[0080] A3. The Hilbert-Huang Transform (HHT) method was used to perform time-frequency analysis on the above independent component signals related to heart rate, and the above third time-frequency analysis results were obtained.

[0081] The HWT mentioned above is an adaptive, nonlinear, non-stationary signal time-frequency analysis method, which consists of two steps: Empirical Mode Decomposition (EMD) and Hilbert spectral analysis.

[0082] The above-mentioned Hilbert-Huang transform was used to perform time-frequency analysis on the independent component signals related to heart rate, resulting in the third time-frequency analysis result, including: The above-mentioned independent component signal related to heart rate is adaptively decomposed into several intrinsic mode functions through the above-mentioned empirical mode decomposition; the above-mentioned Hilbert spectrum analysis is performed on each of the above-mentioned intrinsic mode functions to extract the instantaneous frequency of each of the above-mentioned intrinsic mode functions, and the above-mentioned instantaneous frequency is used as the result of the above-mentioned third time-frequency analysis.

[0083] In this embodiment, the heart rate-related independent component signal is expanded layer by layer into several Intrinsic Mode Functions (IMFs) according to local feature scale using EMD: s_HR(t) = Σ(i=1 to K) IMF_i(t) + r_K(t). Where s_HR(t) is the heart rate-related independent component signal, IMF_i(t) is the i-th IMF, representing an oscillation mode of the signal, K is the number of IMFs (which can be 8), and r_K(t) is the residual component (i.e., the trend term). The aforementioned IMF_i(t) can be obtained through the following filtering method: Initialization: h_0(t) = s_HR(t), For i = 1 to K: Repeat the following steps until h satisfies the IMF condition: Step a: Identify all local maxima and minima of h(t); Step b: Connect all maxima using cubic spline interpolation to form the upper envelope e_max(t); connect all minima using cubic spline interpolation to form the lower envelope e_min(t); Step c: Calculate the average envelope: m(t) = [e_max(t) + e_min(t)] / 2; Step d: Subtract the average envelope from the signal: h'(t) = h(t) - m(t); Step e: Determine whether h'(t) satisfies the IMF condition: Condition 1: The number of extreme points is equal to or differs by at most 1 from the number of zero-crossing points; Condition 2: At any given time, the average value of the upper and lower envelopes is 0; If IMF_i(t) = h'(t) is satisfied, proceed to step f; If not satisfied: Let h(t) = h'(t), return to step a; Step f: Calculate the residuals: r_i(t) = s_HR(t) - Σ(j=1 to i) IMF_j(t); Step g: If r_i(t) is monotonic or has fewer than 3 extreme points, stop the decomposition; Otherwise, let h_0(t) = r_i(t) and continue to decompose the next IMF.

[0084] After obtaining several IMFs, a Hilbert transform can be performed on each IMF component to obtain the instantaneous frequency and instantaneous amplitude, ultimately constructing a Hilbert spectrum with high time-frequency resolution. HWT, through EMD adaptive decomposition, is not limited by preset window functions and can capture non-stationary and nonlinear heart rate changes.

[0085] In some embodiments, considering that some of the decomposed IMFs may not contain heart rate signals, after obtaining several IMFs, it is also possible to first determine which IMFs contain heart rate signals, and then perform Hilbert transform on the IMFs containing heart rate signals. That is, before performing the Hilbert spectral analysis on each of the aforementioned intrinsic mode functions, the process further includes: B1. Calculate the dominant frequency of each of the above intrinsic mode functions.

[0086] B2. Based on the dominant frequency of each of the above intrinsic mode functions and the frequency range corresponding to the above heart rate, the above intrinsic mode functions with dominant frequencies within the frequency range corresponding to the above heart rate are selected to obtain the heart rate intrinsic mode functions.

[0087] Correspondingly, the Hilbert spectral analysis performed on each of the aforementioned intrinsic mode functions includes: The Hilbert spectral analysis was performed on each of the aforementioned intrinsic heart rate modes.

[0088] For the i-th IMF (i.e., IMF_i), calculate its dominant frequency: f_dominant_i = arg max_f PSD_i(f); The power spectral density of IMF_i is obtained by calculating the power spectral density of PSD_i using FFT.

[0089] If f_dominant_i ∈ [0.67Hz, 3.33Hz], meaning f_dominant_i is within the heart rate range, then IMF_i is labeled as a candidate heart rate component. Considering that heart rate information at any given moment is carried only by the "cleanest and most energy-concentrated" IMF, and other candidate IMFs are either low-energy or contain noise / harmonics / motion artifacts, the IMF with the highest energy or the highest signal-to-noise ratio is usually chosen as the intrinsic heart rate mode function IMF_HR to suppress the interference of noise modes on the heart rate trajectory and reduce computational load. IMF_HR = arg max_{IMF_i ∈ candidate heart rate components} Energy(IMF_i); Among them, Energy(IMF_i) = Σ_t [IMF_i(t)]².

[0090] In this embodiment, since the IMFs (i.e., intrinsic heart rate mode functions) within the heart rate range are first screened out by the main frequency, and then Hilbert spectrum analysis is performed on the screened intrinsic heart rate mode functions, it is beneficial to avoid motion artifacts, power frequency interference and other modes from entering the Hilbert spectrum, so that the instantaneous frequency curve no longer shows multi-value jumps. At the same time, since the number of intrinsic mode functions for Hilbert spectrum analysis is reduced, it is also beneficial to improve the calculation speed.

[0091] In some embodiments, considering that the instantaneous frequencies extracted after Hilbert spectral analysis may include frequencies outside the heart rate frequency range, these outliers can be filtered out. That is, performing the aforementioned Hilbert spectral analysis on each of the aforementioned intrinsic mode functions (or heart rate intrinsic mode functions) to extract the instantaneous frequencies of each of the aforementioned intrinsic mode functions includes: The Hilbert spectral analysis was performed on each of the aforementioned intrinsic mode functions, and outliers in the instantaneous frequencies of each of the aforementioned intrinsic mode functions extracted after the Hilbert spectral analysis were filtered out according to the frequency range corresponding to the aforementioned heart rate.

[0092] In this embodiment of the application, a Hilbert transform is performed on each IMF to extract the corresponding instantaneous frequency. Specifically: for the i-th (i from 1 to K) IMF, let z_i(t) = IMF_i(t) + j·H[IMF_i(t)], and then extract the instantaneous amplitude and instantaneous phase from the analytic signal z_i(t): A_i(t) = |z_i(t)| = √[IMF_i(t)² + H[IMF_i(t)]²]; φ_i(t) = arctan[H[IMF_i(t)] / IMF_i(t)].

[0093] Where z_i(t) is the analytic signal obtained after Hilbert transform, and H[·] is the Hilbert transform: H[IMF_i(t)] = (1 / π) × PV ∫_{-∞}^{∞} IMF_i(τ) / (t-τ) dτ, where PV represents the Cauchy principal value. Optionally, H[IMF_i(t)] is calculated using the FFT fast algorithm: H[IMF_i(t)] = IFFT[j·sign(f)·FFT(IMF_i(t))].

[0094] Where A_i(t) is the instantaneous amplitude and φ_i(t) is the instantaneous phase.

[0095] After obtaining the instantaneous phase, calculate the instantaneous frequency based on the instantaneous phase: f_i(t) = (1 / 2π) × dφ_i(t) / dt.

[0096] Of course, in actual calculations, numerical differentiation can also be used to obtain the result: f_i(t) = (1 / 2π) × [φ_i(t+Δt) - φ_i(t-Δt)] / (2Δt).

[0097] In this embodiment of the application, the frequency range corresponding to the heart rate can be set to [0.5Hz, 4Hz]. Assuming the instantaneous frequency is f_HR(t), if it is determined that f_HR(t) < 0.5Hz, or if it is determined that f_HR(t) > 4Hz, then f_HR(t) = NaN, that is, f_HR(t) is marked as invalid, which means that f_HR(t) is an outlier of the instantaneous frequency, and the outlier is filtered out.

[0098] Optionally, after filtering outliers, the corresponding normal value can be estimated to replace the outlier. This further includes, after filtering outliers of the instantaneous frequencies of the intrinsic mode functions extracted after the Hilbert spectral analysis based on the frequency range corresponding to the heart rate, the method also includes: Estimate the normal value corresponding to the above-mentioned outlier instantaneous frequency; replace the above-mentioned outlier instantaneous frequency with the estimated normal value to obtain the replaced instantaneous frequency sequence; and smooth the above-mentioned replaced instantaneous frequency sequence.

[0099] In this embodiment, NaN can be filled using cubic spline interpolation or linear interpolation, and then the filled instantaneous frequency sequence is smoothed. For example, a 1-second moving average window can be used to smooth the instantaneous frequency sequence: f_HR_smooth(t) = MovingAverage[f_HR(t), window=1 second]. Since the positions of the filtered outliers are filled after filtering out the outliers, it is equivalent to correcting the instantaneous frequency jumps, spikes, and blank segments in the Hilbert spectrum caused by noise, low perfusion, or motion artifacts, so that the final heart rate sequence maintains the details of physiological changes while having both continuity and low noise characteristics.

[0100] In this embodiment, true multi-resolution time-frequency analysis is achieved by fusing STFT, CWT, and HHT methods, which has the following advantages: a) Full-band coverage: STFT: Fixed resolution, unified processing across the entire frequency band; CWT: Low-frequency high-frequency resolution, high-frequency high-time resolution, adaptive coverage; HHT: Instantaneous frequency per sampling point, highest time resolution; By combining these three methods, at least one approach can provide the optimal representation of the heart rate, regardless of its frequency band.

[0101] b) Time-varying characteristic capture: From resting to exercise: The heart rate rises rapidly from 60 bpm to 120 bpm, and the instantaneous frequency of HHT can accurately track the rise process; Heart rate variability: subtle changes in the intervals between heartbeats (milliseconds), which can be captured by the high temporal resolution of HHT and CWT; Respiratory sinus arrhythmia: Heart rate changes periodically with respiration (0.2Hz modulation), and STFT can identify low-frequency modulation.

[0102] c) Robustness improvement: Noise interference: When one method is affected by noise, other methods may not be affected (e.g., STFT is sensitive to high-frequency noise, while HHT can be filtered out by EMD decomposition). Motion artifacts: Rapid motion causes signal non-stationarity; HHT and CWT have better adaptability than STFT. Sudden changes in illumination: These changes cause signal jumps, which can be suppressed by temporal smoothing and multi-method fusion.

[0103] d) Balance between accuracy and speed: STFT calculation is the fastest (10-20ms) and is preferred for preliminary estimation; CWT calculation is moderate (50-100ms), providing fine-grained analysis; HHT calculations are relatively slow (100-200ms), but they provide the most accurate estimates when needed (in motion scenarios). The three processes are computed in parallel (total time is about 200ms), which still meets the real-time requirements (outputs heart rate once per second).

[0104] S15, estimate the heart rate of the subject based on the above time-frequency analysis results.

[0105] In this embodiment of the application, the heart rate can be estimated based on each time-frequency analysis result, and then the heart rate of the subject to be tested can be obtained by fusing the various estimation results.

[0106] When performing time-frequency analysis using the three methods—short-time Fourier transform, continuous wavelet transform, and Hilbert-Huang transform—the estimation of the heart rate of the subject based on the results of each of these time-frequency analyses includes: C1. Estimate the candidate heart rate of the subject based on the first time-frequency analysis results to obtain the first candidate heart rate.

[0107] In this embodiment of the application, the first time-frequency analysis result includes the short-time Fourier transform coefficients STFT(f, t) of the independent component signal related to heart rate at time t and frequency f. Then, the first candidate heart rate (or the first candidate heart rate sequence) HR_STFT(t) can be: HR_STFT(t) = arg max_f STFT(f,t)×60, where arg max_f STFT(f,t) represents the frequency f that makes STFT(f,t) reach its maximum value among all frequencies at a fixed time t.

[0108] C2. Estimate the candidate heart rate of the subject based on the above second time-frequency analysis results to obtain the second candidate heart rate.

[0109] In this embodiment of the application, the second time-frequency analysis result includes the wavelet power spectrum |W(a,t)|², and the second candidate heart rate (or second candidate heart rate sequence) HR_CWT(t) obtained from the wavelet power spectrum is: HR_CWT(t) = arg max_f CWT(f,t)×60.

[0110] C3. Based on the results of the third time-frequency analysis above, estimate the candidate heart rate of the subject to be tested and obtain the third candidate heart rate.

[0111] When the third time-frequency analysis result is a smoothed instantaneous frequency sequence: f_HR_smooth(t), then the third candidate heart rate HR_HHT(t) is: HR_HHT(t) = f_HR_smooth(t) × 60.

[0112] C4. Estimate the heart rate of the subject under test based on the first candidate heart rate, the second candidate heart rate, the third candidate heart rate, and the pre-set candidate heart rate weights for each candidate heart rate, wherein each candidate heart rate weight is related to the current scene of the subject under test.

[0113] In this embodiment, the heart rate of the subject is obtained by multiplying each candidate heart rate by its corresponding candidate heart rate weight and then summing the results. Since the accuracy of the heart rate determined by the three methods—Short Time Fourier Transform (STFT), Continuous Wavelet Transform (CWT), and Hilbert-Huang Transform (HHT)—is related to the current scenario of the subject, in a resting scenario, STFT and CWT can provide high-precision heart rate, with HHT mainly serving as verification and supplement. However, in a motion scenario, the heart rate changes rapidly, and the instantaneous frequency tracking capability of HHT significantly improves measurement accuracy. Therefore, setting the weights of each candidate heart rate according to the current scenario of the subject can maintain optimal estimation under both "steady-state" and "transient" extreme physiological conditions. For example, increasing the weights of STFT / CWT at rest utilizes their high signal-to-noise ratio and computational efficiency to obtain a smooth, low-latency heart rate baseline; increasing the weights of HHT during sudden motion changes leverages its instantaneous frequency tracking advantage to quickly capture steep changes, avoiding lag and jumps.

[0114] In some embodiments, to further improve the accuracy of the set candidate heart rate weights, the candidate heart rate weights can also be calculated in the following manner: Assuming the candidate heart rate weights include a first candidate heart rate weight, a second candidate heart rate weight, and a third candidate heart rate weight, then the first candidate heart rate weight, the second candidate heart rate weight, and the third candidate heart rate weight can be determined using the following method: D1. Calculate the normalized peak energy based on the above first time-frequency analysis results, and calculate the weight of the first candidate heart rate based on the above normalized peak energy, the above first candidate heart rate, and the preset first base weight.

[0115] Assuming the weight of the first candidate heart rate is w_STFT(t), then w_STFT(t) can be: w_STFT(t) = α_STFT × ω_peak(t) × ω_stable(t).

[0116] Wherein, α_STFT represents the first basic weight, which is related to the current scene of the object under test. For example, it can be set to 0.5 in a resting scene and 0.3 in a moving scene; ω_peak(t) = E_peak(t) represents the normalized peak energy (E_peak / Σ_f STFT(f,t)), and the higher the value, the more reliable the frequency estimation; ω_stable(t) = exp(-|HR_STFT(t) -HR_STFT(t-1)| / σ_HR) represents the time stability, and σ_HR represents the smoothing parameter, which can be set to 10 bpm.

[0117] D2, and / or, determine the wavelet ridge frequency based on the second time-frequency analysis result, determine the ridge duration and matching degree based on the wavelet ridge frequency, and calculate the second candidate heart rate weight based on the ridge duration, the matching degree and the preset second basic weight, wherein the matching degree is the matching degree between the wavelet coefficient at the ridge and the wavelet basis function of the continuous wavelet transform generated based on the heart rate information.

[0118] Assuming the weight of the second candidate heart rate is w_CWT(t), then w_CWT(t) can be: w_CWT(t) = α_CWT × ω_ridge(t) × ω_match(t).

[0119] Wherein, α_CWT represents the second basic weight, which can be set to 0.4 in the resting scenario and 0.3 in the moving scenario; ω_ridge(t) = Persistence(t) represents the ridge persistence, which is used to measure the stability of the wavelet ridge frequency over time. The wavelet ridge frequency is the trajectory formed by tracking the peak along the time axis; ω_match(t) = Matching(t) represents the matching degree with the wavelet basis function (normalized correlation coefficient).

[0120] D3, and / or, calculate the energy proportion and instantaneous frequency change rate of the heart rate intrinsic mode function based on the above third time-frequency analysis results, and calculate the above third candidate heart rate weight based on the energy proportion of the above heart rate intrinsic mode function, the above instantaneous frequency change rate, and the preset third basic weight. The above heart rate intrinsic mode function is: the above intrinsic mode function whose main frequency is within the frequency range corresponding to the above heart rate.

[0121] In this embodiment of the application, assuming the weight of the third candidate heart rate is w_HHT(t), then w_HHT(t) can be: w_HHT(t) = α_HHT × ω_IMF(t) × ω_smooth(t).

[0122] Wherein, α_HHT represents the third basic weight, which can be set to 0.1 in the resting scenario and 0.4 in the exercise scenario; ω_IMF(t) = IMF_energy(t) represents the energy proportion of the intrinsic heart rate mode function; ω_smooth(t) = exp(-|df_HHT / dt| / σ_f) represents the instantaneous rate of change of frequency; σ_f represents the smoothing parameter, which can be 0.5Hz / s.

[0123] After calculating the weights of the first, second, and third candidate heart rates separately, normalization can be performed before calculating the fused heart rate. Assume the normalized weights are as follows: W_STFT(t) = w_STFT(t) / [w_STFT(t) + w_CWT(t) + w_HHT(t)]; W_CWT(t) = w_CWT(t) / [w_STFT(t) + w_CWT(t) + w_HHT(t)]; W_HHT(t) = w_HHT(t) / [w_STFT(t) + w_CWT(t) + w_HHT(t)]; The fused heart rate HR_fused(t) can then be: HR_fused(t)=W_STFT(t)×HR_STFT(t) + W_CWT(t)×HR_CWT(t) + W_HHT(t)×HR_HHT(t).

[0124] Among them, HR_STFT(t) is the first candidate heart rate, HR_CWT(t) is the second candidate heart rate, and HR_HHT(t) is the third candidate heart rate.

[0125] Optionally, in a resting scenario: the weight of the first candidate heart rate (W_STFT(t)) can be set to 0.5, the weight of the second candidate heart rate (W_CWT(t)) can be set to 0.4, and the weight of the third candidate heart rate (W_HHT(t)) can be set to 0.1. In an exercise scenario: the weight of the first candidate heart rate can be set to 0.3, the weight of the second candidate heart rate can be set to 0.3, and the weight of the third candidate heart rate can be set to 0.4.

[0126] In this embodiment of the application, after obtaining the fused heart rate, the confidence level of the fused heart rate can be evaluated. Assume the confidence level of the fused heart rate is Confidence(t): Confidence(t)=β1×Consistency(t) + β2×Agreement(t) + β3×Quality(t); Consistency(t) represents the consistency among the three methods.

[0127] Consistency(t) = 1 - [Std(HR_STFT, HR_CWT, HR_HHT) / Mean(HR_STFT, HR_CWT, HR_HHT)], where higher consistency (i.e., smaller standard deviation) indicates higher confidence.

[0128] Here, Agreement(t) represents the agreement between pairs of terms.

[0129] Agreement(t) = (1 / 3) × [I(|HR_STFT-HR_CWT|<5) + I(|HR_STFT-HR_HHT|<5) + I(|HR_CWT-HR_HHT|<5)]; The I(·) above is an indicator function. If the pairwise difference between the three methods is less than 5 bpm, then the consistency is 1. Where Quality(t) represents the signal quality.

[0130] Quality(t)=(SNR(t) - SNR_min) / (SNR_max - SNR_min); Normalized signal-to-noise ratio, ranging from [0,1]; The β1, β2, and β3 mentioned above represent weighting coefficients, which can be set as: β1=0.5, β2=0.3, β3=0.2.

[0131] After calculating the confidence level of the fused heart rate, the final fused heart rate can be determined by combining this confidence level.

[0132] Optionally, the final fused heart rate is obtained by time-domain filtering based on confidence level: If Confidence(t) < 0.5 (low confidence), then: HR_fused(t) = MedianFilter[HR_fused(t-2:t+2)] (replace with median filter); Then, Kalman filtering is used for smoothing: HR_smooth(t)=KalmanFilter[HR_fused(t)]; The HR_smooth(t) mentioned above is the final fused heart rate.

[0133] Optionally, the final fused heart rate is displayed, along with the corresponding confidence level, so that the user can determine the reliability of the final fused heart rate based on the confidence level.

[0134] It should be understood that the sequence number of each step in the above embodiments does not imply the order of execution. The execution order of each process should be determined by its function and internal logic, and should not constitute any limitation on the implementation process of the embodiments of this application.

[0135] Corresponding to the heart rate estimation method described in the above embodiments, Figure 2 A structural block diagram of a heart rate estimation device provided in an embodiment of this application is shown. For ease of explanation, only the parts related to the embodiment of this application are shown.

[0136] Reference Figure 2 The heart rate estimation device 2 is applied to an electronic device and includes: a video image sequence acquisition module 21, a time-series signal acquisition module 22, an independent component signal determination module 23, a time-frequency analysis module 24, and a heart rate determination module 25. Wherein: The video image sequence acquisition module 21 is used to acquire a continuous video image sequence of the face of the object under test.

[0137] The timing signal acquisition module 22 is used to determine the RGB timing signal that changes over time based on the above video image sequence.

[0138] The independent component signal determination module 23 is used to perform blind source separation on the above RGB time-series signal using an independent component algorithm to obtain independent component signals related to heart rate.

[0139] The time-frequency analysis module 24 is used to perform time-frequency analysis on the independent component signal related to heart rate using at least two time-frequency analysis methods to obtain the corresponding time-frequency analysis results. Each time-frequency analysis method is based on a different principle when performing time-frequency analysis.

[0140] The heart rate determination module 25 is used to estimate the heart rate of the subject under test based on the above time-frequency analysis results.

[0141] Since each time-frequency analysis method is based on a different principle, the errors in the time-frequency analysis results obtained after the analysis are independent of each other. This helps to reduce the error in the heart rate of the subject estimated based on each time-frequency analysis result, thereby improving the accuracy of the obtained heart rate.

[0142] Optionally, the above time-frequency analysis results include a first time-frequency analysis result, a second time-frequency analysis result, and a third time-frequency analysis result. The above-mentioned time-frequency analysis of independent component signals related to heart rate using at least two time-frequency analysis methods yields corresponding time-frequency analysis results, including: The short-time Fourier transform method was used to perform time-frequency analysis on the above-mentioned independent component signals related to heart rate, and the first time-frequency analysis result was obtained. The continuous wavelet transform method was used to perform time-frequency analysis on the above independent component signals related to heart rate, and the second time-frequency analysis result was obtained. The Hilbert-Huang transform method was used to perform time-frequency analysis on the above-mentioned independent component signals related to heart rate, and the third time-frequency analysis result was obtained.

[0143] Optionally, before performing time-frequency analysis on the heart rate-related independent component signal using the short-time Fourier transform method, the method further includes: Based on the physiological characteristics of heart rate and the requirements for time-frequency resolution and / or statistical stability, the window length range of the above short-time Fourier transform is determined. The window length for the short-time Fourier transform is determined based on the aforementioned window length range.

[0144] Optionally, determining the window length of the short-time Fourier transform based on the aforementioned window length range includes: Determine the current scene of the subject being tested. The scene includes: a resting scene, a moving scene, and a transition scene from resting to moving. Within the aforementioned window length range, the corresponding window length for the short-time Fourier transform is determined based on the aforementioned scenario.

[0145] Optionally, the heart rate estimation device 2 further includes: The wavelet basis function generation module is used to generate the wavelet basis function of the continuous wavelet transform based on the heart rate information before performing time-frequency analysis on the independent component signal related to heart rate using the continuous wavelet transform method.

[0146] Optionally, the wavelet basis function for generating the continuous wavelet transform based on heart rate information includes: The wavelet basis function of the continuous wavelet transform is generated based on at least one of the following information about the heart rate: frequency, morphology, the heart rate corresponding to the heart rate, and the systolic peak value of the heart rate corresponding to the heart rate.

[0147] Optionally, the above-mentioned continuous wavelet transform method is used to perform time-frequency analysis on the above-mentioned independent component signals related to heart rate to obtain the above-mentioned second time-frequency analysis result, including: The range of scale parameters is determined based on the frequency range corresponding to the heart rate mentioned above. For each of the above-mentioned scale parameters within the above-mentioned scale parameter range, wavelet coefficients are calculated based on the above-mentioned scale parameter, the preset time shift, the wavelet basis function of the above-mentioned continuous wavelet transform, and the above-mentioned independent component signal related to heart rate. The wavelet power spectrum is calculated based on each of the above wavelet coefficients, and the wavelet power spectrum is used as the result of the second time-frequency analysis.

[0148] Optionally, the Hilbert-Huang transform method includes empirical mode decomposition and Hilbert spectral analysis. The above-mentioned time-frequency analysis of the independent component signals related to heart rate using the Hilbert-Huang transform yields the third time-frequency analysis result, including: The above empirical mode decomposition adaptively decomposes the independent component signals related to heart rate into several intrinsic mode functions. The Hilbert spectral analysis is performed on each of the aforementioned intrinsic mode functions to extract the instantaneous frequency of each of the aforementioned intrinsic mode functions. The instantaneous frequency is used as the result of the aforementioned third time-frequency analysis.

[0149] Optionally, the heart rate estimation device 2 further includes: The main frequency calculation module is used to calculate the main frequency of each of the aforementioned intrinsic mode functions before performing the Hilbert spectral analysis on each of the aforementioned intrinsic mode functions. The heart rate intrinsic mode function determination module is used to filter out the intrinsic mode functions whose main frequency is within the frequency range corresponding to the heart rate based on the main frequency of each of the above intrinsic mode functions and the frequency range corresponding to the heart rate, so as to obtain the heart rate intrinsic mode functions. Correspondingly, the Hilbert spectral analysis performed on each of the aforementioned intrinsic mode functions includes: The Hilbert spectral analysis was performed on each of the aforementioned intrinsic heart rate modes.

[0150] Optionally, the above-mentioned Hilbert spectral analysis of each of the aforementioned intrinsic mode functions to extract the instantaneous frequencies of each of the aforementioned intrinsic mode functions includes: The Hilbert spectral analysis was performed on each of the aforementioned intrinsic mode functions, and outliers in the instantaneous frequencies of each of the aforementioned intrinsic mode functions extracted after the Hilbert spectral analysis were filtered out according to the frequency range corresponding to the aforementioned heart rate.

[0151] Optionally, after filtering out the outliers of the instantaneous frequencies of each of the intrinsic mode functions extracted after the Hilbert spectral analysis based on the frequency range corresponding to the heart rate, the method further includes: Estimate the normal values ​​corresponding to the abnormal values ​​of the above instantaneous frequencies; The abnormal values ​​of the instantaneous frequencies are replaced with the estimated normal values ​​to obtain the replaced instantaneous frequency sequence. The instantaneous frequency sequence after the above replacement is smoothed.

[0152] Optionally, estimating the heart rate of the subject based on the aforementioned time-frequency analysis results includes: Based on the results of the first time-frequency analysis, the candidate heart rate of the subject to be tested is estimated to obtain the first candidate heart rate; Based on the results of the second time-frequency analysis, the candidate heart rate of the subject under test is estimated to obtain the second candidate heart rate. Based on the above third time-frequency analysis results, the candidate heart rate of the above-mentioned subject is estimated to obtain the third candidate heart rate; Based on the first candidate heart rate, the second candidate heart rate, the third candidate heart rate, and the pre-set candidate heart rate weights for each candidate heart rate, the heart rate of the subject to be tested is estimated, wherein each of the above candidate heart rate weights is related to the current scene of the subject to be tested.

[0153] Optionally, the aforementioned candidate heart rate weights include a first candidate heart rate weight, a second candidate heart rate weight, and a third candidate heart rate weight. Before estimating the heart rate of the subject based on the first candidate heart rate, the second candidate heart rate, the third candidate heart rate, and the candidate heart rate weights pre-set for each candidate heart rate, the method further includes: The normalized peak energy is calculated based on the first time-frequency analysis results. The weight of the first candidate heart rate is calculated based on the normalized peak energy, the first candidate heart rate, and the preset first baseline weight. And / or, determine the wavelet ridge frequency based on the above second time-frequency analysis results, determine the ridge duration and matching degree based on the above wavelet ridge frequency, and calculate the above second candidate heart rate weight based on the above ridge duration, the above matching degree and the preset second basic weight, wherein the above matching degree is the matching degree between the wavelet coefficient at the ridge and the wavelet basis function of the above continuous wavelet transform generated based on the heart rate information; And / or, calculate the energy proportion and instantaneous frequency change rate of the heart rate intrinsic mode function based on the above third time-frequency analysis results, and calculate the above third candidate heart rate weight based on the energy proportion of the above heart rate intrinsic mode function, the above instantaneous frequency change rate and the preset third basic weight. The above heart rate intrinsic mode function is: the above intrinsic mode function with the main frequency in the frequency range corresponding to the above heart rate. The aforementioned first basic weight and / or the aforementioned second basic weight and / or the aforementioned third basic weight are related to the current scenario of the aforementioned test object.

[0154] It should be noted that the information interaction and execution process between the above-mentioned devices / units are based on the same concept as the method embodiments of this application. For details on their specific functions and technical effects, please refer to the method embodiments section, and they will not be repeated here.

[0155] Figure 3 This is a schematic diagram of the structure of an electronic device provided in an embodiment of this application. Figure 3As shown, the electronic device 3 of this embodiment includes: at least one processor 30 ( Figure 3 The diagram shows only one processor, a memory 31, and a computer program 32 stored in the memory 31 and executable on the at least one processor 30, wherein the processor 30 executes the computer program 32 to implement the steps in any of the above method embodiments.

[0156] The electronic device 3 can be a desktop computer, laptop, handheld computer, or cloud server, etc. This electronic device may include, but is not limited to, a processor 30 and a memory 31. Those skilled in the art will understand that... Figure 3 This is merely an example of electronic device 3 and does not constitute a limitation on electronic device 3. It may include more or fewer components than shown in the figure, or combine certain components, or different components. For example, it may also include input / output devices, network access devices, etc.

[0157] The processor 30 may be a Central Processing Unit (CPU), or it may be other general-purpose processors, digital signal processors (DSPs), application-specific integrated circuits (ASICs), field-programmable gate arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc. A general-purpose processor may be a microprocessor or any conventional processor.

[0158] In some embodiments, the memory 31 may be an internal storage unit of the electronic device 3, such as a hard disk or memory of the electronic device 3. In other embodiments, the memory 31 may be an external storage device of the electronic device 3, such as a plug-in hard disk, smart media card (SMC), secure digital (SD) card, flash card, etc., equipped on the electronic device 3. Furthermore, the memory 31 may include both internal and external storage units of the electronic device 3. The memory 31 is used to store the operating system, applications, bootloader, data, and other programs, such as the program code of the computer program. The memory 31 can also be used to temporarily store data that has been output or will be output.

[0159] Those skilled in the art will clearly understand that, for the sake of convenience and brevity, the above-described division of functional units and modules is merely an example. In practical applications, the above functions can be assigned to different functional units and modules as needed, that is, the internal structure of the device can be divided into different functional units or modules to complete all or part of the functions described above. The functional units and modules in the embodiments can be integrated into one processing unit, or each unit can exist physically separately, or two or more units can be integrated into one unit. The integrated unit can be implemented in hardware or as a software functional unit. Furthermore, the specific names of the functional units and modules are only for easy differentiation and are not intended to limit the scope of protection of this application. The specific working process of the units and modules in the above system can be referred to the corresponding process in the foregoing method embodiments, and will not be repeated here.

[0160] This application also provides a network device, which includes: at least one processor, a memory, and a computer program stored in the memory and executable on the at least one processor, wherein the processor executes the computer program to implement the steps in any of the above method embodiments.

[0161] This application also provides a computer-readable storage medium storing a computer program that, when executed by a processor, can implement the steps in the above-described method embodiments.

[0162] This application provides a computer program product that, when run on an electronic device, enables the electronic device to implement the steps described in the various method embodiments above.

[0163] If the integrated unit is implemented as a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, all or part of the processes in the methods of the above embodiments of this application can be implemented by a computer program instructing related hardware. The computer program can be stored in a computer-readable storage medium, and when executed by a processor, it can implement the steps of the various method embodiments described above. The computer program includes computer program code, which can be in the form of source code, object code, executable files, or certain intermediate forms. The computer-readable medium can include at least: any entity or device capable of carrying computer program code to a photographic device / electronic device, a recording medium, a computer memory, a read-only memory (ROM), a random access memory (RAM), an electrical carrier signal, a telecommunication signal, and a software distribution medium. Examples include USB flash drives, portable hard drives, magnetic disks, or optical disks. In some jurisdictions, according to legislation and patent practice, computer-readable media cannot be electrical carrier signals or telecommunication signals.

[0164] In the above embodiments, the descriptions of each embodiment have different focuses. For parts that are not described in detail or recorded in a certain embodiment, please refer to the relevant descriptions of other embodiments.

[0165] Those skilled in the art will recognize that the units and algorithm steps of the various examples described in conjunction with the embodiments disclosed herein can be implemented in electronic hardware, or a combination of computer software and electronic hardware. Whether these functions are implemented in hardware or software depends on the specific application and design constraints of the technical solution. Those skilled in the art can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of this application.

[0166] In the embodiments provided in this application, it should be understood that the disclosed apparatus / network devices and methods can be implemented in other ways. For example, the apparatus / network device embodiments described above are merely illustrative. For instance, the division of modules or units is only a logical functional division, and in actual implementation, there may be other division methods. For example, multiple units or components may be combined or integrated into another system, or some features may be ignored or not executed. Furthermore, the coupling or direct coupling or communication connection shown or discussed may be through some interfaces; the indirect coupling or communication connection between devices or units may be electrical, mechanical, or other forms.

[0167] The units described as separate components may or may not be physically separate. The components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the units can be selected to achieve the purpose of this embodiment according to actual needs.

[0168] It should be noted that the information collection process (such as the facial image collection process, fingerprint information collection process, etc.) / feature extraction process involved in this application is carried out with the user's knowledge and permission. That is, the information collection process / feature extraction process complies with the requirements of laws and regulations and does not constitute an act that harms the public interest.

Claims

1. A heart rate estimation method, characterized by, include: Acquire a continuous video image sequence of the face of the subject under test; The RGB timing signal that varies over time is determined based on the video image sequence; The independent component algorithm is used to perform blind source separation on the RGB time-series signal to obtain independent component signals related to heart rate; At least two time-frequency analysis methods are used to perform time-frequency analysis on the independent component signal related to heart rate to obtain corresponding time-frequency analysis results, wherein each of the time-frequency analysis methods is based on a different principle when performing time-frequency analysis; The heart rate of the subject under test is estimated based on the time-frequency analysis results.

2. The heart rate estimation method of claim 1, wherein, The time-frequency analysis results include a first time-frequency analysis result, a second time-frequency analysis result, and a third time-frequency analysis result. The step of performing time-frequency analysis on the independent component signal related to heart rate using at least two time-frequency analysis methods to obtain the corresponding time-frequency analysis results includes: The independent component signal related to heart rate is subjected to time-frequency analysis using the short-time Fourier transform method to obtain the first time-frequency analysis result; The independent component signal related to heart rate was subjected to time-frequency analysis using the continuous wavelet transform method to obtain the second time-frequency analysis result. The Hilbert-Huang transform method was used to perform time-frequency analysis on the independent component signal related to heart rate, and the third time-frequency analysis result was obtained.

3. The heart rate estimation method as described in claim 2, characterized in that, Before performing time-frequency analysis on the heart rate-related independent component signal using the short-time Fourier transform method, the method further includes: The window length range of the short-time Fourier transform is determined based on the requirements for time-frequency resolution and / or statistical stability of the physiological characteristics of heart rate. The window length of the short-time Fourier transform is determined based on the window length range.

4. The heart rate estimation method as described in claim 3, characterized in that, Determining the window length of the short-time Fourier transform based on the window length range includes: Determine the current scene of the object under test, including: a resting scene, a moving scene, and a transition scene from resting to moving; Within the window length range, the corresponding window length of the short-time Fourier transform is determined according to the scenario.

5. The heart rate estimation method of claim 2, wherein, Before performing time-frequency analysis on the heart rate-related independent component signal using the continuous wavelet transform method, the method further includes: The wavelet basis functions of the continuous wavelet transform are generated based on the heart rate information.

6. The heart rate estimation method of claim 5, wherein, The step of generating the wavelet basis function of the continuous wavelet transform based on heart rate information includes: The wavelet basis function of the continuous wavelet transform is generated based on at least one of the following information about the heart rate: frequency, morphology, the heart rate corresponding to the heart rate, and the systolic peak of the heart rate corresponding to the heart rate.

7. The heart rate estimation method of claim 5, wherein, The step of performing time-frequency analysis on the independent component signal related to heart rate using the continuous wavelet transform method to obtain the second time-frequency analysis result includes: The range of the scaling parameters is determined based on the frequency range corresponding to the heart rate. For each scale parameter within the scale parameter range, wavelet coefficients are calculated based on the scale parameter, a preset time shift, the wavelet basis function of the continuous wavelet transform, and the independent component signal related to heart rate. The wavelet power spectrum is calculated based on each of the wavelet coefficients, and the wavelet power spectrum is used as the result of the second time-frequency analysis.

8. The heart rate estimation method of claim 2, wherein, The Hilbert-Huang transform method includes empirical mode decomposition and Hilbert spectral analysis. The use of the Hilbert-Huang transform to perform time-frequency analysis on the independent component signals related to heart rate, yielding the third time-frequency analysis result, includes: The heart rate-related independent component signal is adaptively decomposed into several intrinsic mode functions through the empirical mode decomposition. The Hilbert spectral analysis is performed on each of the intrinsic mode functions to extract the instantaneous frequency of each of the intrinsic mode functions, and the instantaneous frequency is used as the result of the third time-frequency analysis.

9. The heart rate estimation method of claim 8, wherein, Before performing the Hilbert spectral analysis on each of the intrinsic mode functions, the method further includes: Calculate the dominant frequency of each of the intrinsic mode functions; Based on the dominant frequency of each intrinsic mode function and the frequency range corresponding to the heart rate, the intrinsic mode functions with the dominant frequency within the frequency range corresponding to the heart rate are selected to obtain the heart rate intrinsic mode function; The Hilbert spectral analysis of each of the intrinsic mode functions includes: The Hilbert spectral analysis was performed on each of the aforementioned heart rate intrinsic mode functions.

10. The heart rate estimation method of claim 8, wherein, The step of performing Hilbert spectral analysis on each of the intrinsic mode functions to extract the instantaneous frequencies of each of the intrinsic mode functions includes: The Hilbert spectrum analysis is performed on each of the intrinsic mode functions, and outliers in the instantaneous frequencies of each of the intrinsic mode functions extracted after the Hilbert spectrum analysis are filtered out according to the frequency range corresponding to the heart rate.

11. The heart rate estimation method of claim 10, wherein, After filtering out outliers in the instantaneous frequencies of each intrinsic mode function extracted after the Hilbert spectral analysis based on the frequency range corresponding to the heart rate, the method further includes: The estimated normal value corresponding to the abnormal value of the instantaneous frequency; The abnormal values ​​of the instantaneous frequencies are replaced with the estimated normal values ​​to obtain the replaced instantaneous frequency sequence; The replaced instantaneous frequency sequence is then smoothed.

12. Heart rate estimation method according to any of claims 2 to 11, characterized in that, The step of estimating the heart rate of the subject based on the various time-frequency analysis results includes: Based on the first time-frequency analysis result, the candidate heart rate of the test subject is estimated to obtain the first candidate heart rate; Based on the second time-frequency analysis results, the candidate heart rate of the test subject is estimated to obtain the second candidate heart rate; Based on the third time-frequency analysis results, the candidate heart rate of the test subject is estimated to obtain the third candidate heart rate; The heart rate of the test subject is estimated based on the first candidate heart rate, the second candidate heart rate, the third candidate heart rate, and the candidate heart rate weights pre-set for each candidate heart rate, wherein each candidate heart rate weight is related to the current scene of the test subject.

13. The heart rate estimation method of claim 12, wherein, The candidate heart rate weights include a first candidate heart rate weight, a second candidate heart rate weight, and a third candidate heart rate weight. Before estimating the heart rate of the subject based on the first candidate heart rate, the second candidate heart rate, the third candidate heart rate, and the candidate heart rate weights pre-set for each candidate heart rate, the method further includes: The normalized peak energy is calculated based on the first time-frequency analysis result, and the weight of the first candidate heart rate is calculated based on the normalized peak energy, the first candidate heart rate, and the preset first base weight. And / or, determine the wavelet ridge frequency based on the second time-frequency analysis result, determine the ridge duration and matching degree based on the wavelet ridge frequency, and calculate the second candidate heart rate weight based on the ridge duration, the matching degree and the preset second basic weight, wherein the matching degree is the matching degree between the wavelet coefficient at the ridge and the wavelet basis function of the continuous wavelet transform generated based on the heart rate information; And / or, calculate the energy proportion and instantaneous frequency change rate of the heart rate intrinsic mode function based on the third time-frequency analysis result, and calculate the third candidate heart rate weight based on the energy proportion of the heart rate intrinsic mode function, the instantaneous frequency change rate, and the preset third basic weight. The heart rate intrinsic mode function is: the intrinsic mode function with the main frequency in the frequency range corresponding to the heart rate. The first basic weight and / or the second basic weight and / or the third basic weight are related to the current scene of the object under test.

14. A heart rate estimation apparatus, characterized by include: The video image sequence acquisition module is used to acquire a continuous video image sequence of the face of the object under test; A timing signal acquisition module is used to determine the RGB timing signal that changes over time based on the video image sequence; The independent component signal determination module is used to perform blind source separation on the RGB time-series signal using an independent component algorithm to obtain independent component signals related to heart rate. The time-frequency analysis module is used to perform time-frequency analysis on the independent component signal related to heart rate using at least two time-frequency analysis methods to obtain the corresponding time-frequency analysis results, wherein each of the time-frequency analysis methods is based on a different principle when performing the time-frequency analysis; The heart rate determination module is used to estimate the heart rate of the subject under test based on the various time-frequency analysis results.

15. An electronic device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the computer program, it implements the method as described in any one of claims 1 to 13.

16. A computer-readable storage medium, the computer-readable storage medium storing a computer program, characterized in that, When the computer program is executed by a processor, it implements the method as described in any one of claims 1 to 13.

17. A computer program product, characterised in that, Includes a computer program, which, when executed, causes the electronic device to perform the method according to any one of claims 1 to 13.

Images