Preprocessing of HRV data

Accurate preprocessing of HRV data is crucial before any HRV analysis. Artifacts in RR interval time series, whether technical or physiological, can significantly disrupt signal analysis. Technical artifacts may arise from issues like noise or algorithm inaccuracies, manifesting as missing or misaligned beats. Physiological artifacts include ectopic beats and arrhythmias. To ensure data integrity, manually check ECG recordings for artifacts, using only clean sections for analysis ^[10]. If artifact-free data is insufficient, interpolation methods may mitigate these effects ^{[2], [3], [6]}.

Furthermore, slow changes in mean heart rate during recordings can adversely affect some HRV analysis parameters. These nonstationarities, common in HRV signals, should be removed before analysis. Detrending can be achieved using polynomial models ^[5], [7] or the advanced smoothness priors method included in Kubios HRV analysis software ^[9].

Noise detection in Kubios HRV software

Kubios HRV software automatically detects noise using raw ECG or interbeat interval data. This feature is crucial for analyzing long-term ambulatory recordings where signal quality varies. To ensure the reliability of HRV analysis results, periods of noisy data are excluded so that HRV parameters remain unaffected by noise.

An example of automatic noise detection for ECG data is illustrated in the figure below. Please note that the purpose of noise detection is to identify segments where noise distorts several consecutive beat detections, potentially affecting HRV analysis accuracy. It’s important to understand that individual intermittent abnormal beat intervals, such as ectopic beats, are not typically marked as noise since they can be accurately corrected by the automatic beat correction described next.

Figure 1: Illustration of automatic noise detection (screenshot taken from Kubios HRV Scientific software).

Beat correction in Kubios HRV software

Kubios HRV data analysis software offers robust tools for HRV preprocessing. Initially, beat detections can be manually corrected when raw ECG (or PPG) data is available. This manual correction ensures the accuracy of beat intervals but is feasible only with a limited number of misdetections. Additionally, the software includes two algorithms for beat correction: 1) Automatic beat correction, and 2) Threshold-based beat correction. Below, we detail these methods.

Automatic beat correction algorithm

The automatic beat correction algorithm, detailed in [4], detects artifacts from the dRR series—a time series of differences between successive RR intervals. This series helps robustly differentiate ectopic and misplaced beats from normal sinus rhythm using a time-varying threshold (Th). The threshold is adapted to varying HRV levels by estimating it from the time-varying distribution of dRR values. A quartile deviation of the 90 surrounding beats is calculated and multiplied by a factor of 5.2 to cover 99.95% of beats, assuming normal distribution. However, since RR intervals often deviate from normality, some normal beats may exceed this threshold, necessitating a decision algorithm to identify artifact beats.

Table 1: Performance of the automatic artefact correction algorithm in detecting normal sinus rhythm beats, and simulated missed, extra and misaligned (with three different levels of displacement) beats.

Kubios HRV beat correction algorithm accuracy

Ectopic beats typically show negative-positive-negative (NPN) or positive-negative-positive (PNP) patterns in the dRR series, while long beats show positive-negative (PN) and short beats show negative-positive (NP) patterns. Only dRR segments with these patterns are classified as artifacts. Missed or extra beats are detected by comparing the current RR value with the median of the surrounding 10 RR intervals (medRR). A beat is flagged as missed if the current RR interval, RR(i), satisfies the condition

$\begin{equation*} \left|\frac{RR(i)}{2}-medRR(i)\right|<2 Th \end{equation*}$

and a beat is flagged as extra if two successive RR intervals, RR(i) and RR(i+1), meet the criteria

$\begin{equation*} \left|RR(i)+RR(i+1)-medRR(i)\right|<2 Th. \end{equation*}$

Correction of detected artefacts: Detected ectopic beats are corrected by replacing the corrupted RR times with interpolated values. Similarly, beats that are too long or too short are adjusted by interpolating new values into the RR time series. Missed beats are corrected by inserting a new R-wave occurrence time, and extra beats are addressed by removing the superfluous R-wave detection and recalculating the RR interval series.

Threshold based beat correction algorithm

The threshold-based beat correction algorithm evaluates each RR interval against a median-filtered local average, which excludes outliers from affecting the average. If an RR interval deviates from this local average by more than a specified threshold, it is marked for correction. Threshold levels range from very low to very strong , and can be custom set:

Very low: 0.45 sec
Low: 0.35 sec
Medium: 0.25 sec
Strong: 0.15 sec
Very strong: 0.05 sec
Custom threshold in seconds

For instance, the Medium level (0.25 sec) targets RR intervals significantly deviating from the local average. These are corrected using cubic spline interpolation. It’s important to note that these thresholds are adjusted based on the mean heart rate; higher heart rates have smaller thresholds due to expected decreased variability, and vice versa.

Because artifacts are identified through simple thresholding and heart rate variability levels are highly individual, the correction level must be tailored to each case. Initially, determine if any beat intervals require correction. If so, select the lowest possible correction level that effectively corrects abnormal beats without causing overcorrection.

HRV detrending in Kubios HRV software

To avoid the effect of slow nonstationary trends (slow changes in mean HR basically) over HRV analysis results, such trends by default are removed in Kubios HRV software using a smoothness priors method. The theory behind the smoothness priors detrending method is described shortly here. For more details, see [9].

Let $z\in\R^N$ denote the RR interval time series which can be considered to consist of two components

$\begin{equation*} z = z_{\rm{stat}}+z_{\rm{trend}} \end{equation*}$

where $z_{\rm{stat}}$ is the nearly stationary RR interval series of interest, $z_{\rm{trend}}$ is the low frequency aperiodic trend component, and $N$ is the number of RR intervals. Suppose that the trend component can be modeled with a linear observation model as

$\begin{equation*} z_{\rm{trend}} = H\theta + e \end{equation*}$

where $H\in\R^{N\times p}$ is the observation matrix, $\theta\in\R^p$ are the regression parameters, and $e$ is the observation error. The task is then to estimate the parameters by some fitting procedure so that $\hat{z}_{\rm{trend}} = H\hat{\theta}$ can be used as the estimate of the trend. The properties of the estimate depend strongly on the properties of the basis vectors (columns of the matrix $H$ ) in the fitting. A widely used method for the solution of the estimate $\hat{\theta}$ is the least squares method. However, a more general approach for the estimation of $\hat{\theta}$ is used here. That is, the so-called regularized least squares solution

$\begin{equation*} \hat{\theta}_\lambda = \arg \min_\theta \left\{ \Vert z-H\theta \Vert^2 + \lambda^2\Vert D_2 (H\theta)\Vert^2 \right\} \end{equation*}$

where $\lambda$ is the regularization parameter and $D_2$ indicates the discrete approximation of 2nd order derivative. This is a modification of the ordinary least squares solution, dragging the solution to the direction in which the norm $\Vert D_2(H\theta)\Vert$ gets smaller. In this way, prior information about the predicted trend $H\theta$ can be implemented to the estimation. The solution of previous equation can be written in the form

$\begin{equation*} \hat{\theta}_\lambda = (H^TH + \lambda^2H^TD_2^TD_2H)^{-1} H^T z \end{equation*}$

and the estimate for the trend which is to be removed can be written

$\begin{equation*} \hat{z}_{\rm{trend}} = H\hat{\theta}_\lambda. \end{equation*}$

The selection of the observation matrix $H$ can be implemented according to some known properties of the data $z$ . However, here we use a trivial choice of identity matrix $H = I\in\R^{N\times N}$ . In this case, the regularization part of the minimized function simply draws the solution towards a curve for which the 2nd order difference is zero, i.e. towards first order linear curve. With these specific choices, the detrended nearly stationary RR series can be written as

$\begin{equation*} \hat{z}_{\rm{stat}} = z- H\hat{\theta}_\lambda = (I-(I+\lambda^2D_2^TD_2)^{-1}) z. \end{equation*}$

This detrending method works like a time-varying highpass filter, where the cutoff frequency can be adjusted by changing the smoothing parameter $\lambda$ . In Kubios HRV Premium software, the cutoff frequency related to the given smoothing parameter is presented. The smoothing parameter is by default selected to remove only VLF frequency components from the RR time series.

Frequently Asked Question (FAQ)

What are the common sources of noise in HRV data?

The most prevalent source of noise and artifacts in HRV analysis is motion artifacts, which distort ECG or PPG waveform data, leading to incorrect beat detections and abnormal interbeat interval values. Other significant sources include poor electrode-skin contact in ECG recordings and electrical activity from nearby muscles. In PPG recordings, changes in peripheral circulation or variations in ambient lighting can also introduce noise.

What role does beat correction play in HRV analysis?

Accurate HRV analysis requires interbeat interval data to reflect normal sinus rhythm and be free from artifacts. Even a single abnormal beat interval or ectopic beat, if not corrected, can significantly distort HRV metrics. This distortion is particularly pronounced in parameters that reflect beat-to-beat variability, such as RMSSD or high frequency (HF) power, which are highly sensitive to such artifacts.

How do you analyze HRV data that is corrupted by noise or artifacts?

When analyzing HRV data affected by noise or artifacts, you have two options. First, you can analyze only the periods of data that are free from noise and artifacts. If this is not feasible, you must preprocess the HRV data before analysis by excluding noisy periods and applying beat correction to address any remaining intermittent artifacts and abnormalities. Kubios HRV software products provide excellent tools to accomplish this effectively.

How do nonstationarities affect HRV analysis, and why should they be removed?

Nonstationarities, such as slow changes in heart rate that fall into the very low frequency range (below 0.04 Hz), are common in HRV signals. These changes are typically not directly related to the autonomic nervous system’s regulation of heart function. Thus, for accurate short-term HRV analysis assessing ANS function and the balance between the sympathetic and parasympathetic branches, removing these nonstationarities is crucial. They primarily affect HRV parameters that reflect overall variability, such as SDNN or total spectral power.

References

G.G. Berntson, J.T. Bigger Jr., D.L. Eckberg, P. Grossman, P.G. Kaufmann, M. Malik, H.N. Nagaraja, S.W. Porges, J.P. Saul, P.H. Stone, and M.W. Van Der Molen. Heart rate variability: Origins, methods, and interpretive caveats. Psychophysiol, 34:623–648, 1997.
N. Lippman, K.M. Stein, and B.B. Lerman. Nonlinear predictive interpolation: a new method for the correction of ectopic beats for heart rate variability analysis. J Electrocardiol, 26:S14–S19, 1993.
N. Lippman, K.M. Stein, and B.B. Lerman. Comparison of methods for removal of ectopy in measurement of heart rate variability. Am J Physiol, 267(1):H411–H418, July 1994.
J.A. Lipponen and M.P. Tarvainen. A robust algorithm for heart rate variability time series artefact correction using novel beat classification. J Med Eng Technol, 43(3):173-181, 2019.
D.A. Litvack, T.F. Oberlander, L.H. Carney, and J.P. Saul. Time and frequency domain methods for heart rate variability analysis: a methodological comparison. Psychophysiol, 32:492–504, 1995.
J. Mateo and P. Laguna. Analysis of heart rate variability in the presence of ectopic beats using the heart timing signal. IEEE Trans Biomed Eng, 50(3):334–343, March 2003.
I.P. Mitov. A method for assessment and processing of biomedical signals containing trend and periodic components. Med Eng Phys, 20(9):660–668, November-December 1998.
S.W. Porges and R.E. Bohrer. The analysis of periodic processes in psychophysiological research. In J.T. Cacioppo and L.G. Tassinary, editors, Principles of Psychophysiology: Physical Social and Inferential Elements, pages 708–753. Cambridge University Press, 1990.
M.P. Tarvainen, P.O. Ranta-aho, and P.A. Karjalainen. An advanced detrending method with application to HRV analysis. IEEE Trans Biomed Eng, 49(2):172–175, February 2002.
Task force of the European society of cardiology and the North American society of pacing and electrophysiology. Heart rate variability – standards of measurement, physiological interpretation, and clinical use. Circulation, 93(5):1043–1065, March 1996.