Finding Ellipse Axes

Question

Consider the following ellipse, generated by the bounding region of the following points

ps = {{-11, 5}, {-12, 4}, {-10, 4}, {-9, 5}, {-10, 6}};
rec = N@BoundingRegion[ps, "FastEllipse"];
Graphics[{rec, Red, Point@ps}]

We have that the ellipse 'rec' is given in the form

Ellipsoid[{-10.4, 4.8}, {{2.77333, 0.853333}, {0.853333, 1.49333}}]

How can I retrieve the lengths of the two main axes of such ellipsoid? Following this representation and Mathematica's general definition of Ellipsoid

I tried the following, using the eigenvalues of rec[[2]]

eigs = Eigenvectors[Inverse[rec[[2]]]]
eigv = Eigenvalues[Inverse[rec[[2]]]];
lens = 2/Sqrt[eigv]

Out[]= {2.06559, 3.57771}

where the 2 factor comes from the fact what what I retrieve from the eigenvalues is actually half the length of the main axis. Indeed we get

Graphics[{rec, Red, Point@ps,
  Blue, Line[
   RegionCentroid@rec + # & /@ {-(lens[[1]] eigs[[1]])/
      2, (lens[[1]] eigs[[1]])/2}],
  Line[RegionCentroid@rec + # & /@ {-(lens[[2]] eigs[[2]])/
      2, (lens[[2]] eigs[[2]])/2}]}]

Is this correct? Is there a quicker way of doing this?

Answer 1

Stolen from @J.M.'s answer, https://mathematica.stackexchange.com/a/239797/4999, with one correction (is that enough to make it not a duplicate?):

Nodes = ps;
ellipsoidBR = BoundingRegion[Nodes, "FastEllipse"]; (* not "FastEllipsoid" *)
center = ellipsoidBR[[1]];
{vals, vecs} = Eigensystem[ellipsoidBR[[2]]];
{a, b} = Sqrt[vals];
major = N@{center - a vecs[[1]], center + a vecs[[1]]}
minor = N@{center - b vecs[[2]], center + b vecs[[2]]}
Graphics[{ellipsoidBR, Red, Point@ps, Green, Point@center, 
  Line@{major, minor}}, Frame -> True]

Answer 2

Here's the alternate way using SingularValueDecomposition:

pt = rec[[1]]; mat = rec[[2]];
{u, s, v} = SingularValueDecomposition[(mat + Transpose[mat])/2]; 
func = Composition[AffineTransform[{u, pt}], ScalingTransform[Sqrt[Diagonal[s]]]];

and the length:

EuclideanDistance @@@ {func@{{0, -1}, {0, 1}}, func@{{-1, 0}, {1, 0}}}

{2.06559, 3.57771}

Graphics[{rec, Blue, Line[func@{{-1, 0}, {1, 0}}], 
  Line[ func@{{0, -1}, {0, 1}}]}]

Answer 3

I think the help of Ellipsoid is incomplete because it does not explain the input: Ellipsoid[p,[CapitalSigma]], where the second argument is called the "weight matrix".

You will remember that an ellipse (for simplicity I am explaining the 2D case, nD is similar and assume the ellipse is centered at the origin) can be written by:

x^2/rx^2 + y^2/ry^2 == 1

We may write this as:

{x,y}.{{1/rx^2,0},{0,1/ry^2}}.{x,y} == r.mat.r == 1

Note that the inverse Sqrt of the eigenvalues of m0 are the half axes of the ellipse.

If we now rotate the coordinate system by a rotation matrix: rot (r'=rot.r where r' are the new coordinates) the ellipse will be rotated (in the inverse sense) in the new coordinates:

r'. Transpose[rot].mat. rot .r' == r' . mat' . r' == 1

Therefore a rotated ellipse may be represented by a an symmetric (positive definite) matrix: mat'. This is called "weigh matrix" in the help.

Note that the eigenvalues of the matrix are not changed by a rotation. The eigenvectors point in the directions of the half axes.

Here is an example: Let

rx=2;
ry=1;
m0=DiagonalMatrix[{1/rx^1,1/ry^2}]
{x,y}.m0.{x,y}==1

This represents an axis aligned ellipse with half axes rx and ry:

Region[ImplicitRegion[{x, y} . m0 . {x, y} == 1, {x, y}], 
 Axes -> True]

If we now rotate the matrix m0:

rot = RotationMatrix[-Pi/4];
m= Transpose[rot].m0.rot;

we get a rotated ellipse:

Region[ImplicitRegion[{x, y} . m . {x, y} == 1, {x, y}], Axes -> True]

Threrfore, the half axes are obtained by the inverse Sqrt of the eigenvalues of the weight matrix. And the directions of the axes are given by the eigenvectors.