-->stacksize(20000000) -->load signals; -->ecg = signals(12);We can see that the sampling rate was 100 Hz
-->duration
duration =
10.
-->nsamples(12)
ans =
1000.
-->duration/nsamples(12)
ans =
0.01
From an EDF file we can extract information to scale the signal with the next code
-->physmin(12)
ans =
- 1920.2
-->physmax(12)
ans =
1926.
-->digmin(12)
ans =
- 2048.
-->digmax(12)
ans =
2047.
-->physdim(12)
ans =
uV
The result means that a value of -2048 corresponds to -1926 uV and a value of 2047 to 1926 uV. From these values, the amplitude of the signal can be corrected. But now we are going to calculate the power spectrum from the original values.
We are going to establish an epoch of one second corresponding to 100 values. We are going to build a vector containing the frequency scale
-->inct = 0.01
inct =
0.01
-->t = [1:100]*inct;
-->incf = 1/(100*inct)
incf =
1.
-->fre = [0:99]*incf;
Now we are going to consider three epochs and calculate the spectrum of each one of them
-->ndx = 1:100; // first epoch -->x1 = ecg(ndx); -->xsetech([0,0,1/2,1/3]) -->plot2d(t,x1) -->xsetech([1/2,0,1/2,1/3]) -->ffx1 = fft(x1,-1); -->plot2d(fre,abs(ffx1)^2 /100); -->ndx = 101:200; // second epoch -->x2 = ecg(ndx); -->xsetech([0,1/3,1/2,1/3]) -->plot2d(t,x2) -->xsetech([1/2,1/3,1/2,1/3]) -->ffx2 = fft(x2,-1); -->plot2d(fre,abs(ffx2)^2 /100); -->ndx = 201:300; // third epoch -->x3 = ecg(ndx); -->xsetech([0,2/3,1/2,1/3]) -->plot2d(t,x3) -->xsetech([1/2,2/3,1/2,1/3]) -->ffx3 = fft(x3,-1); -->plot2d(fre,abs(ffx3)^2 /100);This result shows that even although the signals are very similar, the power spectrum varies considerably. Notice that frequencies above 50 Hz are above Nyquist frequency and that the spectra are symmetrical
Now we are going to calculate the spectrum by using `pspect' with the next code
-->xbasc() -->xsetech([0,0,1,1]) -->Pxx = pspect(50,100,'re',ecg); -->plot2d(fre,Pxx);`Pspect' calculates the power spectrum of the whole signal after dividing it into epochs of 100 values. Notice that this spectrum is smoother than the individual spectra of the previous figure.
At this moment, we could ask whether the use of a longer epoch would have produced the same result. The next figure clearly shows that it is not the case. With the next code we are going to calculate the spectrum with epochs progressively longer
-->xbasc()
-->n = 200; // epoch of 200 values
-->fre = [0:n-1] * (1/(n*inct));
-->ffx1 = fft (ecg(1:n),-1);
-->xsetech([0,0,1,1/3])
-->plot2d(fre, (abs(ffx1)^2)/n);
-->xtitle('n = 200')
-->n = 500; // epoch of 500 values
-->fre = [0:n-1] * (1/(n*inct));
-->ffx2 = fft (ecg(1:n),-1);
-->xsetech([0,1/3,1,1/3])
-->plot2d(fre, (abs(ffx2)^2)/n);
-->xtitle('n = 500')
-->n = 1000; // epoch of 1000 values
-->fre = [0:n-1] * (1/(n*inct));
-->ffx3 = fft (ecg(1:n),-1);
-->xsetech([0,2/3,1,1/3])
-->plot2d(fre, (abs(ffx3)^2)/n);
-->xtitle('n = 1000')
Notice that, by increasing the length of the epoch, we increase the definition in frequency but we do not eliminate the variability among frequencies (increment of frequency becomes smaller). We need to average the epochs to get a clean, smooth power spectrum (alternatively we could average the values inside the spectrum).
![]() |