A Trigonometric Projection Method for Overcoming High Intensity Heart Rate Caused Waveform Deformation in Electrocardiogram Biometrics

Although electrocardiogram (ECG) has been proven as a biometric for human identification, applying biometric technology remains challenging with diverse heart rate circumstances in which high intensity heart rate caused waveform deformation may not be known in advance when ECG templates are registered. A calibration method that calculates the ratio of the length of an unidentified electrocardiogram signal to the length of an electrocardiogram template is proposed in this paper. Next, the R peak is used as an axis anchor point of a trigonometric projection (TP) to attain the displacement value. Finally, the unidentified ECG signal is calibrated according to the generated trigonometric value, which corresponds to the trigonometric projection degree of the ratio and the attained displacement measurement. The results reveal that the proposed method provides superior overall performance compared with that of the conventional downsampling method, based on the percentage root mean square difference (PRD), correlation coefficients, and mean square error (MSE). The curve fitting equation directly maps from the heart rate levels to the TP degree without prior registration information. The proposed ECG calibration method offers a more robust system against heart rate interference when conducting ECG identification.


Introduction
Electrocardiogram (ECG), an electrical signal from the heart, has been used for medical diagnosis purposes for more than a century. Many studies have demonstrated the use of one-lead, resting ECG signals to identify people [1][2][3]. ECG characteristics are usually used in biometric systems, such as ECG waveform morphology [1,4,5], amplitude and interval features from fiducial points [5][6][7][8][9], frequency components [10,11], eigenvalues [12], mathematical models [13], and chaos [14]. However, heart-rhythm-related variations cause instability in ECG biometric systems, especially with highintensity heart rate (HR). HR normalization is crucial for conducting accurate diagnosis, even in clinical applications. Bazett's and Framingham's formulas are commonly used for normalizing QT time intervals. These equations adjust the QT intervals by dividing them by the square and cube roots of the RR intervals. However, Bazett's formula has been criticized for being inaccurate for measuring rapid heart rates (HRs). Moreover, Framingham's formula may not be accurate for low HRs [15]. Because of variations in human ECG signals, normalizing the factors of an individual is challenging if no prior data are collected for certain individuals.
Therefore, here we report the state of the art of ECG. Several studies on ECG identification have been conducted on low-and moderate-intensity HR changes, such as in the time domain [5][6][7][8][9], frequency domain [10,11], time-frequency domain [16][17][18][19], and chaos-related projections [14]. The HR variation range is 10-30 beats/min increased from resting HR. Awal et al. [13] used mathematical model for fitting the ECG model using the nonlinear least square technique, which can potentially simulate different HR conditions. Khalil et al. [8] found that the QRS complex of an ECG signal has the most unique signature, and applied high-order Legendre polynomials to HR without the use of any ECG calibration methods. Chen et al. [20] applied a complexity-based approach to handle abnormal ECG signals while conducting biometric identification. Irvine et al. [12] designed a robust eigenPulse method to downsample (DS) and normalize ECG features through principal component analysis (PCA) to enhance the overall identification performance. Wang et al. [10] proposed a model that combines autocorrelation and discrete cosine transform to overcome HR using frequency-domain-related features. However, most of the aforementioned methods are suitable for low-and moderate-intensity HR applications or for those in which the HR range is unspecified.
Therefore, the aim of this study was to provide a timedomain ECG calibration method that overcomes HR interference in ECG morphology through comparison with the traditional downsampling method for human identification. We hypothesized that a calibration method only stores healthy resting ECG signal morphologies as prior information (without any high-intensity HR ECG waveform in gallery), which enables comparison possible among ECG waveforms at different HR levels. Therefore, a general formula is provided for normalizing ECG waveforms on the basis of HR or RR intervals.

Experimental Setup
Instead of healthy resting or arrhythmia [21] ECG databases, various HR ECG signals obtained from healthy individuals on bicycles with similar demographics were employed under the regulations of Tzu Chi Hospital, Taiwan (IRB099-125). A total of 50 randomly selected volunteers with similar demographic information, including 33 men and 17 women with no heart disease were included in this study. Their demographic information was similar for age (23.22 ± 2.46 years), height (168.3 ± 8.61 cm), weight (63.94 ± 12.22 kg), and body mass index (BMI; 22.5 ± 3.57 kg/m 2 ). For the experiment depicted in Fig. 1, the volunteers were asked to rest for 2 min before using an exercise bike (818E Ergomedic Fitness Testing Bike, Monark) and to bike for at least 3 min until their heart rate reached at least 130 beats per minute (bpm) or 60 bps above their resting heart rate (rHR). The resistance or tension of the exercise bike was controlled at 2 kg-force (kgf). Once the required HR was achieved, the volunteers were asked to remain in their seats for at least 10 min until their HR stabilized to their rHR, which was computed as one over their RR interval. Each session was monitored using a data acquisition system (MP35, BioPack, Inc.) at 500 samples per second (sps) to obtain one-lead ECG signals. ECG electrodes were placed on the inner sides of participants' front arms to avoid interference with leg movements while riding the bike. Instead of the standard ECG recording on the right leg, the ground lead was placed on the inner left arm, which creates more interference challenges than the chest ECG. Hence, the signal bandwidth reduction is a trade-off filter design to handle such noise. In addition, the inner arm placement causes slightly lower amplitude of the entire ECG from the standard lead-one ECG measurement because the arm distance can be regarded as resistance in an equivalent circuit [1,5].
High-intensity interval exercising comprises short workout intervals (70% to 90% maximum HR) and low-intensity intervals (55% to 65% maximum HR) [22]. The filtered ECG was segmented heartbeat-by-heartbeat and categorized into either gallery G (identified data) or probe P (unidentified data) sets. No overlap was observed between the randomly selected heartbeats from the probe and gallery data sets. The  Fig. 1 a HR ranges of 50 individuals in the probe set at six different HR levels (300 randomly selected beats); b Experimental setup gallery set included a randomly selected resting heartbeat from each individual. This set was used as a template database to form a system development data set. The probe set included an unidentified ECG set as the input signals of the proposed method, which comprised six subsets-p0, p1, p2, p3, p4, and p5-representing six levels of rhythm strengths without including the rHR. This set was used to form an unidentified ECG data set. The six levels of the probe set p n , where n = 0, … , 5 , included randomly selected ECG signals within HR ranges from [(rHR − 5) + 10 × (n + 1)] to [(rHR + 5) + 10 × (n + 1)] bpm, where n = 0, … , 5 . Here, rHR ± 5 describes the selected ECG beats obtained by varying the rHR (± 5). For example, p0 contained the HR range of all randomly selected beats from (rHR -5) + 10 to (rHR + 5) + 10 bpm. The high intensity HR caused waveform deformation groups were categorized as p4 and p5; the moderate-intensity, HR-caused waveform deformation groups were categorized as p2 and p3; the low-intensity, HR-caused waveform deformation groups were categorized as p0 and p1. In total, the probe set including six subsets had 300 (six levels × 50 people) randomly selected beats. The overall HR ranges of each subset of P are listed as follows:

Preprocessing and Heartbeat Segmentation
Many ECG biometric characteristics are not for medical diagnosis purposes. Tan's biometric system in mobile devices showed a total 99.52% subject verification accuracy. Their bandpass filter was designed using 50 Hz as the high cut-off frequency to retain as much of the ECG signal energy as possible while removing the power line (60 Hz) and other high frequency interferences [23,24]. Other ECG biometric studies applied bandwidths of 1-40 Hz [25,26] and 0.5-40 Hz [27,28]. We also chose a low pass filter with a cut-off of 50 Hz after preprocessing [1,23,24]. The baseline wander interferences were eliminated using a median filter. The median filter includes both large-scale (600 points) and small-scale (200 points) windows to extract baseline wander waveforms [29,30]. A detailed definition of segmenting the waveform of an ECG signal is provided in a previous study [31].

Trigonometric Projection Method
The trigonometric projection (TP) method was developed to provide an alternative downsampling process for the maintenance of more biometric details in some ECG intervals with a trade-off of unimportant sample points because not all intervals among the ECG fiducial points are equally important in ECG biometrics [8]. The proposed TP calibration method projects ECG heartbeats at all HR levels to the same sampling length by adjusting the projection angle.
In general, resting ECG heartbeats have the longest sample sequences and the sample points reduce continuously when the HR level increases. In the TP method, various ECG lengths from probe and gallery data sets are projected to target distance length L targ . The downsampling target distance length L targ is set as a constant that is generally smaller than the original distance length of the fastest HR heartbeat because of its limitation on ECG cut-off frequency (50 Hz), aliasing (< 100 Hz), and HR limitation (200 points for a 150 bpm ECG heartbeat). In this study, the target distance lengths L targ were set to 50, 100, 150, and 170 sample points for evaluation. The original distance length L orig was defined as the length of an ECG heartbeat in the probe and gallery data sets.
As shown in Fig. 2, the calibration is initiated by the receipt of the uncalibrated ECG signal as an input. The coordinates of the unidentified ECG signals u with various numbers of sample points are multiplied by a generated trigonometric value that corresponds to a TP degree used for projecting ECG signal u to ECG signal v at a determined length on a plane. An ECG gallery signal is first converted to attain a calibrated ECG template with the same length. The detail steps are described below.
First, all ECG templates x g with original distance length L orig (x g ) obtained from the gallery data set are projected to the target distance lengths L targ . Then, an unidentified ECG signal x p (obtained from the probe set) is received to calculate the ratio of the original distance length L orig (x p ) of the unidentified ECG signal x p to the target distance length L targ for a comparison with the calibrated ECG template x g obtained from the gallery data set. Hence, the ratio Eq. (1) is as follows: where the distance length L indicates the number of sample points on a segmented waveform, such as L targ (x) and L orig . This ratio is less than 1 in general.
In Fig. 2, L orig , located on the S1 plane, is generally greater than or equal to L targ , located on the S2 plane. The TP degree θ was calculated on the basis of the cosine of the arc of the ratio presented in Eq. (2): where θ represents the ratio degree of the distance length of an uncalibrated ECG signal to the new target distance length of the ECG heartbeats. The calibrated ECG template x g is presented on the x-axis (unit: number of points), and its R peak (R g ) is set as the anchor point. The R-peak (R p ) x-axis section of the unidentified ECG signal x p is aligned on the x-axis portion of R g . The axis anchor point attains a displacement measurement to shift the coordinates of all ECG sample points for further TP, and extra samples outside the certain fix length are truncated. The generated projection angle θ and its cosine values are applied for projecting the distance length of the uncalibrated ECG signal on a plane S2 with the required width. Finally, the unidentified ECG is calibrated according to the calculated trigonometric angle, which corresponds to the TP degree and the attained displacement measurement. However, the amplitudes of the ECG signals remained unchanged in this study. The projection formula of the coordinate positions in the x-axis can be summarized in Eq. (3) as: Based on point-to-point transform, the new x-axis coordinate is equal to the distance from the processing sample point to R peak ( L toR ) times cos plus the original x-axis coordinate ( L R ) of an R peak and one edge point, where the function fix() represents rounding toward zero and the R peak location is set as the new origin (0, 0). L toR is positive when the processing points are on the right-hand side of the R point; L toR is negative when the processing points are on the left-hand side of the R point.
In digital signal processing, downsampling, decimation, and resampling are the most common methods to reduce samples in the time domain to form a fixed-length signal [32]. Therefore, the proposed TP method was compared with a conventional downsampling method represented as x(Mn) and a resampling [33] method. The downsampling factor (commonly denoted by M) is usually an integer or a rational (3) x new = fix L toR * cos + L R + 1.
fraction (as previous mentioned, y/x). An integer factor to reduce samples is not always found for ECG decimation. Hence, downsampling with a rational factor decimation was applied for comparison in this study. Unlike downsampling, the resampling method samples the original data sequence based on the interpolation points, which may not have amplitude values overlapping between the original and resampled points. The filtering effect on resampling method risks that the identification features related to fictional ECG points may change due to the algorithm passing the sharp fictional points, for example, the R points. For comparison, TP method maintains the R point as the anchor points without any change in R amplitude.
In this study, the projection or downsampling process was not used to reduce the data rate or the amount of the data but to reduce the morphological differences due to the intensity of the HRs. Finally, all six levels of the ECG heartbeats were compared with the resting template of the same person at a preset length according to the PRD, correlation coefficients, and MSE [34,35].

Results
The advantage of the proposed method is that various lengths of ECG signals can be normalized without requiring any preregistrations. Figure 3 illustrates the PRD (%) values of the proposed method and the downsampling method, which are 44.37% and 45.32%, respectively. The lower the PRD values, the higher the similarity between two waveforms. This result in Fig. 3 indicates that the TP method yielded not only lower PRD values but also less distortion in the waveforms.  Table 1 shows the comparison results among high-, moderate-, and low-intensity HR subsets using three methods: PRD, correlation coefficients, and MSE. The TP method produced the best performance in terms of correlation coefficients, MSE, and PRD for the low-intensity HR subsets. The result showed that the resampling filter effect increased the performance in terms of PRD, correlation coefficients, and MSE compared with pure downsampling. However, this is an unfair comparison because the TP method retains the original values of samples with no antialiasing filter effect, which is the same as pure downsampling. Hence, we focused on a comparison between the TP and DS methods.
The overall PRD results indicated that PRD increased when HR increased for all target lengths. The overall PRD of the proposed method was lower than that of the down-sampling method. A higher number of sample points provided better PRD performance on average. A detailed comparison is provided in Fig. 4a.
The PRD values of high-intensity HR subsets were more than twice as high as those of the low-intensity HR subsets. The TP method demonstrated better performance than the downsampling method. Nevertheless, the difference in PRD values decreased when the number of sample points increased.
The overall correlation results indicated that the correlation coefficients decreased when HR increased for all target lengths. The overall correlation coefficients of the proposed method were higher than those of the downsampling method. A higher number of sample points produced better correlation coefficients on average. A detailed comparison is displayed in Fig. 4b.  Table 1 Comparison among high-, moderate-and low-intensity HR subsets by using three methods in terms of PRD (%), correlation coefficients and MSE (10 −3 mV 2 ) The bold font indicates the best performance among three methods  The average correlation coefficients are listed in Table 1 for the high-, moderate-, and low-intensity HR subsets. A strong correlation between the probe and template signals was observed even for the high-intensity HR subsets. The correlation coefficients obtained using the trigonometric method were greater than 0.88, which is the suggested prescreening indicator before further feature classification in the ECG identification process [5]. Unlike the PRD, the correlation coefficients still maintained high correlation levels in all subsets. The trigonometric method exhibited better performance than the downsampling method. Comparisons of the TP and DS methods are presented in terms of MSE in Fig. 4c, which shows that the downsampling method produced a larger MSE when the target length was small. However, the TP method remained relatively stable at all four target lengths. Finally, the TP method exhibited smaller errors or distances than the traditional downsampling and resample methods for those parameters under all conditions. However, the differences between the parameters of the two methods decreased as the target length L targ increased.
First rank accuracy is often used to evaluate identification systems. In Fig. 5, the first rank accuracy is used to compare ECG identification performance on multiple time-domain systems using the TP or downsampling methods. The red and blue lines in Fig. 5 illustrate that the norm distances of the features remain stable after the TP method. If the TP method is used as a prescreening process and combined with the norm distances of 17 feature vectors (RQ amplitude, QS duration, RS amplitude, ST amplitude, QT duration, RS slope, QRS triangular area, RS amplitude/TS amplitude, RS 2 amplitude, PQ amplitude, QS amplitude, RP amplitude, RT amplitude, ST slope Angle Q, R, and S) [1], then the overall performance can be further improved. As shown by the green line, due to template matching being based on correlation coefficients, the TP method more improved template matching in terms of average accuracy from p0 to p5 than the downsampling and no-resampling methods (the brown and yellow lines, respectively) with 8.3% and 1.3% increased accuracy, respectively. For the pink line, QRS matching shows the best performance at HR level 5 and 6 without any resampling process, which confirms previous study results [33]. However, QRS matching was not well performed at HR levels 1-3. Overall, the TP calibration method can be applied to the ECG waveform to reduce differences between the six HR levels regardless of the methods subsequently used for identification.
TP degrees were obtained using the ratio of the original ECG length to the target ECG length. Our database contains 300 ECG heartbeats with various HRs, which was used to find the optimal fit in Fig. 6 for mapping the TP degree to the mean HR level. Our results revealed that the nonlinear cubic curve fitting with the norm of the residuals is 0.58 in Eq. (4): where r denotes the HR level and represents the TP degree. The norm of the residuals is a measure of goodness of fit, where a smaller value indicates a better fit.
This formula and the obtained residual plots revealed the goodness of fit among the p0, p2, p3, and p5 subsets. This information may aid in the development of a real-time ECG biometric system through combination with Eq. (4) without any HR-related preregistration as pre-knowledge. Only RR intervals can potentially determine the TP degree without the use of whole ECG beats.

Discussion
We proposed a TP method as an alternative to the downsampling method to overcome ECG waveform deformation caused by high-intensity HR when an ECG biometric system matches multiple ECG waveforms with different HRs. The method successfully reduces the overall PRD results compared with no process or the downsampling method. We suggest projecting one ECG heartbeat to 150 points to increase the overall accuracy of ECG identification and reduce the waveform deformation effect.
(4) θ = 0.13r 3 − 1.4r 2 + 2.1r + 53 Unlike the traditional downsampling method, which truncates samples with a fixed distance, the TP method nonlinearly truncates samples, which produces the best overall performance for ECG identification in the time domain compared with downsampling and resampling. This has less influence on the HR-related intervals such as the PQ and ST intervals. However, this method employs R peaks as the center for truncation, leading to more distortion in the QRS complex and less distortion in other areas. However, compared with the traditional downsampling method, the main limitation of the TP method is that the QRS complexes may exhibit greater distortion than that in other areas of the ECG waveforms.
The above results explain what happens at low HRs when the PRD and correlation coefficients are more similar (no alterations in the ECG due to exercise or higher HRs), and remarkably improved performance can be achieved when the proposed method is used. However, this performance is not observed in HR levels 4 and 6. The equations of the PRD and correlation coefficients reveal whether distortion occurs in the larger amplitude area, namely the QRS complex. The denominators at these points increase, and the PRD and correlation coefficients decrease. For comparison, the denominators of the MSE are not related to signal amplitudes, and are insignificant at lower resolutions.
Technically, it is hard for a biometric system to collect ECG templates with a full range of HRs. Equation (4) provides an approximate model to predict the projection angles from resting ECG to nonresting ECG, so the different heart rhythm ECG waveforms can be normalized for comparison in the biometric database. We suggest that the sampling frequency of the ECG signals should be increased before applying the TP method. When the formula is applied to the ECG biometric system, our algorithm may provide an appropriate calibration procedure for increasing the screening capability in ECG biometric systems.

Conclusions
The proposed calibration method is an effective normalization method for ECG biometrics, especially for use in cases of high intensity HR caused waveform deformation. The calibration was conducted on an unidentified electrocardiogram signal based on the generated trigonometric value, which corresponded with the TP degree of the ratio and the attained displacement measurement. The results revealed that the PRD and correlation coefficients were within the ranges of 22.60-51.49% and 0.8921-0.9737, respectively. The results also indicated that the TP method is favorable in terms of the MSE. The curve fitting equation directly maps the HR levels to the TP degree and contributes to the future development of real-time ECG biometric systems without prior registration information. We presented an HR calibration method for eliminating high intensity HR caused waveform deformation when conducting ECG identification under various HR conditions.