Hello,

I wanted to calculate psd according to the Welch Method and the Lomb Method. Therefore, I would like to ask the following two questions:

1. Function: lomb_psd
   Why is the power at f=0 fixed to 2 (powers[0] = 2)? My psd figure has always a to high y-Axis such that the the results won't be plotted.

2. Just a general question: Why are the max psd values so different between Welch and Lomb?

Thank you in advance

In the time domain documentation (time.rst)
The output keys mentioned for triangular_index are incorrect. They are given as 'tinn_histogram','tinn'....etc
While they should be
'tri_histogram', 'tri_index'

Goal:
[ ] Fixing Documentation

Hi Pedro,

just wanted to verify my own Monkey C implementation of Pointcareé Plot with your data and results. Found major differences, as seen in the following screenshots legends (sd1 and sd2):

in nonlinear.py there is the poincareé plot with the additional ellipse drawn from SD1 and SD2 standard deviation of the points distance to the two axes.

# SD1 & SD2 Computation
	sd1 = np.std(np.subtract(x1, x2) / np.sqrt(2))
	sd2 = np.std(np.add(x1, x2) / np.sqrt(2))

Both arrays are having issues and therefore both standard deviations are based on wrong values.
sd1: when x1-x2 is calculated there can be negative numbers which are not to be treated as distance to the 45 degree axe
sd2: adding x1+x2 is not at all giving the distance to the -45 degree axe since ignoring the position of the centroid point

Correct implementation as follows:

# Load NNI sample series
nni = pyhrv.utils.load_sample_nni()
nn = pyhrv.utils.check_input(nni)

# Prepare Poincaré data
x1 = np.asarray(nn[:-1])
x2 = np.asarray(nn[1:])

# SD1 & SD2 Computation
# sd1 = np.std(np.subtract(x1, x2) / np.sqrt(2))
sd1 = np.std(np.abs(np.subtract(x1, x2)) / np.sqrt(2))

# Initialization of an empty vector for sd2
sd2_values = []
mean_val = np.mean(nn)
x_center = mean_val
y_center = mean_val

# Loop through all elements of x2
for i in range(len(x2)):
    x_beatSample = x1[i]
    y_beatSample = x2[i]
    x_beatMirrorSample = x_center - (y_beatSample - y_center)
    y_beatMirrorSample = y_center - (x_beatSample - x_center)
    x_diff = x_beatMirrorSample - x_beatSample
    y_diff = y_beatMirrorSample - y_beatSample
    distancePoint_i = 0.5 * np.sqrt(x_diff**2 + y_diff**2)
    # Add the sd2 value to the list
    sd2_values.append(distancePoint_i)

# Output results
print("nn", nn)
print("distances used for SD1:", np.abs(np.subtract(x1, x2)) / np.sqrt(2))
print("distances used for SD2:", sd2_values)

sd2 = np.std(sd2_values)

print("SD1:", sd1);
print("SD2:", sd2);

Hope you or someone else can comment if my suggested implementation is correct?

Regards - Klaus.

P.s.: My Garmin Watch app is online since last december: --> https://apps.garmin.com/en-US/apps/1c02ac65-70c5-440e-a6e6-f4a706c8deb0

In dfa, the NN interval argument is nn:

However, the documentation says the NN interval argument is nni: https://pyhrv.readthedocs.io/en/latest/_pages/api/nonlinear.html#detrended-fluctuation-analysis-dfa

python3 setup.py install
ImportError: cannot import name 'ImageTk' from 'PIL' (/usr/lib/python3/dist-packages/PIL/init.py)

having issue to run overall code in colab and also in jupyter notebook? how can i solve this issue ? is python version matter? can anybody tell how can i solve this issue?

The Polar OH1 optical heart rate sensor has ECG accuracy:

https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0217288

and R-R peaks in the sensor stream:

https://towardsdatascience.com/creating-a-data-stream-with-polar-device-a5c93c9ccc59

Furthermore, the Polar OH1 is comfortable to wear and affordable. Sensor data can be obtained with the HR & HRV Logger per Polar H6, H7, H10 e OH1 app:

https://play.google.com/store/apps/details?id=hrmonitor.com.hr_rri_monitor&gl=CH

but does this app also work with the Python HRV toolbox?

Hi Pedro, first of all - thanks for the great project!

I'm currently trying to check how well the biosppy+pyHRV analysis stacks up vs Kubios for a dataset of mine, and I'm having some trouble to reach equal results for FFT and time domain features.

I have seen Table 6 - Parameter Settings from your thesis, and tried to implement all the settings as described there. I'm not sure if I can match them exactly, as FFT - Interpolation frequency, FFT - Window width, and FFT - Window overlap are not exposed through kwargs_welch. However, I checked the source code, and since noverlap=nperseg//2 per default and the sampling frequency is 4, only the window width needs to be set to the same value. I tried with nperseg=300 as set in the source code, as well as nperseg=1200, as I thought that with a sampling frequency of 4Hz you'd need more values in the segments to match 300 seconds, but still, my results are off.

Can you help me? I've kept the other parameters on their default, as I assume they are at the value you used for the verification in the thesis. Could this be due to the R-Peak extraction by biosppy?

I am using python 3.10.2 under windows 10. I am trying to change the frequency bands of the pyhrv welch method and I am getting an UnbondLocalError.

import pyhrv
rr = [1,2,3,4,5,6,7,8,9,10]
pyhrv.frequency_domain.welch_psd(rr, fbands={'vlf': (0.003, 0.04), 'lf': (0.04, 0.15), 'hf': (0.15, 0.4)})

Error obtained:

>>> pyhrv.frequency_domain.welch_psd(rr, fbands={'vlf': (0.003, 0.04), 'lf': (0.04, 0.15), 'hf': (0.15, 0.4)})
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Users\ggaleote\Documents\git\sinus-rhythms-pathologies\.venv\lib\site-packages\pyhrv\frequency_domain.py", line 159, in welch_psd
    fbands = _check_freq_bands(fbands)
  File "C:\Users\ggaleote\Documents\git\sinus-rhythms-pathologies\.venv\lib\site-packages\pyhrv\frequency_domain.py", line 700, in _check_freq_bands
    args_.append(ulf)
UnboundLocalError: local variable 'ulf' referenced before assignment

The currently available tinn() function (pyhrv.time_domain.tinn()) computes wrong results and is time consuming.

Revision goals:
[ ] Fix bug that causes the computation of wrong parameters
[ ] Improve performance/decrease computational time

Hi,

I am using pyhrv to extract features from (my) rr-intervalls. Now I am getting this error often:

path/to/file/frequency_domain.py:947: UserWarning: Attempting to set identical low and high ylims makes transformation singular; automatically expanding.
ax.axis([0, fbands['hf'][1], 0, max(power) * (1 + 0.05)]).

I am a bit lost if this means but I get this error mostly when I am sitting still. If I move and my HR goes up the error vanishes.

in triangular_index() and tinn(), the default value for binsize was set to 7.815, however the correct value should be 7.8125(1000/128).

Thanks for your great work. I got the following error, running
results = pyhrv.hrv(nni=nni_base) on Python3.7

Would you please help me.

KeyError Traceback (most recent call last)
in
----> 1 results = pyhrv.hrv(nni=nni_base)

~/venv/lib/python3.7/site-packages/pyhrv/hrv.py in hrv(nni, rpeaks, signal, sampling_rate, interval, plot_ecg, plot_tachogram, show, fbands, kwargs_ecg_plot, kwargs_tachogram, kwargs_time, kwargs_nonlinear, kwargs_welch, kwargs_lomb, kwargs_ar)
284 # COMPUTE FREQUENCY DOMAIN RESULTS (kwargs are verified by the frequency_domain() function)
285 f_results = fd.frequency_domain(nni=nn, fbands=fbands, kwargs_welch=kwargs_welch, kwargs_lomb=kwargs_lomb,
--> 286 kwargs_ar=kwargs_ar, show=False)
287
288 # COMPUTE NONLINEAR PARAMETERS

~/venv/lib/python3.7/site-packages/pyhrv/frequency_domain.py in frequency_domain(nni, rpeaks, signal, sampling_rate, fbands, show, show_param, legend, kwargs_welch, kwargs_lomb, kwargs_ar)
1125 else:
1126 # Compute Autoregressive PSD with default values
-> 1127 ar_results = ar_psd(nn, show=False, fbands=fbands, legend=legend, show_param=show_param)
1128
1129 # If plots should be shown (show all plots at once)

~/venv/lib/python3.7/site-packages/pyhrv/frequency_domain.py in ar_psd(nni, rpeaks, fbands, nfft, order, show, show_param, legend, mode)
497
498 # Plot PSD
--> 499 figure = _plot_psd('ar', frequencies, powers, freq_i, params, show, show_param, legend)
500 figure = utils.ReturnTuple((figure, ), ('ar_plot', ))
501

~/venv/lib/python3.7/site-packages/pyhrv/frequency_domain.py in _plot_psd(method, freq, power, freq_indices, parameters, show, show_param, legend)
904 ax.set_title("PSD - Welch's Method")
905 elif method == 'ar':
--> 906 ax.set_title("PSD - Autoregressive (Order %i)" % parameters['ar_order'])
907 elif method == 'lomb':
908 ax.set_title("PSD - Lomb-Scargle Periodogram")

~/venv/lib/python3.7/site-packages/biosppy/utils.py in getitem(self, key)
405 if isinstance(key, six.string_types):
406 if key not in self._names:
--> 407 raise KeyError("Unknown key: %r." % key)
408
409 key = self._names.index(key)

KeyError: "Unknown key: 'ar_order'."

Hi, I've tried to use the frequency_domain.frequency_domain() function as follow:
test = np.array([...]) res = frequency_domain.frequency_domain(signal=test)
But it raised this error: AttributeError: module 'biosppy.inter_plotting.ecg' has no attribute 'ecg'. Could you please help me? :)

Hi - I'm working with a Polar H10 heart rate monitor, and I'm wondering if it's possible to use this device with PyHRV. I can pipe the H10 into LSL and save as an XDF, but trying to understand if you've seen anyone analyze data off a Polar monitor using this tool. Thanks!

** Update **
So it seems as if this is not possible. Trying to use PyHRV in a real-time context appears challenging because the module relies on an OpenSignals ECG txt file as input. Converting a real-time stream into an OpenSignals text file and then porting that back into PyHRV probably too much overhead for real-time use.

Which brings me to another question, is there any way to use the PyHRV functions on numpy arrays? The data I read off a polar monitor is pretty simple (timestamp and microvoltage), and it would be great to be able to feed the voltage data +/- the sampling frequency (130hz) into PyHRV for real-time analysis without needing to convert each data chunk into a text file for processing.

Hello Pedro, congrats for your work!

I think there is something wrong with the time_domain.heart_rate function when I use NNI as input. When I use the biosppy.ecg function and calculate the mean of my input signal the result is 80.24485892889525. With your time_domain function the mean result is 160.4897178577905 (exactly twice).

Also, I think that the documentation of this function (heart_rate) is not right.
Currently return is: hr (array): Series of NN intervals in [ms].
While I think should be: hr (array): Series of Heart Rate in [bpm].

Thank you again Pedro!

When I use Welch's periodogram method, which suggested window width and overlapping rate in frequency analysis for 5 minutes of HRV?

We compare ECG results of pyHRV and medilog DARWIN2 (https://www.schiller.ch/en/products/medilog-darwin2-p198). The time domain is almost the same. But frequency domain LF/HF ratio has 20% difference. Do you have any suggestion? We use 5 minutes duration now.

As the title says, the nonlinear.dfa function does always show plots, irrespective of the argument given to show.
This is at least true in an ipyhton environment.
I fixed it locally, but I'm not sure whether the current behavior is desired or not. If it is not, I can submit a fix.

D:\pyhrv>python3 setup.py bdist_wheel
Library mydpss not found
Traceback (most recent call last):
File "D:\pyhrv\setup.py", line 8, in
long_description = fh.read()
UnicodeDecodeError: 'cp950' codec can't decode byte 0xe2 in position 3044: illegal multibyte sequence

When I try to plot tachogram I get this error, please help.

tools.tachogram(signal, interval=[0, 15])

TypeError Traceback (most recent call last)
in
1 # Plot the Tachogram of the ECG signal (interval: 0s to 15s)
----> 2 tools.tachogram(signal, interval=[0, 15])

/opt/miniconda3/lib/python3.7/site-packages/pyhrv/tools.py in tachogram(nni, signal, rpeaks, sampling_rate, hr, interval, title, figsize, show)
373 interval = [0, t[-1]]
374 else:
--> 375 interval = pyhrv.utils.check_interval(interval, limits=[0, t[-1]], default=[0, 10])
376
377 # Prepare figure

/opt/miniconda3/lib/python3.7/site-packages/pyhrv/utils.py in check_interval(interval, limits, default)
228 limits = default
229 elif limits is not None:
--> 230 limits = _check_limits(limits, 'limits')
231
232 # Check interval limits

/opt/miniconda3/lib/python3.7/site-packages/pyhrv/utils.py in _check_limits(interval, name)
283 interval[0], interval[1] = interval[1], interval[0]
284 vals = (name, name, interval[0], interval[1])
--> 285 warnings.warn("Corrected invalid '%s' limits (lower limit > upper limit).'%s' set to: %s" % vals)
286 if interval[0] == interval[1]:
287 raise ValueError("'%f': Invalid interval limits as they are equal." % name)

TypeError: not all arguments converted during string formatting

Hi!
When I am processing my signal I need the returned frequency data, but not necessarily the figures.
This code plots a nice figure in my notebook:

result = fd.welch_psd(intervals, show = False)

Can you help with this?

Thanks a lot for this package!
Orsi

Thank you for the helpful toolkit. When I ran the function "frequency_domain" in a loop , I got a error reported as "Fail to allocate bitmap". However, I still have a lot of available memory on my machine. Does anyone know why?

pip install pyhrv
Defaulting to user installation because normal site-packages is not writeable
Collecting pyhrv
Using cached pyhrv-0.4.1-py3-none-any.whl (3.2 MB)
Collecting biosppy (from pyhrv)
Using cached biosppy-1.0.0-py2.py3-none-any.whl (106 kB)
Requirement already satisfied: matplotlib in /home/andyueng/.local/lib/python3.10/site-packages (from pyhrv) (3.7.1)
Requirement already satisfied: numpy in /home/andyueng/.local/lib/python3.10/site-packages (from pyhrv) (1.24.4)
Requirement already satisfied: scipy in /home/andyueng/.local/lib/python3.10/site-packages (from pyhrv) (1.11.1)
Collecting nolds (from pyhrv)
Using cached nolds-0.5.2-py2.py3-none-any.whl (39 kB)
Collecting spectrum (from pyhrv)
Using cached spectrum-0.8.1.tar.gz (230 kB)
Preparing metadata (setup.py) ... error
error: subprocess-exited-with-error

× python setup.py egg_info did not run successfully.
│ exit code: 1
╰─> [12 lines of output]
Traceback (most recent call last):
File "", line 2, in
File "", line 14, in
File "/home/andyueng/.local/lib/python3.10/site-packages/setuptools/init.py", line 20, in
from setuptools.dist import Distribution, Feature
File "/home/andyueng/.local/lib/python3.10/site-packages/setuptools/dist.py", line 35, in
from setuptools.depends import Require
File "/home/andyueng/.local/lib/python3.10/site-packages/setuptools/depends.py", line 7, in
from .py33compat import Bytecode
File "/home/andyueng/.local/lib/python3.10/site-packages/setuptools/py33compat.py", line 55, in
unescape = getattr(html, 'unescape', html_parser.HTMLParser().unescape)
AttributeError: 'HTMLParser' object has no attribute 'unescape'
[end of output]

note: This error originates from a subprocess, and is likely not a problem with pip.
error: metadata-generation-failed

× Encountered error while generating package metadata.
╰─> See above for output.

note: This is an issue with the package mentioned above, not pip.
hint: See above for details.

import biosppy
import numpy as np
import pyhrv.tools as tools
from opensignalsreader import OpenSignalsReader 

# Load sample ECG signal & extract R-peaks using BioSppy
signal = OpenSignalsReader('./samples/SampleECG.txt').signal('ECG')
signal, rpeaks = biosppy.signals.ecg.ecg(signal, show=False)[1:3]

# Compute NNI
nni = tools.nn_intervals(rpeaks)

pyhrv.tools.nn_intervals() takes in rpeaks (array): R-peak times in [ms] or [s].

The problem is that the rpeaks parameter returned by biosppy.signals.ecg.ecg() does NOT contain R-peak times in [ms] or [s], instead it contains the integer indices of the R-peaks.

Our research team had been working with faulty HRV calculations for half a year before we just now noticed the issue...

thanks for the great work. I found your repo while look for the SD1/SD2 implementation of poincare plot for my heart rate data in python.

It is great to see that you are making standard lib for the hrv calculation for python3. Can you give some hints about the min sample rate of heart rate for your api, since i am going to use it for the data optained by smartwatches, they have different sample rate als clinic ecg.

It will be really great if you can added the Lyapunov Exponent (LE) to the non-linear feature. It is always a great hiccup to get the LE correctly implemented.

Thanks again for keeping up this great work.

FYI - the samples folder does not contain any files.

Thank you.

Hey devs,

Throughout the nonlinear.py codes, import biosppy (line 45) and then calling biosppy.utils.ReturnTuple(args, names) throws the following error:
AttributeError: module 'biosppy' has no attribute 'utils'

Installed libraries

biosppy==2.2.1
pyhrv==0.4.1

Proposed fix

I am circumventing this by directly importing the ReturnTuple class from biospyy.utils:
from biosppy.utils import ReturnTuple
ReturnTuple(args, names)

Please find enclosed PDF file which compares some frequency domain plots for pyHRV vs. Kubios Standard.
I am using the "SampleECG.txt" data from pyHRV.

Do you have some insight regarding the different appearance of the plots between the two softwares?
What could I tune in pyHRV functions to make a similar looking plot as displayed in Kubios Standard?

Thank you,
Ben DeKoven
pyHRV_27May2021.pdf

Hi ,

I am trying to stop plotting from pyhrv.hrv(show=False) but it does not work and shows all the plots! I am using 0.4 version

I have gotten some feedback from some users that there is an issue with the nn_intervals() function which leads to the following error message:

"Incompatible data type in list or numpy array (only int or float allowed)"

This happened because of a logic flaw I had while trying to make pyHRV overly robust and this will be fixed with the next update.

If you are looking for an intermediate solution, please delete the lines 97 & 98 from the tools.py module.

Apologies for the issues this may have been causing to some users!

Hi Understand how to create a plot as:

results = nl.poincare(nni, show=False)

which returns a matlibplot figure...but how do I then either reposition the plot and add more plots to it?

Check condition in line 719:

File "/usr/local/lib/python2.7/dist-packages/pyhrv/time_domain.py", line 255, in sdnn_index
segments, seg = tools.segmentation(nn, full=full, overlap=overlap, duration=duration)
File "/usr/local/lib/python2.7/dist-packages/pyhrv/tools.py", line 661, in segmentation
csegment.append(nn[cindex])
IndexError: index 427 is out of bounds for axis 0 with size 427

Hi, I am working with PPG signal and I wonder why the interpolation of nni is made with 4Hz and what kind of interpolation does frequency domain functions use. I didn´t find in the docs any information about this 4 Hz value.

Hello Pedro,

When I try to run the frequency domain Welch PSD example with SampleECG.txt file I got a different plot figure

Example code:

`
import biosppy
import pyhrv.frequency_domain as fd
from opensignalsreader import OpenSignalsReader

signal = OpenSignalsReader('SampleECG.txt').signal('ECG')
rpeaks = biosppy.signals.ecg.ecg(signal)[2]
result = fd.welch_psd(rpeaks=rpeaks)`

When the result should be:

Hi and thank you very much for your extremely helpful toolbox. I have data from a smartwatch (Gear S3) obtained at 5Hz. The data consists of the Timestamp of each sample, the beats per minute, and the RRinterval of the last detected pulse (Image)

Due to the fact that the signal is sampled at 5Hz, if a new pulse has not been detected then we get the same value multiple times as can be seen in the Figure. One quick way to get the true RRI series would be to drop consecutive values, however we have seen many cases where the output of the watch is the same for a long period (e.g. 5 seconds). We do not know if its because the watch hasn't detected a new pulse, or if the rRInterval detected is exactly the same with the previous.

Do you know if there is a better way/algorithm to extract/approximate the true sequence (apart from dropping consecutive values or some of them) ? Thanks!

Hi, thanks for the great repo which makes life much easier.
I guess there is one import line missing here so that it will throw error when the signal is passed:
from biosppy.signals.ecg import ecg

pyhrv.nonlinear.poincare(nni = data, show=False) seems to not disable the visual output.
Is there any other way to disable the visualization?

pyhrv.report.hrv_report(results, path='C:\ecgt\pyhrvwork', rfile='SampleReport', file_format='csv')

AttributeError: module 'pyhrv' has no attribute 'report'

pyhrv.tools.hrv_report(results, rfile='SampleReport', path='/my/favorite/path/')

same error found in both type of syntax

help me out

Hello,
I have been trying to get the same output as kubios outputs when using the pyhrv pyhrv.frequency_domain.welch_psd() method. I have the default setting on kubios with nothing changed and for the Welch method I left everything on default except for detrend which I set to false. I wondering if you have any tips on producing similar results. I have copied the method and change "fs" to 1 which gave me results close to kubios, but different enough for hf% and lf% not to be close.

thanks!