はじめに
点と点の最短距離は、ユークリッド空間（一般的な場合）では直線距離になる。これを厳密に証明するためには、微分に似た分野である汎関数を用いた変分法を用いる必要性がある。しかし、シミュレーションだけならば、そこまで難しい数学を使用する必要性はないと考えられる。そこで今回は、両端を固定した10個の点の距離の総和を考え、それらが最小になるのは、以下のように同一直線状にあるということを感覚的に理解することを目的とする。
また、次に示すグラフの通り、最急降下法で最小値を探索していった
場合、最小値は直線距離に漸近することが分かった。
各点の間の距離の総和
$xy$平面において$y=0$上の$(0,0)$から$(1,0)$の範囲に$x$方向において、等間隔に$n$点を配置する。ただし、$(0,0)$と$(1,0)$は不変とし、それ以外の点の$y$座標は、0から1のランダムの値を初期値とする。
原点から数えて、$i$番目の点の座標を$(x_i,y_i)$とおくと、$i,i+1$番目の点のユークリッド距離$l$は、以下のように表すことができる。
l_{i}=\sqrt{(x_{i+1}-x_{i})^2+(y_{i+1}-y_{i})^2}
したがって、それらの総距離を$error$とすると以下の式が成立する。
error(y_2,y_3\cdot\cdot\cdot y_{n-2},y_{n-1})=\sum_{i=1}^{n-1} l_{i}=\sum_{i=1}^{n-1}\sqrt{(x_{i+1}-x_{i})^2+(y_{i+1}-y_{i})^2}
ただし、ここで気を付けて欲しいのは、両端の点を固定してあるので、$y_1,y_n$は定数であり0である。
上記の$error$関数は2次関数がベースとなっているようなので、最小値を1つだけ持つと今回は天下りに仮定する。そうすると以下に示す、最急降下法を使用することができる。
具体的には、$error(y_2,y_3\cdot\cdot\cdot y_{n-2},y_{n-1})$の最小値を最急降下法で推定する。
最急降下法
以下、2次元空間における最急降下法を紹介する。これは、次元が増えても同様の考え方で拡張できる。
解空間の最小値を取る極値が一つ解空間の谷が一つの場合に有効な最小値の推定手法である。$(a,b)$に初期値を与え、以下のような式を用いて$(a,b)$を更新していく。つまり、$(a_k,b_k)$を算出する。
a_{k+1}=a_{k}-\alpha(\frac{\delta }{\delta a_{i}}Error(a_{k},b_{k}))
b_{k+1}=b_{k}-\beta(\frac{\delta }{\delta b_{k}}Error(a_k,b_k))
ただし$\alpha,\beta$は$(a,b)$の学習率である。
今回は、学習率を全て、$\alpha$とした。
プログラム
目的関数の推移のグラフとアニメーションを描写するプログラムを以下に示す。
目的関数の最適化の推移
import numpy as np
import matplotlib.pyplot as plt
import japanize_matplotlib
import math
import matplotlib.animation as animation
# ## アニメーション用の設定
# fig = plt.figure()
# plt.xlim(0,1)
# plt.ylim(0,1)
ims = []
m=300
h=10e-5
alpha=0.01
n=10
x=np.linspace(0,1, n)
y1=np.random.rand(n-2)
y0=np.array([0])
y=np.concatenate([y0,y1], axis=0)
y=np.concatenate([y,y0], axis=0)
#print(y)
# plt.plot(x, y, marker='o', linestyle='-', color='b', label='データ')
# plt.show()
def error(y):
sum=0
for i in range(1,n):
dx=x[i]-x[i-1]
dy=y[i]-y[i-1]
d=math.sqrt(dx**2+dy**2)
sum=sum+d
return sum
k_ary=[]
results=[]
for k in range(m):
d_errs=[]
err=error(y)
for i in range(1,n-1):
y[i]=y[i]+h
err2=error(y)
d_err=(err2-err)/h
d_errs.append(d_err)
y[i]=y[i]-h
#print(d_errs)
d_errs=np.concatenate([y0,d_errs,y0], axis=0)
y=y-np.array(d_errs)*alpha
# ## アニメーション用のプロット
# im = plt.plot(x, y, marker='o', linestyle='-', color='b')
# ims.append(im)
results.append(error(y))
k_ary.append(k)
# 目的関数の推移のグラフ化
plt.plot(k_ary,results)
plt.xlabel('試行回数')
plt.ylabel('目的関数（総距離）')
plt.savefig("目的関数.png")
plt.show()
# # 10枚のプロットを 100ms ごとに表示
# ani = animation.ArtistAnimation(fig, ims, interval=100)
# #plt.show()
# ani.save("output.gif", writer="imagemagick")
これを実行すると以下のような、グラフが出力される。
このように、試行回数を増やすにつれて、目的関数は減少し、直線距離である1に漸近する。
アニメーション
import numpy as np
import matplotlib.pyplot as plt
import japanize_matplotlib
import math
import matplotlib.animation as animation
## アニメーション用の設定
fig = plt.figure()
plt.xlim(0,1)
plt.ylim(0,1)
ims = []
m=300
h=10e-5
alpha=0.01
n=10
x=np.linspace(0,1, n)
y1=np.random.rand(n-2)
y0=np.array([0])
y=np.concatenate([y0,y1], axis=0)
y=np.concatenate([y,y0], axis=0)
#print(y)
# plt.plot(x, y, marker='o', linestyle='-', color='b', label='データ')
# plt.show()
def error(y):
sum=0
for i in range(1,n):
dx=x[i]-x[i-1]
dy=y[i]-y[i-1]
d=math.sqrt(dx**2+dy**2)
sum=sum+d
return sum
k_ary=[]
results=[]
for k in range(m):
d_errs=[]
err=error(y)
for i in range(1,n-1):
y[i]=y[i]+h
err2=error(y)
d_err=(err2-err)/h
d_errs.append(d_err)
y[i]=y[i]-h
#print(d_errs)
d_errs=np.concatenate([y0,d_errs,y0], axis=0)
y=y-np.array(d_errs)*alpha
## アニメーション用のプロット
im = plt.plot(x, y, marker='o', linestyle='-', color='b')
ims.append(im)
results.append(error(y))
k_ary.append(k)
## 目的関数の推移のグラフ化
# plt.plot(k_ary,results)
# plt.xlabel('試行回数')
# plt.ylabel('目的関数（総距離）')
# plt.savefig("目的関数.png")
# plt.show()
# 10枚のプロットを 100ms ごとに表示
ani = animation.ArtistAnimation(fig, ims, interval=100)
#plt.show()
ani.save("output.gif", writer="imagemagick")
これを実行すると以下のようなアニメが出力される。
このように、距離の総和が視覚的にも減少し、直線に漸近するのが分かる。
まとめ
今回は、変分法ではなく、その考え方に似た最急降下法を用いた方法で、視覚的に、点と点の最短距離が直線距離になることをシミュレーションした。厳密考えると、循環論法や距離の定義などいろいろな課題があるが、数値計算的には面白い題材だと思う。このように、目的関数が単純な関数ならば、変数が多くても、最急降下法は有効である。
参考文献
変分法について