I've used sometimes a simple approach, though probably could be some more better or more efficient. Assuming that you have a sum of delayed deltas in time, with unit area, you can build its spectrum as the sum of frequency exponentials (with exponent coefficient given by the delay-time of each impulse, according to continuous-time Fourier Transform of a delayed delta). As far as this is complex, you should just sample in frequency (for getting discrete-Fourier) and take the modulus. Alternatively, you could be thinking in using Lomb periodogram (just for comparing HRV interpolated, spectrum of counts and Lomb, and see the similar and the difference). More sophisticated non-uniform interpolation and spectrum algorithms exist, which are seldom used in HRV.

I agree the Jose sir method which use Lomb periodogram to estimate the spectrum. Lomb periodogrm did not consider sampling rate problem.

Thank you for your advice. I used also the Lomb periodogram, but actually when I use it I need to have values of the samples of the signal. The one values of the Dirac impulses are no good for this. So I need to calculate the differences between the times of the impulses and technicality I get an Interval function ( a series of differences of the times allocated in the time of the second R wave of the given pair of subtracted R waves) whit this . I already analyzed the interval function and I just try to analyze the event series directly for comparison.

I also tried by now a direct approach calculating the equation:

D(jω)=Σ (cos(ωti) - j*sin(ωti))

Which is derived from the direct FT of the series of Dirac impulses knowing that the integral has a meaning just in the times of the Dirac impulses and in this times is equal to e-jωt= cos(ωt) - j*sin(ωt). Is this what you meant Mr. Rojo-Álvarez?

It does show the frequency of the cosine wave which I used, but whit a lot of disturbance. Also the spectrum of the unidirectional component of my signal ( 1+0.5*cos(2*pi*f ) ) is calculated as a spectrum of a rectangular signal for some unknown reason.

The code in matlab is:

post=0;

for frek=0:0.001:0.5

post=post+1;

PP(1,post)=sum(cos(2*pi*frek.*time)-i*sin(2*pi*frek.*time));

end

It is like what you said, D(jw) is not enough for giving you a spectral description, though the equation is OK (you can find some corrections to the basic ecuation you're using in the literature). With respect to the code lines, I'm afraid I'm not sure to understand them. Maybe the following script could give you a comparison of Lomb with "alone" spectrum of counts:

@@@@@@@@@@@@@@@@@@@@@@@@@@@

function ExampleLombCounts

%

% Generate random data

nsamples = 50;

tt1 = rand(nsamples,1); % in ms

tt = cumsum(tt1); % in sec, irregular sampling

f = 0.25; % Hz

RR = (1+.5*sin(2*pi*f*tt));

ttRR = cumsum(RR);

figure(1), subplot(311), stem(ttRR,RR), xlabel('s'), ylabel('s')

% Point spectrum

ww = 2*pi*linspace(0,2,1024);

S = zeros(size(ww));

for k=1:length(tt)

S = S + exp(-j*ww*(tt(k)-25));

end

subplot(312), plot(ww/2/pi,abs(S)); xlabel('Hz')

% with Lomb

[P,f,alpha] = lomb(RR,tt);

subplot(313), plot(f,P), xlabel('Hz'), axis tight

@@@@@@@@@@@@@@@@@@@@@@@@@@@

(you should use lomb.m, the normalized Lomb from matlab file exchange site).

You'll see that the spectrum of counts alone is not giving you a good estimate. Why? Maybe you can find useful the following reference (if you're not using it already):

http://elib.uni-stuttgart.de/opus/volltexte/2010/5239/pdf/nag23.pdf

Where there is detailed equations on the relationship between spectrum of counts and RR spectrum.

Best.

Thank you for your help.

To the code I send:

post=0; % just an index variable

for frek=0:0.001:0.5 % changing the frequency for which I count the summation from 0 to 0.5 by step of 0.001.

post=post+1; % increasing the index variable

PP(1,post)=sum(cos(2*pi*frek.*time)-i*sin(2*pi*frek.*time)); %counting the summation from the theoretical equation for selected frequency

end

The variable time is hereby a vector of times were the Dirac impulses are in the series.

Thank you also for the code you send. Could you please give me a link on the lomb.m function you’re using?

Thank you.

Yes I basically am using the same article just published in a different Journal.

Thank you for the code. I just don`t understand why there is the -25 in:

S = S + exp(-j*ww*(tt(k)-25));

You're right, it's spureus code :) just was making previously some examples with the same code, and tried to make about-symmetric in time (from 0 to 50 secs, compared to -25 to 25 secs). But will work the same just removing -25. Sorry for the lapsus.

Ocurrence of rectangular signal spectrum menntioned above (sinc function) can be explained by time truncation (finite duration) effect. Spectrum of counts of time unlimited signal has Dirac pulse at zero frequency. After time truncation, spectrum is convolved in freq. domain with Fourier transform of rectangular function. I discussed this effect in

https://www.researchgate.net/publication/234040105_Spectrum_of_counts_computation_for_HRV_analysis

Regarding spectrum of conts versus Lomb, Sp. of counts is closely related to IPFM model and shows spectrum of modulating signal. This is not true for Lomb, but Lomb spectrum can look "better"- with less disturbances.

Conference Paper Spectrum of counts computation for HRV analysis

I had the same problem when I worked with ECG signals and a decided to compute the Fourier transform of a comb. I published a paper that I can send you

[email protected]

About the Lomb periodogram, I do not like it, because it is a sub-optimum algorithm. It is preferable to use the Fourier transform directly, by putting dirac deltas at the sample positions.