Java SVD with JAMA or else

Question

I have a cloud of Points and I need the best fitting Line. I'm using JAMA but I don't know why, something is not working. Probably it's me who doesn't get how it works. I have a Nx3 Matrix (this is what JAMA svd supports) and I get the right Matrix V from the Svd. The vector I need is the right singular vector (row or coloumn ? ) corresponding to the greatest singular value. This vector is supposed to represent the main direction.

Every vector is supposed to have y as the greatest positive value, x can be both positive or negative and z is supposed to be negative. But sometimes the vector I get has a negative y value, or anyway it's pointing in the wrong direction.

My cloud of Points are pretty regular, the are all postioned almost along the y axis (with z small and negative). So the main direction is supposed to be really easy to find. But it's still not working properly.

In this case I'm getting the row vector (I tried the coloumn vector too) of the right Matrix V. I already substracted the centroid from "pointSet".

public static Matrix bestDirection(Matrix pointSet){

    Matrix bestFittingLine = new Matrix(3,1);
    SingularValueDecomposition svd = pointSet.svd();

    bestFittingLine.set(0, 0, svd.getV().get(0, 0));
    bestFittingLine.set(1, 0, svd.getV().get(0, 1));
    bestFittingLine.set(2, 0, svd.getV().get(0, 2));

    return bestFittingLine;
}

I guess maybe I'm not considering something. Idk maybe I should use another technique or another library.

Answer 1

From Wikipedia on SVD:

Non-degenerate singular values always have unique left and right singular vectors, up to multiplication by a unit phase factor (for the real case up to sign).

Simply put, you can't rely on the sign of the output singular vectors.

You may also need to center the data before SVD.

Why won't you perform a regression?

Answer 2

If your equation takes this form:

z = a0 + a1*x + a2*y

Your matrix equation looks like this for N points:

z(i) = a0 + a1*x(i) + a2*y(i)  i = 1, N

The left hand side is an Nx1 vector; the right hand side is an Nx3 matrix multiplying an unknown vector that's 3x1.

Multiply both sides by A(transpose) and you end up with a 3x3 matrix multiplying a 3x1 vector of unknown coefficients that's equal to a 3x1 vector. Use standard matrix solutions to solve for the unknown coefficients. It'd be easy to do, even in closed form.

That's the linear least squares solution simplified. Here's a scribd document that spells it out in more detail.