Heart Rate Variability Analysis with the HRV Toolkit (2024)

I. Background

Heart rate variability (HRV) analysis attempts to assess cardiacautonomic regulation through quantification of sinus rhythmvariability. The sinus rhythm times series is derived from the QRS toQRS (RR) interval sequence of the electrocardiogram (ECG), byextracting only normal sinus to normal sinus (NN) interbeatintervals. Relatively high frequency variations in sinus rhythmreflect parasympathetic (vagal) modulation, and slower variationsreflect a combination of both parasympathetic and sympatheticmodulation and non-autonomic factors [1-5].

Traditional heart rate variability (HRV) measures are usually dividedinto two broad categories: time domain measures and frequency domainmeasures [3,4]. The time domain heart rate variability statisticscommonly calculated are defined in Table 1. Note,however, that computing pNNx with x< 50 ms in both long- and short-term recordings may provide a morerobust index of fluctuations due to vagal tone than the standard pNN50statistic [6]. Commonly used frequency domain measures are definedin Table 2.

The low frequency band (0.04 - 0.15 Hz) includes physiologicoscillations associated with baroreceptor reflexes and the highfrequency band (0.15 - 0.40 Hz) encompasses respiratory sinusarrhythmia. The powers in these bands has been used to provide indexesof autonomic function. Such measures must be interpreted withcaution, however. As noted, oscillations in the "low" frequency bands appear tobe mediated by parasympathetic and sympathetic components, while the"high" frequency power is mediated exclusively by the vagus.

Traditionally, frequency domain measures are calculated by resamplingthe original NN interval series and then applying the fast Fouriertransform or autoregressive spectral estimation (the maximum entropymethod). This resampling, however, can cause an attenuation in thehigh frequency components. If discontinuities exist in the NNinterval series, either because of the presence of abnormal beats orbecause of gaps or extreme noise in the original ECG recording,traditional approaches require either discarding the data or guessworkto estimate the locations of missing normal beats. To eliminate theneed for evenly sampled data required by Fourier or maximum entropymethods, frequency domain spectra can be calculated usingthe Lomb periodogram forunevenly sampled data [7,8] (the method used in this toolkit).

Although the long term (24-hour) statistics of SDANN, SDNNIDX and ULFpower can be calculated for shorter data lengths, they will becomeincreasingly unreliable. For short-term data (less than 15 minutes inlength), only the time domain measures of AVNN, SDNN, rMSSD and pNN50and the frequency domain measures of total power, VLF power, HF powerand LF/HF ratio should be used.

A number of the HRV measures are highly correlated with eachother. These include SDNN, SDANN, total power and ULF power; SDNNIDX,VLF power and LF power; and rMSSD, pNN50 and HF power. The LF/HF ratiodoes not correlate strongly with any other HRV measures [4].

Heart rate variability (HRV) has been widely applied in basic andclinical research studies. Its clinical application is verylimited at present, however. These limitations are due to lack ofstandardization of methodology and application to differentnon-comparable subsets of subjects, as well as to the confoundingeffects of age, gender, drugs, health status, and chronobiologicvariations, among others. Furthermore, outliers due to ectopy andartifact can have major effects on computed HRV values. In elderlysubjects, especially, a spuriously high value of certain measures maybe due to the effects of "erratic supraventricular rhythm" [9] due tosubtle atrial ectopy, wandering atrial pacemaker, or sinus nodeconduction abnormalities. Additional information on heart ratedynamics and analysis techniques, including non-linear and complexitybased measures, can be found in theHRV 2006 coursenotes and elsewhere on PhysioNet (see forexample: Detrended FluctuationAnalysis, Multiscale EntropyAnalysis, and Information-BasedSimilarity, among others).

PhysioNet's HRV Toolkit, available here, is a rigorously validatedpackage of open source software for HRV analysis, includingvisualization of NN interval time series, automated outlier removal,and calculation of the basic time- and frequency-domain HRV statisticswidely used in the literature, including all of those listed in thetables below.

Several other high-quality, freely available HRV toolkits may also beof interest to researchers; links to them are provided at the end ofthis page.

Selected References:

[1] Wolf MM, Varigos GA, Hunt D, Sloman JG. Sinus arrhythmia in acutemyocardial infarction. Med J Aust 1978;2:52-53.

[2] Kleiger RE, Miller JP, Bigger JT, Moss AJ, and the MulticenterPost-Infarction Research Group. Decreased heart rate variability and itsassociation with increased mortality after acute myocardial infarction.Am J Cardiol 1987;59:256-262.

[3] Task Force of the European Society of Cardiology and the NorthAmerican Society of Pacing and Electrophysiology. Heart rate variability:Standards of measurement, physiological interpretation, and clinical use. Circulation 1996; 93:1043-1065.

[4] Mietus JE. Timedomain measures: from variance to pNNx.http://physionet.org/events/hrv-2006/mietus-1.pdf

[5] Parati G, Mancia G, Di Rienzo M, Castiglioni P, Taylor JA, Studinger P.Point-Counterpoint: Cardiovascular variability is/is not an index ofautonomic control of circulation. J Appl Physiol 2006; 101: 676 - 682.

[6] Mietus JE, Peng C-K, Henry I, Goldsmith RL, Goldberger AL.The pNNx-files: Reexamining a widely-used heart rate variability measure.Heart 2002;88:378-380.

[7] Press WH, Teukolsky SA, Vetterling WT, Flannery BP. Numerical Recipes in C: The Art of Scientific Computing,2nd ed. Cambridge Univ. Press, 1992, pp. 575-584.

[8] Moody GB. Spectral analysisof heart rate without resampling. Computers in Cardiology 1993;715-718.

[9] Stein PK, Yanez D, Domitrovich PP, Gottdiener J, Chaves P, Kronmal R,Rautaharju P. Heart rate variability is confounded by the presence oferratic sinus rhythm. Computers in Cardiology 2002; 669-72.

II. Obtaining the HRV Toolkit

The HRV Toolkit: Contents

The HRV Toolkit consists of a Usages file,a Makefile, and the following programs:

plt_rrs and get_hrv are the main scripts used in calculatingHRV as illustrated below. These two scripts call the various C programs toaccomplish their calculations. The C programs can be used separately and theirusage and various options can be found in Usages,or by running any of these programs with the -h option.

Downloading and Installing the HRV Toolkit

Installing the HRV toolkit is easy:

1. On Windows, install the free Cygwin software first; see ourtutorial for details. Cygwin includesthe gcc C compiler and all other utilities needed to build and run thecomponents of the HRV Toolkit on Windows. Start a Cygwin (terminalemulator) window and use it for all of the remaining steps.
2. Install the freeWFDB andplt software packages.
3. The HRV toolkit is available as a tarball ofsources (for all platforms) or as tarballs of prebuilt binariesfor GNU/Linux (x86),Mac OS X (x86),Solaris (SPARC),or Windows/Cygwin. Download theversion of your choice.
4. Unpack the HRV toolkit tarball you downloaded. (See the PhysioNet FAQfor information on unpacking tarballs.)
5. If you downloaded the sources, enter the source directory (HRV.src)and compile and install the toolkit by typing:
   

    make install

   If you downloaded the binaries, move the contents of the HRV directoryinto some directory in your PATH (or add the HRV directoryto your PATH). The binaries require the same additional packages asthe source distribution (WFDB, plt, and (on Windows) Cygwin).

III. Using the HRV toolkit

User interface

Input data format

Both plt_rrs (for plotting the RR/NN interval series)and get_hrv (for calculating the HRV statistics) can take asinput either a PhysioBank-compatiblebeat annotation file or a text file containing an RR intervallist. RR interval lists can be in any of four formats:

• 3 columns (T, RR, A)
• 2 columns (RR, A)
• 2 columns (T, RR)
• 1 column (RR)

Although T is assumed to be expressed in seconds by default,the -I option of plt_rrs and get_hrv (seebelow) can be used if the RR interval list contains T values expressedas clock time (hh:mm:ss.xxx), hours (hh.xxxxxxx) or minutes(mm.xxxxx). Similarly, although RR is assumed to be expressed inseconds, intervals in milliseconds can be input by using the `-m'option. (See details below or in Usages).

Beat annotation files are available for most of the PhysioBank records thatinclude ECGs. For information about record and annotation conventions see thePhysioNet FAQ. If you wish tostudy a recording for which no beat annotation file or RR intervallisting is available, you may be able tocreate an annotation fileusing software from the WFDBsoftware package. Additional information on RR intervals, heart rate,and HRV can be found at the RR/HR/HRVHowto.

Input used in this tutorial

The examples below read the ecg annotations fromrecord chf03 of the BIDMCCongestive Heart Failure Database. To reproduce the results shownbelow, it is not necessary to download this annotation file, becauseapplications that use the WFDB library (including those in the HRVToolkit that read annotation files) can locate and read the copy fromthe PhysioNet web server if no local copy exists. In order to do sosuccessfully, it is necessary to specify the path to the record fromthe PhysioBank archive directory within the record name, as shown inthese examples (use chfdb/chf03 rather thansimply chf03). Examples that illustrate input from RRinterval lists assume that the input file is named chf03.rr,and that it is in the current (local) directory.

plt_rrs: plotting the RR/NN interval time series

NN interval outliers due to missed or false beat detections can seriouslycorrupt HRV statistics. Most frequency domain measures are especiallysusceptible to outliers, particularly LF and HF power which can be inerror by greater than 1000%. Most time domain measures are less affectedby outliers, but still can give errors in excess of 100%. AVNN, pNN50 andULF power are least affected, generally with errors less than 10% (seeTime domain Measures:From Variance to pNNx).

To reduce the errors due to outliers, the first step in HRV analysis is to visually examine the RR/NN interval data and determine the necessity ofoutlier filtering. This can be done using plt_rrs.

plt_rrs can be used with either an annotation file or an RR intervallisting (see below) and has the following options:

 plt_rrs [options] -R rrfile | record annotator [start [end]] Plot RR intervals or RR interval heart rates options : [-P 2|4|8|16|24|32] : plot 2, 4, 8, 16, 24 or 32 hours per page (default: scale page length to data length) [-R rrfile] : RR interval file : time (sec), interval [-N] : plot NN intervals instead of RR intervals [-H] : plot RR/NN interval heart rate [-F "filt hwin"] : filter intervals, plot filtered data [-f "filt hwin"] : filter intervals [-p] : plot points [-I c|h|m] : input time format: hh::mm:ss, hours, minutes (default: seconds) [-m] : RR intervals in msec [-y "ymin ymax"] : y axis limits [-o] : output postscript

 plt_rrs -o other arguments | lpr

Using an annotation file, like this:

 plt_rrs chfdb/chf03 ecg

To plot selected portions of the data, specify a start time and optionally an endtime. For example,

 plt_rrs -P 8 chfdb/chf03 ecg 00:00:00 01:00:00

With the -p option, as in

 plt_rrs -P 8 -p chfdb/chf03 ecg 00:00:00 01:00:00

With the -N option, as in

 plt_rrs -P 8 -N chfdb/chf03 ecg 00:00:00 01:00:00

The "NN : RR = 3823 : 4203 = 0.910 [380 non-NN]" in the titleindicates the total number of NN intervals, the total number of RRintervals, the fraction of RR intervals that are NN intervals and thenumber of non-NN intervals.

To filter outliers, use the -F or -f options, as in

 plt_rrs -P 8 -N -F "0.2 20 -x 0.4 2.0" chfdb/chf03 ecg 00:00:00 01:00:00

Using -F allows plt_rrs to plot the excluded intervals assmall filled circles:

The -f option suppresses output of excluded intervals:

The "Filt : NN : RR = 3762 : 3823 : 4203 = 984 : 0.910 = 0.895[61 Filtered, 380 non-NN]" in the title of the two plots above gives thetotal number of NN intervals remaining after filtering, the total number of NNintervals before filtering, the total number of RR intervals, the fraction of NNintervals remaining after filtering, the fraction of RR intervals that are NNintervals, the fraction of the total number of RR intervals that are NN intervalsremaining after filtering, and the number of NN intervals filtered out togetherwith the number of non-NN intervals.

Filtering the NN interval data may not be necessary if there are no extremeoutliers.

To plot the RR interval series of an RR interval file containing threecolumns of data : time in seconds, interval in seconds and beat label,use the -R option followed by the name of the file, rather thanspecifying the record and annotator names. For example, if the filechf03.rr contains:

1.356 0.952 N2.312 0.956 N3.252 0.940 N4.212 0.960 N5.156 0.944 N6.112 0.956 N7.048 0.936 N7.996 0.948 N8.960 0.964 N9.900 0.940 N...

The above RR interval list was generated using the command

 rrlist ecg chfdb/chf03 -s >chf03.rr

If an RR interval listing contains only time in seconds and intervalin seconds, e.g.,

1.356 0.9522.312 0.9563.252 0.9404.212 0.9605.156 0.9446.112 0.9567.048 0.9367.996 0.9488.960 0.9649.900 0.940...

To plot an RR interval listing containing only RR intervals in secondsand beat labels, e.g.,

0.952 N0.956 N0.940 N0.960 N0.944 N0.956 N0.936 N0.948 N0.964 N0.940 N...

Finally, if the RR interval listing contains only RR intervals in seconds,e.g.,

0.9520.9560.9400.9600.9440.9560.9360.9480.9640.940...

If the RR intervals in any of the above formats are listed in millisecondsrather than seconds, use plt_rrs with the -m option.

get_hrv: calculating the HRV statistics

To calculate the time and frequency domain HRV statistics use get_hrv,where the options are:

get_hrv [options] -R rrfile | record annotator [start [end]]Get HRV statistics : REC : NN/RR AVNN SDNN SDANN SDNNINDX RMSSD PNN : TOTPWR ULF VLF LF HF LF/HFoptions : [-R rrfile] : RR interval file : time (sec), interval [-f "filt hwin"] : filter outliers [-p "nndiff ..."] : nn difference for pnn (default: 50 msec) [-P "lo1 hi1 lo2 hi2 lo3 hi3 lo4 hi4"] : power bands(default : 0 0.0033 0.0033 0.04 0.04 0.15 0.15 0.4) [-s] : short term stats of REC : NN/RR AVNN SDNN RMSSD PNN : TOTPWR VLF LF HF LF/HF [-I c|h|m] : input time format: hh::mm:ss, hours, minutes (default: seconds) [-m] : RR intervals in msec [-M] : output statistics in msec rather than sec [-L] : output statistics on one line [-S] : plot HRV results on screenplotting options : [-F "filt hwin"] : filter outliers, plot filtered data [-y "ymin ymax"] : time series y-axis limits ("- -" for self-scaling) [-X maxfreq] : fft maximum frequency (default : 0.4 Hz) [-Y fftmax] : fft maximum ("-" for self-scaling) [-o] : output plot in postscript

The command

 get_hrv -f "0.2 20 -x 0.4 2.0" -p "20 50" chfdb/chf03 ecg

 chfdb/chf03 : NN/RR = 0.942899 AVNN = 0.892206 SDNN = 0.054712 SDANN = 0.0466773 SDNNIDX = 0.0241124 rMSSD = 0.0177694 pNN20 = 0.0688698 pNN50 = 0.0252389 TOT PWR = 0.0033882 ULF PWR = 0.00270489 VLF PWR = 0.00039018 LF PWR = 0.000124104 HF PWR = 0.000169027 LF/HF = 0.734226

By default, get_hrv outputs the HRV statistics in seconds. To outputresults in milliseconds use the -M option, as in

 get_hrv -M -f "0.2 20 -x 0.4 2.0" -p "20 50" chfdb/chf03 ecg

 chfdb/chf03 : NN/RR = 0.942899 AVNN = 892.206 SDNN = 54.712 SDANN = 46.6773 SDNNIDX = 24.1124 rMSSD = 17.7694 pNN20 = 6.88698 pNN50 = 2.52389 TOT PWR = 3388.65 ULF PWR = 2705.19 VLF PWR = 390.25 LF PWR = 124.156 HF PWR = 169.057 LF/HF = 0.734403

To output the HRV statistics all on one line, use the -L option:

 get_hrv -L -M -f "0.2 20 -x 0.4 2.0" -p "20 50" chfdb/chf03 ecg

 chfdb/chf03 : 0.942899 892.206 54.712 46.6773 24.1124 17.7694 6.88698 2.52389 : 3388.65 2705.19 390.25 124.156 169.057 0.734403

 REC : NN/RR AVNN SDNN SDANN SDNNINDX RMSSD PNN20 PNN50 : TOTPWR ULF VLF LF HF LF/HF

To limit the output to the statistics named as short-term HRV statistics intables 1 and 2 above, use the -s option, as in

 get_hrv -s -M -f "0.2 20 -x 0.4 2.0" -p "20 50" chfdb/chf03 ecg 0:00:00 1:00:00

 chfdb/chf03 : NN/RR = 0.90173 AVNN = 867.958 SDNN = 39.7191 rMSSD = 20.8659 pNN20 = 6.26553 pNN50 = 2.20811 TOT PWR = 1702.96 VLF PWR = 1423.01 LF PWR = 103.941 HF PWR = 176.008 LF/HF = 0.590547

To output the short-term HRV statistics all on one line, use both the-s and the -L options:

 get_hrv -L -s -M -f "0.2 20 -x 0.4 2.0" -p "20 50" chfdb/chf03 ecg 0:00:00 1:00:00

 chfdb/chf03 : 0.90173 867.958 39.7191 20.8659 6.26553 2.20811 : 1702.96 1423.01 103.941 176.008 0.590547

 REC : NN/RR AVNN SDNN RMSSD PNN20 PNN50 : TOTPWR VLF LF HF LF/HF

As with plt_rrs, get_hrv with the -Roption can accept a text file containing an RR interval listing in anyof the formats described above. If the interval times are notlisted, get_hrv will calculate the time of the RR intervalfrom the RR interval sequence, and if there are no beat labels,get_hrv -R will treat all intervals as normal.

If the RR intervals are listed in milliseconds rather than seconds, useget_hrv with the -m option.

To plot the HRV results in a window including the NN intervaltime series, the NN interval histogram and the NN interval power spectrumuse the -S option to get_hrv.As with plt_rrs the -o will output PostScript, and the-F filtering option will plot the excluded intervals of the NNinterval time series: the non-normal beats as small open circles andoutliers as small filled circles.For example

 get_hrv -S -s -M -F "0.2 20 -x 0.4 2.0" -p "20 50" chf03 ecg 0:00:00 1:00:00

Note that plots made this way show the calculated HRV statistics with reducedprecision.

Additional get_hrv plotting options include -y for scalingthe time series y-axis limits, the -X option for scaling the powerspectrum x-axis limit and the -Y option for scaling the power spectrumy-axis limit.

IV. Representative measurements of HRV

Table 3 gives representative values of HRV measurements, obtainedusing the HRV toolkit from the 72 subjects in theMIT-BIH Normal Sinus RhythmDatabase and the Normal SinusRhythm RR Interval Database. These representative values(comparable to those previously reported in healthy subjects) are not,however, intended as a standard normative database.

HRV measurements obtained for each of these 72 subjects are availablehere.

Links to other HRV toolkits

Several other free HRV toolkits exist. Among them are:

• Physionet Cardiovascular Signal Toolbox (https://physionet.org/physiotools/physionet-cardiovascular-signal-toolbox/)

  An open-source modular program for calculating heart rate variability (HRV) implemented in Matlab with evidence-based algorithms and output formats. The Toolbox is compatible with 64-bit MATLAB on GNU/Linux, Mac OS X, and MS-Windows. It was compared to several other open source and proprietary tools including the PhysioNet C HRV Toolkit to create a benchmark for the field. It was shown to be equivalent to the PhysioNet C HRV Toolkit. It contains the most extensive set of tools in any HRV algorithm collection so far published. It has no dependencies outside of Matlab (tested on Matlab R2017a and R2017b).

• KUBIOS-HRV (http://kubios.uku.fi/)

  Authors: JP Niskanen and colleagues at Kuopio University, Finland

  KUBIOS-HRV 2.0 is available for Windows and Linux. It is astandalone Matlab application (only the Matlab runtime compiler andlibraries, both included in the package, are needed in order to useit). KUBIOS-HRV includes a well-written 53-page user guide withdetailed descriptions of all of the algorithms used, mathematicalderivations where needed, and an extensive bibliography.

  The package includes automated and interactive methods for outlieridentification and removal; a smoothness-based detrending algorithmfor preprocessing time series; derivation of SDNN, SDSD, RMSSD, pNN50,HRV triangular index, and TINN; FFT and AR power spectrum estimates, VLF,LF, HF, and LF/HF ratio; and Poincaré plots, ApEn, SampEn, DFA,correlation dimension, and recurrence plots. Inputs are RR intervallists in text format.

  KUBIOS-HRV 2.0 and its user guide were released in 2002, and it wasdescribed in detail in:

  Juha-Pekka Niskanen, Mika P Tarvainen, Perttu O Ranta-aho, Pasi A Karjalainen,"Software for advanced HRVanalysis." Computer Methods and Programs in Biomedicine 76(1):73-81(2004 Oct; DOI: 10.1016/j.cmpb.2004.03.004)

• R-HRV (http://cran.r-project.org/web/packages/RHRV/)

  Authors: M Lado, A Mendez, L Rodriguez-Linares, X Vila (Univ. of Vigo, Spain)

  R-HRV is freely available in source form (for all platforms) and asprebuilt binaries for Mac OS X and Windows, together with a 13-pageuser guide. It requires R (a very widely used open-source statisticspackage). R-HRV includes functions for importing PhysioBank(WFDB) annotation files, automated outlier handling, and frequency domain HRVanalysis.

  R-HRV was described in:

  L. Rodriguez-Linares, X. Vila, A. Mendez, M. Lado, D. Olivieri, "RHRV: AnR-based software package for heart rate variability analysis of ECGrecordings." 3rd Iberian Conference in Systems and Information Technologies(CISTI 2008), 21-23 June 2008.

• Software for Heart Rate Variability (http://www.macalester.edu/~kaplan/hrv/doc/)

  Authors: Danny Kaplan and Phil Staffin (Macalester College)

  This software, written as Matlab m-code, includes routines for outlierhandling, SDNN, RMSSD, a generalization of pNN50, spectral analysisusing FFT (VLF, LF, HF, LF/HF), ApEn, and DFA. Inputs are RR intervallists in text format.