How to generate M *almost* mutually orthogonal vectors of dimensionality N, where M>>N ?

Generating N mutually orthogonal vectors in is easy enough. I am looking for M vectors in , where , and the M vectors are almost mutually orthogonal. With large enough N one may find many such vectors, subject to some upper bounded small deviation from orthogonality. The idea not necessarily to get all the pairs of vectors have the same dot product, but all pairs have as small as possible overlap (deviation of the dot product from zero), with some small upper bound on the worst case deviation (left for the algorithm to determine). There are blogs and references to papers in MathOverflow on methods, upper bounds, and such, but before I plunge into implementing anything, I am looking to see if there is an implementation or any suggestions for solutions with Matlab.

For instance, I am looking to generate (any) set of 4096 almost orthogonal vectors in .

Thanks.

[Shlomo Geva]

These figure may assist -

In the first figure I select 48 random unit vectors in and plot the histogram of all mutual dot-products.

In the second figure I choose 2 sets of 24 vectors in , and concatenate to form a set of 48 vectors. Each set of 24 vectors is a randomly oriented orthonormal basis set in . I then plot the histogram of all mutual dot-products.

The second set has a larger area under the curve in the middle - meaning that most pairs are closer to orthognality than in the first set. Alas, it also has a few outliers with larger overlap. 'No free lunch' comes to mind. Obviously, as we increase the number of near orthonormal vectors we wish to pack into , the more difficult it becomes to find such vectors. But it seems that the second approach is still on average superior, but does have a cost in having outliers with higher overlap.

(note, the histogram peak at 1is self dot product)

I apologise - this is a real doozy, or so it seems. The problem is not new of course, but I am seeking an implementation before attempting to reinvent the wheel.

[Matt J]

The problem is not new of course, but I am seeking an implementation before attempting to reinvent the wheel.

It seems to me though that you haven't defined an acceptability criteria. It has already been proven to you mathematicaly that you cannot make all overlaps smaller than 1/M simultaneously. What then is "small enough" for your purposes?

[Shlomo Geva]

The acceptability criteria is - as close as possible to the theoretical limit.

Here is the solution in your earlier comment Matt (I missed it, sorry).

Here is what it looks like, as in the figures above (for 48 vectors in :

And here with 4096 vectors in

[Shlomo Geva]

@Matt J the vectors were normalised to unit vectors in the plots above.

Answer 1

[Matt J]

One way to pose that mathematically is to solve the minimization problem,

where is the MxM identity matrix. Here the columns of A are the M vectors that you are looking for. If is the s.v.d. of A, this can be rewritten,

where is the space of vectors in with at most N non-zeros. It is clear, therefore, that the objective function f() can never get smaller than .

A consequence of this is that if the vectors are all required to be unit vectors, there will always be at least two distinct vectors in the set whose dot product is at least .

[Shlomo Geva]

Thanks!

[Matt J]

Note also that an optimally orthogonal solution, under the above criteria can be determined as

[U,~,V]=svd( rand(24,4096) , 'econ' );

A=U*V.';

[Shlomo Geva]

That's a pretty good solution, thanks!

[Matt J]

Note, however, that this solution does not produce unit vectors. You can post-normalize each column of A if you wish, but that will take you further (considerably) from the unconstrained cost function minimum.

[Shlomo Geva]

Yes, I noticed that and normalised the vectors. Normalization does not change the angle between vectors and so it does not matter in terms of orthogonality. It is necessary to normalise when comparing solutions.

[Matt J]

Normalization does not change the angle between vectors and so it does not matter in terms of orthogonality.

It matters. The unconstrained minimization is only concerned about making the majority of the dot products small, but remember that just because the dot product of two vectors is close to zero doesn't mean that their separation angle is close to 90 degrees. Think of the two vectors u=[1,0,0] and v=1e-20*u, whose dot product is the small number 1e-20 but which are perfectly parallel to one another.

You could use fminunc() to minimize subject to an orthonormality condition, as below. It might run better, though, if the SpecifyObjectiveGradient option were used.

[N,M]=deal(24,4096);

[U,~,V]=svd( rand(N,M) , 'econ' );
A0=U*V.';

A=fminunc(@objective,A0);

function val=objective(A)

  [N,M]=size(A);
  A=normalize(A,,1,'n');
  val=norm(A'*A-speye(M),'fro');

end

[Shlomo Geva]

Thanks.

Perhaps I should have been clearer in the problem defnition - looking for unit vectors, i.e near orthonormal, hence the dot product value is constrained to [+1...-1] and only depends on the angle.

I tried the above code, but it went into a black hole (no problem with computational resources - optimisation making no progress though). Have not used fminunc before, so so sure how to use the SpecifyObjectiveGradient option for, perhaps, a more productive execution.

PS: finding M binary vectors b in is of particular interest, and with {+1,-1} encoding their lengths are obviously equal, but that's a somewhat different problem. The vectors are then constrained to lie on the corners of a hypercube (within a unit hypesphere.) At the moment just taking the sign() of the continuous solution. Not sure if there's anything better that can be done.

Answer 2

[John D'Errico]

I'm sorry, but your goal is difficult, virtually impossible. It depends on how close they must be to "orthogonal". In fact, in only 24 dimensions, if you take any 4096 vectors in 24 dimensions, one of them will actually be quite close to at least ONE of the other vectors. Nearly mutually orthogonal will be almost impossible, although you could think about how to make them maximally different. First though, you need to understand what is the issue.

You can view this as a dense packing problem of sorts. To make it easy, think in a LOW number of dimensions. 2-d or 3-d are nice things, because we can visualize what is happening.

Start with 2-d. Suppose I asked you to find a set of two vectors that are maximally orthogonal to each other in the plane? That is easy. In fact, we can visualize it in terms of finding two points on a unit half circle that are 90 degrees apart. This is because a vector will define a line through the origin, and where that line crosses the unit circle will also define the vector. So there is a one-to one correspondence between the set of possible vectors we care about, and the set of points on the perimeter of a circle. And since a line passes through a unit circle in two places, diametrically opposite, then we need only think about the semi-circle in the first and second quadrants.

For example...

% first, I'll draw a semi-circle
theta = linspace(0,pi);
x = cos(theta),y = sin(theta);

x = 1×100

1.0000 0.9995 0.9980 0.9955 0.9920 0.9874 0.9819 0.9754 0.9679 0.9595 0.9501 0.9397 0.9284 0.9161 0.9029 0.8888 0.8738 0.8580 0.8413 0.8237 0.8053 0.7861 0.7660 0.7453 0.7237 0.7015 0.6785 0.6549 0.6306 0.6056

% pick one random "vector" in the first quadrant
t1 = rand(1)*pi/2;
% and a second vector that is exactly 90 degrees away from it.
t2 = t1 + pi/2;
v1 = [cos(t1),sin(t1)];
v2 = [cos(t2),sin(t2)];
plot(x,y,'g-',[0 0;v1(1),v2(1)],[0 0;v1(2),v2(2)],'r-')
axis equal
yline(0)

That was easy, right? So we can think of two vectors in the plane, simply by thinking about where they intersect the unit upper semi-circle. But suppose now I asked you to create 4096 such vectors in the plane? Now we can think of 4096 points on the unit semi-circle. Would you agree the correspondence still works? But would you also agree that any of those vectors will actually be EXTREMELY close to its immediate neighbors on the unit semi-circle?

You simple could not ask to find 4096 nearly orthogonal vectors in the plane. In fact, they are all incredibly close to at least one other vector. Again, we can view this a a packing problem on the unit circle. Find 4096 points on the unit circle that are as far apart from each other as possible. In 2-d, this is easy.

T = linspace(0,pi,4097);
T(end) = []; % throw away the last one, because it is the same as the first.
X = cos(T); Y = sin(T);
plot(x,y,'o')
axis equal

Interestingly, it looks like the plot did not show 4096 points, since they would be packed so denssely, as to become visually indistinguisable. Regardless, you get the idea.

If we wanted to find vectors that are maximally orthogonal to each other, then we want to pack points on the surface of a unit sphere, so they are maximally distant from each other.

Therefore, finding 4096 maximally orthogonal vectors as you wish to do is the fundamental equivalent to packing 4096 maximally equidistant points on the surface of a sphere in 24 dimensions. Yes, a 24 dimensional hyper-sphere is a big place to fit a lot of things. But not that big. Effectively, you are asking to find 4096 points on the surface of that unit 24 dimensional hyper-sphere. They simply will not be very orthogonal. What is worse though, is the problem of sphere packing in higher numbers of dimensions is a difficult one to solve.

You may almost be in luck, in a sense.

https://www.quantamagazine.org/sphere-packing-solved-in-higher-dimensions-20160330/

In that link, we see that a Ukrainian mathematician has solved the sphere packing problem in specifically 8 and 24 dimensions. Does that effort apply to packing on the surface of a sphere? Sadly, probably not, since there you effectively need to pack in 23 dimensions.

Ok. How might I try to most simply solve your problem, as closely as I can, so that I can find n-vectors that are maximally orthogonal? I would start with a large number of randomly scattered points on the surface of the unit hyper-sphere. Then find the closest pair of points, and discard one of them, arbitrarily. Repeat until you have only 4096 points remaining.