平面上の点を近似する直線

プログラムを作っていて、平面上に散らばった点を直線に近づけるように力をかけたいという問題にぶちあたった。

しかし、Google で探してもみつからない (笑) ので、自分で解く羽目になった。せっかくなのでwebに貼っておけば、 誰かが参考にするのではないかと期待してはってみる。間違いがあったら連絡ください。

方向ベクトル $ \vec l = (\cos \theta, \sin \theta)$ の直線 $ l$ の方程式は

$\displaystyle x \sin \theta - y \cos \theta + c = 0$ (1)

点 $ \mathrm{P}_1 = (x_1, y_1)$ から$ l$ におろした垂線の足を $ \mathrm{P}_0 = (x_0,y_0)$ とすると、 式 (1) と $ \overrightarrow{\mathrm{P_1P_0}} \cdot \vec l=0$ より

$\displaystyle (x_0,y_0)$ $\displaystyle = ( x_1 \cos^2 \theta + y_1 \sin \theta \cos \theta - c \sin \theta , x_1 \sin \theta \cos \theta + y_1 \sin^2 \theta + c \cos \theta )$ (2)

直線 $ l$ と点 $ \mathrm{P}_1$ の距離は、公式使って

$\displaystyle d_1 % &= \sqrt{(x_0-x_1)^2+(y_0-y_1)^2} \\ $ $\displaystyle = x_1 \sin \theta - y_1 \cos \theta + c$ (3)

$ \theta$ が与えられているとき、 $ N$ 個の点 $ \mathrm{P}_1 = (x_1,y_1), \ldots, \mathrm{P}_N = (x_N,y_N)$ からの距離の二乗和$ D$ は

$\displaystyle D =$ $\displaystyle \sum d_k^2 = \sum (x_k \sin \theta - y_k \cos \theta + c)^2$ (4)
$\displaystyle =$ $\displaystyle \sum x_k^2 \sin^2 \theta - 2 \sum (x_ky_k) \sin\theta\cos\theta +... ...k^2 \cos^2 \theta + 2 \sum x_k c \sin \theta - 2 \sum y_k c \cos \theta + N c^2$ (5)

式 (5) は $ c$ の二次式なので、$ \theta$ を固定と考えると、$ D$ を最小にする$ c$ は

$\displaystyle c$ $\displaystyle = -\frac{1}{N} \left(\sum x_k \sin \theta - \sum y_k \cos \theta \right)$ (6)

$ D$ を最小にする$ \theta$ を考える。式 (6) を 式 (5) に代入して 倍角の公式とかで整理すると

$\displaystyle D =$ $\displaystyle -\frac{1}{2}\left(\sum x_k^2 - \sum y_k^2 - \frac{1}{N} \left(\Bigl(\sum x_k\Bigr)^2 - \Bigl(\sum y_k\Bigr)^2\right)\right) \cos 2\theta$    
  $\displaystyle + \left( - \sum (x_ky_k) + \frac{1}{N} \Bigl(\sum x_k\Bigr) \Bigl(\sum y_k\Bigr) \right) \sin 2\theta$    
  $\displaystyle + \frac{1}{2}\left(\sum x_k^2 + \sum y_k^2 - \frac{1}{N} \left(\Bigl(\sum x_k\Bigr)^2 + \Bigl(\sum y_k\Bigr)^2\right)\right)$ (7)

ここで

$\displaystyle A$ $\displaystyle = \frac{1}{2} \Bigl( \sum x_k^2 - \frac{1}{N} \Bigl(\sum x_k\Bigr)^2 \Bigr)$ (8)
$\displaystyle B$ $\displaystyle = \frac{1}{2} \Bigl( \sum y_k^2 - \frac{1}{N} \Bigl(\sum y_k\Bigr)^2 \Bigr)$ (9)
$\displaystyle C$ $\displaystyle = - \Bigl( \sum (x_ky_k) - \frac{1}{N} \Bigl(\sum x_k\Bigr) \Bigl(\sum y_k\Bigr) \Bigr)$ (10)

とおくと

$\displaystyle D$ $\displaystyle {} = - (A - B) \cos 2\theta + C \sin 2\theta + (A + B)$ (11)
$\displaystyle \frac{dD}{d\theta} %=& - \left( \left(\sum y_k^2 - \frac{1}{N}\Bi... ... \frac{1}{N} \Bigl(\sum x_k\Bigr) \Bigl(\sum y_k\Bigr)\right) \cos 2\theta \\ $ $\displaystyle {} = 2 (A - B) \sin 2\theta + 2 C \cos 2\theta$ (12)
  $\displaystyle {} = 2 \sqrt{(A-B)^2+C^2} \sin \Bigl(2\theta + \arctan \frac{C}{A-B}\Bigr)$ (13)

$ D$ が最小になるのは $ \frac{dD}{d\theta}=0$ となるときなので、 $ 2\theta + \arctan \frac{C}{A-B}=0$ 、つまり

$\displaystyle \theta$ $\displaystyle = -\frac{1}{2} \arctan \frac{C}{A-B} %\\ $ (14)

このとき

$\displaystyle D$ $\displaystyle = - \sqrt{(A-B)^2+C^2} + A + B$ (15)