2
$\begingroup$

I see here and there that there is a closed form expression of the least-square fit of a sphere of radius $r$ and center $\mathbf{c}$ to $N$ data points $\{\mathbf{x_i}\}_{i\in(1,\cdots,N)}$.

How is this expression obtained?


I considered this energy to be minimized: $$ E(r,\mathbf{c}) = \sum_{i=1}^N (\|\mathbf{x}_i-\mathbf{c}\|^2 - r^2)^2, $$ then I tried to find radius $r^*$ and center $c^*$ that minimize $E$, so that it makes its gradient null. I obtained: $$ \frac{\partial E}{\partial r}(r^*,\mathbf{c}^*) = 0 \ \Rightarrow \ -4\sum_{i=1}^N (\|\mathbf{x}_i - \mathbf{c}^*\|^2 - {r^*}^2) r^* = 0 \ \Rightarrow\ {r^*}^2 = \frac{1}{N}\sum_{i=1}^N \|\mathbf{x}_i - \mathbf{c}^*\|^2 , $$ which match what I found on the cited references above.

But I couldn't find any correct way to obtain the system of equation $A\mathbf{c} = \mathbf{b}$ in order to find the center $\mathbf{c}^* = (A^T A)^{-1}A^T\mathbf{b}$.

Is it by expressing the partial deivative of $E$ with respect to $\mathbf{c}$? $$ \frac{\partial E}{\partial \mathbf{c}} (r^*,\mathbf{c}^*) = 0 \ \Rightarrow\ -4\sum_{i=1}^N(\|\mathbf{x}_i-\mathbf{c}^*\|^2 - {r^*}^2)(\mathbf{x}_i - \mathbf{c}^*) = 0. $$ I tried to develop further this equation by replacing $r^*$ by its expression above but I couldn't put it in the form $A\mathbf{c}^* = \mathbf{b}$

$\endgroup$
2
+100
$\begingroup$

$\def\vec{\boldsymbol}\def\peq{\mathrel{\phantom{=}}{}}$Because\begin{align*} \frac{\partial E}{\partial \vec{c}}(r, \vec{c}) &= -4\sum_{k = 1}^n (\|\vec{x}_k - \vec{c}\|^2 - r^2)(\vec{x}_k - \vec{c})\\ &= 4\left( \sum_{k = 1}^n (\|\vec{x}_k - \vec{c}\|^2 - r^2) \right) \vec{c} - 4\sum_{k = 1}^n (\|\vec{x}_k - \vec{c}\|^2 - r^2) \vec{x}_k\\ &= 4\left( \color{blue}{\sum_{k = 1}^n \|\vec{x}_k - \vec{c}\|^2 - nr^2} \right) \vec{c} + 4r^2 \sum_{k = 1}^n \vec{x}_k - 4\sum_{k = 1}^n \|\vec{x}_k - \vec{c}\|^2 \vec{x}_k, \end{align*} Plugging ${r^*}^2 = \dfrac{1}{n} \sum\limits_{k = 1}^n \|\vec{x}_k - \vec{c}^*\|^2$ into $\dfrac{\partial E}{\partial \vec{c}}(r^*, \vec{c}^*) = \vec{0}$ yields\begin{gather*} \sum_{k = 1}^n \|\vec{x}_k - \vec{c}^*\|^2 \vec{x}_k - \frac{1}{n} \left( \sum_{k = 1}^n \|\vec{x}_k - \vec{c}^*\|^2 \right) \left( \sum_{k = 1}^n \vec{x}_k \right) = \vec{0}.\tag{1} \end{gather*} Since $\|\vec{x}_k - \vec{c}^*\|^2 = \|\vec{x}_k\|^2 + \|\vec{c}^*\|^2 - 2\langle\vec{x}_k, \vec{c}^*\rangle$ for all $1 \leqslant k \leqslant n$, then\begin{gather*} \sum_{k = 1}^n \|\vec{x}_k - \vec{c}^*\|^2 \vec{x}_k = \sum_{k = 1}^n (\|\vec{x}_k\|^2 + \|\vec{c}^*\|^2 - 2\langle\vec{x}_k, \vec{c}^*\rangle) \vec{x}_k\\ = \sum_{k = 1}^n \|\vec{x}_k\|^2 \vec{x}_k + \|\vec{c}\|^2 \sum_{k = 1}^n \vec{x}_k - 2\sum_{k = 1}^n \langle\vec{x}_k, \vec{c}\rangle \vec{x}_k, \end{gather*}\begin{gather*} \left( \sum_{k = 1}^n \|\vec{x}_k - \vec{c}^*\|^2 \right) \left( \sum_{k = 1}^n \vec{x}_k \right) = \left( \sum_{k = 1}^n (\|\vec{x}_k\|^2 + \|\vec{c}^*\|^2 - 2\langle\vec{x}_k, \vec{c}^*\rangle) \right) \left( \sum_{k = 1}^n \vec{x}_k \right)\\ = \left( \sum_{k = 1}^n \|\vec{x}_k\|^2\right) \left( \sum_{k = 1}^n \vec{x}_k \right) + n\|\vec{c}\|^2 \sum_{k = 1}^n \vec{x}_k - 2\left( \sum_{k = 1}^n \langle\vec{x}_k, \vec{c}\rangle \right) \left( \sum_{k = 1}^n \vec{x}_k \right), \end{gather*} and (1) is equivalent to\begin{gather*} 2\sum_{k = 1}^n \left( \langle\vec{x}_k, \vec{c}\rangle - \frac{1}{n} \sum_{j = 1}^n \langle\vec{x}_j, \vec{c}\rangle \right) \vec{x}_k = \sum_{k = 1}^n \|\vec{x}_k\|^2 \vec{x}_k - \frac{1}{n} \left( \sum_{k = 1}^n \|\vec{x}_k\|^2 \right) \left( \sum_{k = 1}^n \vec{x}_k \right).\tag{2} \end{gather*} Now,\begin{align*} &\peq \sum_{k = 1}^n \left( \langle\vec{x}_k, \vec{c}\rangle - \frac{1}{n} \sum_{j = 1}^n \langle\vec{x}_j, \vec{c}\rangle \right) \vec{x}_k = \sum_{k = 1}^n \langle \vec{x}_k - \vec{\overline{x}}, \vec{c} \rangle \vec{x}_k\\ &= \sum_{k = 1}^n ((x_k - \overline{x}) x_c + (y_k - \overline{y}) y_c + (z_k - \overline{z}) z_c) \begin{bmatrix}x_k \\ y_k \\ z_k\end{bmatrix}\\ &= \begin{bmatrix} x_c \sum\limits_{k = 1}^n x_k (x_k - \overline{x}) + y_c \sum\limits_{k = 1}^n x_k (y_k - \overline{y}) + z_c \sum\limits_{k = 1}^n x_k (z_k - \overline{z})\\ x_c \sum\limits_{k = 1}^n y_k (x_k - \overline{x}) + y_c \sum\limits_{k = 1}^n y_k (y_k - \overline{y}) + z_c \sum\limits_{k = 1}^n y_k (z_k - \overline{z})\\ x_c \sum\limits_{k = 1}^n z_k (x_k - \overline{x}) + y_c \sum\limits_{k = 1}^n z_k (y_k - \overline{y}) + z_c \sum\limits_{k = 1}^n z_k (z_k - \overline{z}) \end{bmatrix}\\ &= \begin{bmatrix} \sum\limits_{k = 1}^n x_k (x_k - \overline{x}) & \sum\limits_{k = 1}^n x_k (y_k - \overline{y}) & \sum\limits_{k = 1}^n x_k (z_k - \overline{z})\\ \sum\limits_{k = 1}^n y_k (x_k - \overline{x}) & \sum\limits_{k = 1}^n y_k (y_k - \overline{y}) & \sum\limits_{k = 1}^n y_k (z_k - \overline{z})\\ \sum\limits_{k = 1}^n z_k (x_k - \overline{x}) & \sum\limits_{k = 1}^n z_k (y_k - \overline{y}) & \sum\limits_{k = 1}^n z_k (z_k - \overline{z}) \end{bmatrix} \begin{bmatrix}x_c \\ y_c \\ z_c\end{bmatrix}, \end{align*}\begin{gather*} \sum_{k = 1}^n \|\vec{x}_k\|^2 \vec{x}_k = \sum_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) \begin{bmatrix}x_k \\ y_k \\ z_k\end{bmatrix} = \begin{bmatrix} \sum\limits_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) x_k\\ \sum\limits_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) y_k\\ \sum\limits_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) z_k \end{bmatrix},\\ \frac{1}{n} \left( \sum_{k = 1}^n \|\vec{x}_k\|^2 \right) \left( \sum_{k = 1}^n \vec{x}_k \right) = \sum_{k = 1}^n \|\vec{x}_k\|^2 \vec{\overline{x}}\\ = \sum_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) \begin{bmatrix}\overline{x} \\ \overline{y} \\ \overline{z}\end{bmatrix} = \begin{bmatrix} \sum\limits_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) \overline{x}\\ \sum\limits_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) \overline{y}\\ \sum\limits_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) \overline{z} \end{bmatrix}, \end{gather*} then (2) is equivalent to$$ 2\begin{bmatrix} \sum\limits_{k = 1}^n x_k (x_k - \overline{x}) & \sum\limits_{k = 1}^n x_k (y_k - \overline{y}) & \sum\limits_{k = 1}^n x_k (z_k - \overline{z})\\ \sum\limits_{k = 1}^n y_k (x_k - \overline{x}) & \sum\limits_{k = 1}^n y_k (y_k - \overline{y}) & \sum\limits_{k = 1}^n y_k (z_k - \overline{z})\\ \sum\limits_{k = 1}^n z_k (x_k - \overline{x}) & \sum\limits_{k = 1}^n z_k (y_k - \overline{y}) & \sum\limits_{k = 1}^n z_k (z_k - \overline{z}) \end{bmatrix} \begin{bmatrix}x_c \\ y_c \\ z_c\end{bmatrix} = \begin{bmatrix} \sum\limits_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) (x_k - \overline{x})\\ \sum\limits_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) (y_k - \overline{y})\\ \sum\limits_{k = 1}^n (x_k^2 + y_k^2 + z_k^2) (z_k - \overline{z}) \end{bmatrix}. $$

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ Nicely explained. :) $\endgroup$ – SpiritLevel Apr 14 at 7:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.