Finding the coordinates of points from distance matrix

Question

I have a set of points (with unknow coordinates) and the distance matrix. I need to find the coordinates of these points in order to plot them and show the solution of my algorithm.

I can set one of these points in the coordinate (0,0) to simpify, and find the others. Can anyone tell me if it's possible to find the coordinates of the other points, and if yes, how?

Thanks in advance!

EDIT Forgot to say that I need the coordinates on x-y only

Buttons

That's... going to need a lot of brute-forcing... – Lactone

Consider three points (a triangle). There are two orientations, and an infinite number of rotations that would give the same distance matrix. – Havener

One step further, are we talking a one-dimensional space, or two, or three, or four.... The answer will change in each case. By (0,0), should we accept its two-dimensional? – Fervency

Once you fix an orientation and angle, aren't all points after the third unambiguous? (In 2D space, since the OP gave a 2D zero.) – Stater

@KevinReid In general, yes. But if the first three points lie on one straight line, you have two options (obtained from each other by reflection in the line) until you chose the first point that doesn't lie on that line. – Prebo

Answer 1

Step 1, arbitrarily assign one point P1 as (0,0).

Step 2, arbitrarily assign one point P2 along the positive x axis. (0, Dp1p2)

Step 3, find a point P3 such that

Dp1p2 ~= Dp1p3+Dp2p3
Dp1p3 ~= Dp1p2+Dp2p3
Dp2p3 ~= Dp1p3+Dp1p2

and set that point in the "positive" y domain (if it meets any of these criteria, the point should be placed on the P1P2 axis).
Use the cosine law to determine the distance:

cos (A) = (Dp1p2^2 + Dp1p3^2 - Dp2p3^2)/(2*Dp1p2* Dp1p3)
P3 = (Dp1p3 * cos (A), Dp1p3 * sin(A))

You have now successfully built an orthonormal space and placed three points in that space.

Step 4: To determine all the other points, repeat step 3, to give you a tentative y coordinate. (Xn, Yn).
Compare the distance {(Xn, Yn), (X3, Y3)} to Dp3pn in your matrix. If it is identical, you have successfully identified the coordinate for point n. Otherwise, the point n is at (Xn, -Yn).

Note there is an alternative to step 4, but it is too much math for a Saturday afternoon

Fervency

@BrunoBruck the cosine law give the angle (the first equation) between P1P2 and P1P3. The next part is to get the projection of P3 onto the P1P2 axis. Knowing the distance P1P3, and setting it as the hypotenuse of a triangle, the X and Y values are simply the cos and sine times the hypotenuse, respectively. – Fervency

What you did with P2 is ok, but in the case of P3 I can't select a point of my set, which is not in the same line of P1 and P2, and say for sure that it's on y-axis. – Buttons

Ok, I think I got it. At first we pretend that P3 is in the y-axis to get a right triangle, and in such case we can create the equations for the coordinates. But we know the real distance between P3 and P2, so we are able to get the real angle between P1P2 and P1P3, and using it in the equations for the coordinates we can get the real values for Xp3 and Yp3. Did I understand right? – Buttons

@BrunoBruck sorry for the confusion, P3 shouldn't be on the y-axis, rather it is in the half-plane where y is positive. In other words, build the second dimension so that the third point has a positive y-value. I modified the text to that effect. – Fervency

Oh, now it makes sense and solves my problem. Thank you for the answer! – Buttons

Answer 2

The answers based on angles are cumbersome to implement and can't be easily generalized to data in higher dimensions. A better approach is that mentioned in my and WimC's answers here: given the distance matrix D(i, j), define

M(i, j) = 0.5*(D(1, j)^2 + D(i, 1)^2 - D(i, j)^2)

which should be a positive semi-definite matrix with rank equal to the minimal Euclidean dimension k in which the points can be embedded. The coordinates of the points can then be obtained from the k eigenvectors v(i) of M corresponding to non-zero eigenvalues q(i): place the vectors sqrt(q(i))*v(i) as columns in an n x k matrix X; then each row of X is a point. In other words, sqrt(q(i))*v(i) gives the ith component of all of the points.

The eigenvalues and eigenvectors of a matrix can be obtained easily in most programming languages (e.g., using GSL in C/C++, using the built-in function eig in Matlab, using Numpy in Python, etc.)

Note that this particular method always places the first point at the origin, but any rotation, reflection, or translation of the points will also satisfy the original distance matrix.

Mandell

This should be the answer. There is no need to code it yourself though, Multi-Dimensional Scaling functions can be found in Python or R. – Brocky

I met the same problem recently, could you please help me to see if my derivation process is correct? The link: math.stackexchange.com/questions/3663043/… – Recession

I am also unable to follow up, could you please help me see if I made it correct? math.stackexchange.com/questions/4566927/…, thank you – Bainter

@SamWashington, I commented on that question. Make sure you get your JS syntax right -- in particular, note that the power operator in some languages (including Javascript) is ** not ^. – Mandell

Answer 3

This is a math problem. To derive coordinate matrix X only given by its distance matrix.

However there is an efficient solution to this -- Multidimensional Scaling, that do some linear algebra. Simply put, it requires a pairwise Euclidean distance matrix D, and the output is the estimated coordinate Y (perhaps rotated), which is a proximation to X. For programming reason, just use SciKit.manifold.MDS in Python.

Caswell

Answer 4

If for points p, q, and r you have pq, qr, and rp in your matrix, you have a triangle.

Wherever you have a triangle in your matrix you can compute one of two solutions for that triangle (independent of a euclidean transform of the triangle on the plane). That is, for each triangle you compute, it's mirror image is also a triangle that satisfies the distance constraints on p, q, and r. The fact that there are two solutions even for a triangle leads to the chirality problem: You have to choose the chirality (orientation) of each triangle, and not all choices may lead to a feasible solution to the problem.

Nevertheless, I have some suggestions. If the number entries is small, consider using simulated annealing. You could incorporate chirality into the annealing step. This will be slow for large systems, and it may not converge to a perfect solution, but for some problems it's the best you and do.

The second suggestion will not give you a perfect solution, but it will distribute the error: the method of least squares. In your case the objective function will be the error between the distances in your matrix, and actual distances between your points.

Kelwunn

Thank you for the answer. I don't know if this is the best approach because in some cenarios I have a lot o points and a metaheuristic does not always return the optimal solution, or in this case, a feasible solution. So I could spend a lot of time with it and still doesn't get a feasible answer. – Buttons

@DeepYellow: I like your answer in part because it may help to answer a different, harder question, posted by a different user yesterday. I tried to answer that other question and failed. If the challenge interests you, here is the URL: #10957859 – Kenyatta

@thb: Thanks for pointing that question out. I posted what I think is a correct solution, let me know what you think. – Kelwunn

Answer 5

The "eigenvector" method given by the favourite replies above is very general and automatically outputs a set of coordinates as the OP requested, however I noticed that that algorithm does not even ask for a desired orientation (rotation angle) for the frame of the output points, the algorithm chooses that orientation all by itself! People who use it might want to know at what angle the frame will be tipped before hand so I found an equation which gives the answer for the case of up to three input points, however I have not had time to generalize it to n-points and hope someone will do that and add it to this discussion. Here are the three angles the output sides will form with the x-axis as a function of the input side lengths:

angle side a = arcsin(sqrt(((c+b+a)*(c+b-a)*(c-b+a)*(-c+b+a)*(c^2-b^2)^2)/(a^4*((c^2+b^2-a^2)^2+(c^2-b^2)^2))))*180/Pi/2

angle side b = arcsin(sqrt(((c+b+a)*(c+b-a)*(c-b+a)*(-c+b+a)*(c^2+b^2-a^2)^2)/(4*b^4*((c^2+b^2-a^2)^2+(c^2-b^2)^2))))*180/Pi/2

angle side c = arcsin(sqrt(((c+b+a)*(c+b-a)*(c-b+a)*(-c+b+a)*(c^2+b^2-a^2)^2)/(4*c^4*((c^2+b^2-a^2)^2+(c^2-b^2)^2))))*180/Pi/2

Those equations also lead directly to a solution to the OP's problem of finding the coordinates for each point because: the side lengths are already given from the OP as the input, and my equations give the slope of each side versus the x-axis of the solution, thus revealing the vector for each side of the polygon answer, and summing those sides through vector addition up to a desired vertex will produce the coordinate of that vertex. So if anyone can extend my angle equations to handling beyond three input lengths (but I note: that might be impossible?), it might be a very fast way to the general solution of the OP's question, since slow parts of the algorithms that people gave above like "least square fitting" or "matrix equation solving" might be avoidable.

Callboard