[Math] How to find the minimal axis-parallel ellipse enclosing a set of points.

conic sectionsgeometry

I have a set $X$ of points in $\mathbb{R}^2$ and I'm trying to find the smallest encompassing ellipse which main axes are parallel to the coordinate system's (to put it differently, its both centres share one coordinate). I need the gravitational centre and the vertical and horizontal radius.

Now, I managed to do this for a horizontally aligned rectangle, but that isn't much help (although it's a first approximation, as I can easily draw a rectangle around $X$). I also found some formulas on the internet, but they seem to be wrong as I keep getting $0$ for the radius.

Can this be done algebraically?

Best Answer

There is a very elegant algorithm to find an encompassing ellipse for points arbitrarily positioned in space via PCA (Principal Component Analysis) approach to find axes of the ellipse.

Suppose you have $n$ points of an ellipse stored in the matrix $E$:

$$ E = \begin{bmatrix} x_1 & y_1 \\ x_2 & y_2 \\ \vdots & \vdots \\ x_n & y_n \end{bmatrix} $$

Let us also suppose that the center of the ellipse is in the origin $(0, 0)$ (if not, just subtract corresponding means from each of the coordinates, i.e. $x_i := x_i-\mu_x, y_i := y_i - \mu_y $).

Then, compute the covariance matrix $C$, in our case we can use the following formula as the resulting matrix will be proportional to the real covariance matrix:

$$ C = E^TE $$

Finally, find eigenvectors $v_i$ of $C$, in our particular case, $i = 1, 2$. The resulting eigenvectors will have the same direction as the axes of the ellipse you work on. Also, if you normalize a matrix $V = [v_1\;v_2]$, you may see it is actually a rotation matrix. So, if you want to rotate your ellipse in a way that all its points are parallel to coordinate axes, just do the multiplication: $EV$.

In order to find magnitudes of ellipse axes, use the rotated ellipse $EV$: find the biggest values of $x$ and $y$ in the set of points, and these values are the ones you need.

This algorithm will definitely work if the points lie in a 3D space, and it does not require points to lie exactly on the ellipse.

Sources on PCA, eigenvalues, eigenvectors:

  1. http://www.doc.ic.ac.uk/~dfg/ProbabilisticInference/IDAPILecture15.pdf
  2. http://en.wikipedia.org/wiki/Eigenvalues_and_eigenvectors