How can I add zeros to a vector to make a continuous time axis from a discontinuous one?

I have need to create a continuous time axis in milliseconds to by adding zeros to a vector where data is missing. If you load the file, you can see

time = 
'3225 sec'
'3226 sec'
'3879.9 sec'
'3880.9 sec'

dat =

  4×3 timetable

    Record Time         ECG1               ECG2             ECG1_ECG2   
    ___________    _______________    _______________    _______________

    3225 sec       {1000×1 double}    {1000×1 double}    {1000×1 double}
    3226 sec       {1000×1 double}    {1000×1 double}    {1000×1 double}
    3879.9 sec     {1000×1 double}    {1000×1 double}    {1000×1 double}
    3880.9 sec     {1000×1 double}    {1000×1 double}    {1000×1 double}

Data was recorded at 1kHz but some data was lost. In the end, I need to save the EKG1 variable as a double in milliseconds. I am trying to upsample the time vector to milliseconds and add zeros to the data axis where it's missing

load('sample.mat');
EKG = cat(2,cell2mat(dat{:,1}),cell2mat(dat{:,2}),cell2mat(dat{:,3})); % Convert timetable to matrix
EKG1 = [EKG(1:2000,:); zeros(653900,3); EKG(2001:end,:)]; % Add some zeros to the middle manually

recorded_times = dat.Properties.RowTimes; % extract the record times

ms=milliseconds(recorded_times(1):recorded_times(end))'; % I am trying to convert the axis to milliseconds but it doesn't work

I can see the amount of missing time from recorded_times with diff()

diff(time)

ans = 

  3×1 duration array

        1 sec
   653.85 sec
        1 sec

This seems to work great except for being off by 0.05 seconds, but then

find(diff(time)>1)

ans =

  0×1 empty double column vector

doesn't work seemingly because this is a duration? Not sure. I just need to add 653.9 milliseconds worth of zeros to my EKG variable, which I have done manually but I need to repeat this process for multiple datasets, sometimes with multiple spots of missing data.

Answer 1

[Star Strider]

I‘m not certain what data are in the table, however the first column appears to be the length of the recordings. There do not appear to be any missing values in the two EKG records that I looked at, however since there is no accompanying time vector for the EKG traces, it isn’t possible to determine if any EKG data are missing. And wth no time vector, interpolation is impossible, since there are no gaps to fill.

If a time vector for each record were present., the Signal Processing Toolbox resample function would be the preferred interpolation method, since it includes an anti-aliasing filter. This is important if signal processing is the ultimate goal.

.

[Dustin7]

Sorry for the misunderstanding. I am analyzing EKG data that was sampled at 1 kHz and originally saved in EDF+D due to the discontinuities. Each row contains 1 seconds worth of data. I want to save the data using edfwrite() as EDF+C and therefore need to add zeros when the signal was lost.

The first column of the dat variable are timestamps of when the data was collected starting from zero (I removed the beginning to make the file smaller). From '3226 sec' to '3879.9 sec', there is 653.9 seconds of missing data and I need to add zeros to the EKG data vector to account for this.

This is my solution:

EKG1 = [EKG(1:2000,:); zeros(653900,3); EKG(2001:end,:)]; % Add some zeros to the middle manually

but I need to know how to get 653.9 seconds from the data instead of manually adding it.

[Star Strider]

O.K., none of that was obvious.

To just zero-pad the vector, do something like this —

EKG1 = randn(1,10);                                         % Data VEctor
Fs = 1000;                                                  % Sampling Frequency
t = linspace(0, numel(EKG1)-1, numel(EKG1))/Fs;             % Data Vector Time Vector
EKGv = zeros(1,25);                                         % Desired LEngth 'zeros' Vector
tv = linspace(0, numel(EKGv)-1, numel(EKGv))/Fs;            % Desired Length Time Vector
EKG1v = EKGv;                                               % Define Result Vector
EKG1v(1:numel(EKG1)) = EKG1;                                % Fill Result Vector With Data Vector

figure
subplot(2,1,1)
plot(t, EKG1)
grid
xlim([min(tv) max(tv)])
subplot(2,1,2)
plot(tv, EKG1v)
grid

So the only operation that needs to be done is to define the zeros vector to be the appropriate length. After that, just index into it. The data vector can go anywhhere within the ‘Desired Result’ vector, for example —

EKG1v = EKGv;                                               % Define Result Vector
EKG1v((1:numel(EKG1))+10) = EKG1;                           % Fill Result Vector With Data Vector (Here, Offset By 10 Elements)

figure
subplot(2,1,1)
plot(t, EKG1)
grid
xlim([min(tv) max(tv)])
subplot(2,1,2)
plot(tv, EKG1v)
grid

That random vector looks enough like an EKG that I’ll just go with it!

.

[Dustin7]

Thanks for the quick answers, generally this looks right but it doesn't work for my exact problem. Are you able to load the sample.mat file attached to my original post? I need to upsample the data from the 'time' variable. If you load my data and use this code, you will see missing components of the "time_ms" variable. How can I add something for each millisecond and add zeros to the data vector at the same locations?

Here is my very dirty, brute force method that achieves what it needs to. How can I make it work better and more flexibly?

q = seconds(time) * 1000 % Convert duration to double as milliseconds
t1=linspace(q(1),q(1)+999,999) % Add timepoints of the rest of the first second
t2=linspace(q(2),q(2)+999,999) % 999 remaining timepoints for the second second
t3=linspace(q(3),q(3)+999,999) % 999 remaining timepoints for the third second
t4=linspace(q(4),q(4)+999,999)
time_ms = [q(1), t1, q(2), t2, q(3), t3, q(4), t4]

it = find(diff(time_ms) > 1000) % Find where to begin adding to the time vector

time2 =[time_ms(1:it), linspace(time_ms(it),time_ms(it+1),time_ms(it+1)-time_ms(it)-1), time_ms(it:end)]
% This is close to correct for the time vector but it's very ugly. Then I
% could maybe use 'it' to add data to EKG1 but it seems so tedious.

EKG1 = [EKG(1:2000,:); zeros(time_ms(it+1)-time_ms(it)-1,3); EKG(2001:end,:)];

% How can I get to the below line in a more efficient and automated way?
plot(time2,EKG);

[Star Strider]

I was able to load it and plot it.

The problem is that the time vectors associated with those records are missing, so any sort of time-respective processing is simply impossible. The only option is to use the indices as the time vector and go from there. That’s horribly inaccurate, however it’s the only option.

.

[Dustin7]

Yes, I know the time vectors are missing. All I need to do is replace them with zeros so that the data will be continuous again because I need to open it in another EDF-based program that only accepts EDF+C. I updated my previous post with additional code.

[Star Strider]

The problem is that because the corresponding time vectors are not included, there is no record of which data are missing, so interpolating them is absolutley impossible. If the corresponding time vectors had been saved, this would be straightforward (instead of impossible). This is just math. There is no way around it, and wantinf it to be possible does not make it so.

[Dustin7]

OK, I can make the time vector easily. How can I more easily align the original data with the pieces where it is known and add zeros to the middle of 'time_cont' ? There should be a way I can index the new time vector using the old time vector.

time = [3225 3226 3879.9 3880.9] *1000 % Original time now in ms
time_cont = time(1):0.001:time(end) % Time vector sampled at 1 kHz

[Star Strider]

‘OK, I can make the time vector easily.‘

Sorry, but simply No. The time vector has to correspond to the times the samples were collected, and has to be created at that time the signal was sampled. Thats the whole idea in having it. If there are missing signal elements and there is no corresponding time vector that would indicate when they should have been collected, that informaiton is lost forever and any subsequent processing is essentially impossible.

Again, a trivial example —

x = 0:0.1:20;
y = sin(2*pi*x*0.1) + randn(size(x))/10;

figure
plot(x, y)
grid
title('Original Complete Signal & Time Vector')
text(7.5, 1.3, sprintf('Signal Length = %d', numel(y)))

Lv = randi([0 1], size(x))==1                                   % Logical Index Vector

Lv = 1×201 logical array
   1   0   0   0   1   1   1   1   1   0   1   0   0   0   0   1   0   0   0   1   0   1   1   0   0   1   1   1   1   0   0   0   0   1   1   1   0   0   1   0   0   1   1   0   0   0   1   0   1   1

figure
plot(y(Lv))
grid
title('Signal With Missing Values Without Time Vector')
text(40, 1.3, sprintf('Signal Length = %d', numel(y(Lv))))

Now, figure out what samples are missing in the second plotted signal, and fill them in.

.

[Dustin7]

You would be correct if random samples were missing, but this isn't the case. I know that the data is missing between 3226.999 and 3879.9 seconds and I can make a time axis based on knowing the sampling rate of 1 kHz. Matlab reads data in EDF+D format as in the sample.mat in my orignal post, i.e., without the full time vector, which needs to be created anyway.

[Star Strider]

It would certainly be possible to careat a time vector that spanned the entire time the EKG is supposed to be recorded, and then calculate the indices that were actually recorded and plot them against the full time vector, however filling the missing information would be impossible.

Another trivial example —

t1 = [1:0.1:10 20:0.1:30];
y1 = sin(2*pi*t1*0.1) + randn(size(t1))/10;

figure
plot(y1)
grid
title('Original y_1 Without Time Vector')

t2 = linspace(min(t1), max(t1), t1(end)/((t1(2)-t1(1))));
y2 = interp1(t1, y1, t2);

figure
plot(t2, y2)
grid
title('Interpolated Signal To Fill Missing Region')

.