Cconditional subtraction of 2 vectors of different lengths

I have two very large vectors 412432x1 and 1293843x1 of 2 different signals (the vectors are the arrival times of the signals) that I need to excute only if they satisfy my condition.

Let's say for example I have Vectors A = 5x1 and B = 10x1

and a constant x = number.

All elements of A should be compared with all that of B and if A - B < x a new matrix C should come out that contains the correct elements of A and B which satisfies this condition.

The thing is whenever I use normal subtraction A - B' my matlab stops the code because of the memory, therefore I want to find out the right pairs to subtract, and then execute.

Thank you for your help!

[Dyuman Joshi]

A-B' requires around 3976 GBs of memory, which, I doubt, your PC has that amount of RAM to store the data.

You can run a single (requires storing intermiediate data) through elements of one of the variables or double loop (directly storing the required data) through elements of both the variables, checking the condition and storing data accordingly.

However, it is difficult to pre-allocate in this particular case. And your computer might not have enough RAM to store the resulting data.

[the cyclist]

How frequently do you expect the condition to be satisfied? If it is rare, then you will be able to store the result in a sparse matrix.

As it is, you are also doing a somewhat large number of computations, but this is mainly a memory problem. Looping over the vectors in the dumbest way possible (and doing the condition check and an addition), took 25 minutes on my machine. (This could be made faster by vectorizing one of the loops, but I wanted to see if the dumb way worked.)

a = rand(1,412432);
b = rand(1,1293843);

x = 0.5;

tic
for ii = 1:412432
    if (mod(ii,10000)==0)
        fprintf("iteration %d of 412,432\n",ii)
    end
    for jj = 1:1293843
        if (a(ii) - b(jj)) < x
            c = a(ii)+b(jj); % Note that I am not storing the result. Just looking at timing.
        end
    end
end
toc

One other thing to think about is what is your next step with this (possibly huge) result? Do you really need all of these data points? Maybe some kind of sampling will be adequate? Often when people think they need to solve these huge problems, they really don't.

[Image Analyst]

How exactly do you plan on doing A-B when they have different number of elements, plus A - B' is even worse (see comments above mine)? If A has 5 elements and B has 10 elements, you can subtract the first 5, but what will B(6) be subtracted from??? There is no corresponding A(6) so what are your plans/hopes/dreams there?

If you have any more questions, then attach your data (it's actually not that large), and code to read it in with the paperclip icon after you read this:

TUTORIAL: How to ask a question (on Answers) and get a fast answer

Answer 1

[David Goodmanson]

Hi salah,

You can do quite a bit if the two vectors are sorted. Assume you have two sorted vectors sx and sy of length Nx and Ny. The code below creates an (Nx) x 3 matrix 'sdnupy' with the following properties: The first column is a running row index nx for the sx vector. Columns 2 and 3 are an upper and lower index ndn,nup such that for each nx, sy(ndn:nup) are all the sy points within distance d of sx(nx). There is also an (Ny) x 3 matrix 'sdnupx' such that for each ny, sx(ndn:nup) are all the sx points within distance d of sy(ny). All the required information is contained in each of these matrices.

If x and y have to remain unsorted, then for each nx, the y points are not going to be ndn:nup but could be scattered all over the place and have to be specified individually. If there are a large number of points for each nx, you will get into memory problems. The last part shows for an unsorted situation how to at least obtain all the unsorted y points for a single given x point (and all the unsorted x points for a single given y point).

If the vectors are already sorted you can remove the lines involving sort and invperm and change sx and sy to x and y further on down.

For Nx = 4.2e5, Ny = 1.3e6 this takes less than a second.

y = rand(1,1.3e6);
x = rand(1,4.2e5);
d = .1;

x = x(:);
y = y(:);
Nx = size(x,1);
Ny = size(y,1);
nx = (1:Nx)';
ny = (1:Ny)';
[sx si2ix]  = sort(x);
[sy si2iy]  = sort(y);
i2six = invperm(si2ix);        
i2siy = invperm(si2iy);
% si2ix(sorted index) = unsorted index
% i2six(unsorted index) = sorted index

% find max and min sorted y indices for each sorted x index  
nup = floor(interp1(sy,ny,sx+d));
ndn =  ceil(interp1(sy,ny,sx-d));
fixup = ~isnan(ndn)&isnan(nup);  
fixdn = isnan(ndn)&~isnan(nup);  
nup(fixup) = Ny;
ndn(fixdn) = 1;
sdnupy = [nx ndn nup];

% find max and min sorted x indices for each sorted y index
nup = floor(interp1(sx,nx,sy+d));
ndn =  ceil(interp1(sx,nx,sy-d));
fixup = ~isnan(ndn)&isnan(nup);  
fixdn = isnan(ndn)&~isnan(nup);  
nup(fixup) = Nx;
ndn(fixdn) = 1;
sdnupx = [ny ndn nup];

% unsort if necessary for a single x or y index 

a = 1e5;                            % specify an x index
sa = i2six(a);                      % sorted x index
x(a)                           
sx(sa)                              % same
sindy = sdnupy(sa,2):sdnupy(sa,3);  % max and min sorted y index
indy = si2iy(sindy);                % unsorted y indices
max(abs(y(indy)-x(a)))              % check, should be <= d

b = 2e5;                            % specify a y index
sb = i2siy(b);                     
y(b)
sy(sb)
sindx = sdnupx(sb,2):sdnupx(sb,3);  
indx = si2ix(sindx);                % unsorted x indices  
max(abs(x(indx)-y(b)))              % check, should be <=d


function b = invperm(a)

% inverse permutation for a row or column vector.
% for permutation a, create b such that a(b) = 1:n 
% in which case also b(a) = 1:n

[~, b] = sort(a);