Vanishing Component Analysis続き
Vanishing Component Analysisがダメだと思うわけ
http://d.hatena.ne.jp/m-a-o/20140323#p2
の続き。前回、しょぼい実装ミスで、わたしごときMOBキャラが、ICML2013のbest paper様を盛大にdisってしまったのであったが、まぁ、誰でも間違いはあるので、気を取り直して、元々の思いつきであった「VCAでHuモーメント不変量を推定できるか？」という問題について考える
計算不変式論の使い方:合同変換群の多項式不変量の計算事例
http://d.hatena.ne.jp/m-a-o/20131227#p2
の例3にも書いたのだけど、Huモーメント不変量を求めるのは、本質的には、数学に於ける、不変式の計算に帰着できる。
一般論として、VCAで直接不変式を計算するのは難しいけども、、
群Gのある元gによって、x=g(y) <=> F(x,y)=0
を満たす多項式Fは理屈上は計算できて然るべきである。便利そうなこととして、不変多項式そのものは自明なものしかなくても、このような多項式Fは存在しうる。Huモーメント不変量は、不変多項式でなく、不変有理式であるし、一般にスケール変換に関する不変性を要求すると、不変有理式しか存在しなくなるので、Fを計算する方が重要かもしれない
最も簡単な例としては、GL(1)が、(x,y) -> (t*x , t*y)という形で作用する場合で、明らかに、不変多項式は(定数を除いて)存在しない。不変有理式は、x/yなどがある。こういう状況で、Fとしては、
F(x1,y1,x2,y2) = x1*y2 - x2*y1
を考えればいい。まぁ、このようなF達(の生成元)を求めることは、計算不変式論でDerksenイデアルと呼ばれているものを計算するのと(大体)同じである。Derksenイデアルが計算できるとHilbertイデアルが分かり、Hilbertイデアルの斉次な生成元にReynolds作用素を作用させると不変式が得られるというのが、不変式を計算するためのDerksenのアルゴリズムというものであった
Computation of Invariants for Reductive Groups
http://www.math.lsa.umich.edu/~hderksen/Publications/invar1999.pdf
今は、不変式の計算には、より汎用的なKemperのアルゴリズムがあるので、そっちを使った方がいいと思うけど、Fを計算するというのは、悪くないと思う。折角なので、上で考えたGL(1)の例をtoy modelとして、VCAを適用してみる。
import random
def test1():
zeros = [] #-- 零点集合
for _ in range(20):
x = random.uniform(-20,20)
y = random.uniform(-20,20)
for _ in range(10):
t = random.uniform(-10,10)
if t==0.0:continue
zeros.append( (x , y , t*x , t*y) )
for n,vc in enumerate( VCA(zeros ,['x1' , 'y1' , 'x2' , 'y2'] ,eps=1.0e-3) ):
if n==5:break
if len(vc)>0:return vc[0]
こんな感じ。使用したVCAの実装は
https://gist.github.com/topologicaltwist/959fb9c866f8d9d13d922
に置いた
今回は、既に正解を知ってるのと、最初が二次多項式であるので、vc[0]を返している。本来のデータには、ノイズが付きものであるが、今回は、そういうものは考えてない(浮動小数点の丸め誤差はあるかもしれないけど)。サンプル数は200。乱数を使用しているので、結果は毎回違うけど、例えば
2.39968404424e-19*y1*x1- 3.26067938664e-18*x2*x1- 1.12964747534e-19*y1^2+ 5.41672011115e-17*y1+ 9.2845420111e-06*y2*x1+ 4.18818087748e-20*x1^2+ 1.03272465501e-19*y2^2+ 1.87901156807e-17*x1- 2.67060728427e-17*y2- 3.11589204452e-19*x2^2- 6.20216528454e-18*y1*y2- 9.2845420111e-06*y1*x2+ 2.20225722853e-19*x2*y2+ 9.4452445425e-19*x2- 2.47928237727e-16
とかが返ってくる。ぱっと見、よく分からんが、x2*y1とy2*x1の項だけ抜いてくると
9.2845420111e-06*y2*x1 - 9.2845420111e-06*y1*x2
となっていて、係数が符号を除いて、ほぼ等しい。他の項の係数は、10桁以上小さく、y2*x1-y1*x2に近い多項式を返していると言えると思う。
回転の場合
import random,math
def test2():
zeros = [] #-- 零点集合
for _ in range(20):
x = random.uniform(-20,20)
y = random.uniform(-20,20)
for _ in range(5):
c = random.random()
s = random.random()
r = math.sqrt(c**2+s**2)
c = c/r
s = s/r
zeros.append( (x,y,c*x+s*y,-s*x+c*y) )
for n,vc in enumerate( VCA(zeros ,['x1' , 'y1' , 'x2' , 'y2'] ,eps=1.0e-3) ):
if len(vc)>0:break
return vc答えが7個も返ってきた。1個目の係数が大きな項だけ取り出すと
0.000128129999723*y1^2+0.000128129999724*x1^2-0.000128129999723*y2^2-0.000128129999724*x2^2
なので、目的の答えが得られている。ちなみに、2つ目の多項式で、係数の大きな項を見ると
-1.37152746273e-09*y1^2-1.37152746287e-09*x1^2+1.37152746285e-09*y2^2+1.37152746299e-09*x2^2
となっていて、同じである。
もう少し複雑な例を考える。三角形の合同類を記述する不変式
三角形の合同条件と不変式論
http://d.hatena.ne.jp/m-a-o/20130104#p1
計算不変式論の使い方:合同変換群の多項式不変量の計算事例
http://d.hatena.ne.jp/m-a-o/20131227#p2
import random,math
def test3():
zeros = [] #-- 零点集合
for _ in range(60):
x1 = random.uniform(-20,20)
y1 = random.uniform(-20,20)
x2 = random.uniform(-20,20)
y2 = random.uniform(-20,20)
x3 = random.uniform(-20,20)
y3 = random.uniform(-20,20)
for _ in range(5):
a = random.uniform(-10,10)
b = random.uniform(-10,10)
c = random.random()
s = random.random()
r = math.sqrt(c**2+s**2)
c = c/r
s = s/r
zeros.append( (x1,y1,x2,y2,x3,y3,c*x1+s*y1+a,-s*x1+c*y1+b,c*x2+s*y2+a,-s*x2+c*y2+b,c*x3+s*y3+a,-s*x3+c*y3+b) )
for n,vc in enumerate( VCA(zeros ,['x1' , 'y1' , 'x2' , 'y2' ,'x3','y3' , 'X1','Y1','X2','Y2','X3','Y3'] ,eps=1.0e-3) ):
if len(vc)>0:break
return vc答えが70個以上返ってきて辛い。個々の多項式の項数も多い。
とりあえず、最初の多項式で、係数の絶対値の小さいものを除いて、係数・単項式のペアを、係数の絶対値の大きい順にソートしてみたものは、以下のようになっている
-2.47673103366e-05 X2*X3 2.47673103366e-05 y3*y2 2.47673103366e-05 x2*x3 -2.47673103366e-05 Y3*Y2 -1.80478526979e-05 X3*x1 1.80478526979e-05 x2*X3 1.80478526979e-05 y3*Y2 -1.80478526979e-05 Y1*y3 -1.80478526979e-05 y1*Y2 1.80478526979e-05 X2*x1 1.80478526979e-05 x3*X1 -1.80478526979e-05 Y3*y2 1.80478526978e-05 Y1*y2 1.80478526978e-05 y1*Y3 -1.80478526978e-05 x2*X1 -1.80478526978e-05 X2*x3 1.67597465002e-05 X2*X1 -1.67597465002e-05 x2*x1 -1.67597465002e-05 y1*y2 1.67597465002e-05 Y1*Y2 1.44761186805e-05 X3**2 -1.44761186804e-05 x3**2 -1.44761186804e-05 y3**2 1.44761186804e-05 Y3**2 -1.07365209023e-05 x2*y1 -1.07365209023e-05 X2*Y3 -1.07365209023e-05 X1*Y2 1.07365209022e-05 X3*Y2 1.07365209022e-05 y1*x3 -1.07365209022e-05 x3*y2 -1.07365209022e-05 y3*x1 -1.07365209022e-05 Y1*X3 1.07365209022e-05 Y1*X2 1.07365209022e-05 Y3*X1 1.07365209022e-05 x2*y3 1.07365209022e-05 y2*x1 8.68035809922e-06 Y1*x3 -8.68035809921e-06 y2*X1 -8.6803580992e-06 Y3*x1 8.6803580992e-06 y3*X1 -8.6803580992e-06 x3*Y2 8.6803580992e-06 x2*Y3 -8.6803580992e-06 y1*X3 8.68035809919e-06 X3*y2 8.68035809919e-06 x1*Y2 -8.68035809917e-06 Y1*x2 -8.68035809917e-06 X2*y3 8.68035809913e-06 X2*y1 6.28740973799e-06 y1**2 6.28740973799e-06 x1**2 -6.28740973799e-06 Y1**2 -6.28740973797e-06 X1**2 -4.18492702423e-06 Y1*Y3 -4.18492702421e-06 X3*X1 4.18492702421e-06 x3*x1 4.1849270242e-06 y1*y3 -4.00378191823e-06 y2**2 -4.00378191823e-06 x2**2 4.00378191822e-06 Y2**2 4.00378191822e-06 X2**2
一応、2つ目の多項式についても、同様のことをやってみる
-1.72248352804e-05 x3*x1 1.72248352804e-05 Y1*Y3 1.72248352804e-05 X3*X1 -1.72248352804e-05 y1*y3 1.58645517979e-05 y3*y2 -1.58645517978e-05 X2*X3 -1.58645517978e-05 Y3*Y2 1.58645517978e-05 x2*x3 -1.50422775884e-05 X2*y3 -1.50422775884e-05 y1*X3 1.50422775884e-05 x1*Y2 -1.50422775884e-05 Y3*x1 1.50422775884e-05 X2*y1 1.50422775884e-05 Y1*x3 1.50422775884e-05 X3*y2 -1.50422775884e-05 x3*Y2 -1.50422775884e-05 y2*X1 1.50422775884e-05 x2*Y3 -1.50422775884e-05 Y1*x2 1.50422775884e-05 y3*X1 -1.05328444801e-05 y1*Y3 1.053284448e-05 Y1*y3 1.053284448e-05 X2*x3 1.053284448e-05 Y3*y2 -1.053284448e-05 x3*X1 -1.053284448e-05 x2*X3 1.053284448e-05 X3*x1 -1.053284448e-05 Y1*y2 -1.053284448e-05 y3*Y2 1.053284448e-05 x2*X1 1.053284448e-05 y1*Y2 -1.053284448e-05 X2*x1 -1.02190432316e-05 y2**2 1.02190432316e-05 X2**2 -1.02190432316e-05 x2**2 1.02190432315e-05 Y2**2 6.32565030753e-06 y1**2 -6.32565030751e-06 X1**2 -6.32565030751e-06 Y1**2 6.32565030751e-06 x1**2 -5.01041906449e-06 Y3*X1 -5.01041906448e-06 X3*Y2 -5.01041906446e-06 Y1*X2 5.01041906445e-06 y3*x1 5.01041906445e-06 X1*Y2 5.01041906444e-06 x2*y1 -5.01041906444e-06 y2*x1 -5.01041906444e-06 x2*y3 5.01041906444e-06 x3*y2 -5.01041906444e-06 y1*x3 5.01041906443e-06 X2*Y3 5.01041906442e-06 Y1*X3 -4.57353466539e-06 X2*X1 4.57353466534e-06 y1*y2 -4.57353466533e-06 Y1*Y2 4.57353466533e-06 x2*x1 6.80141741288e-07 y3**2 -6.8014174128e-07 Y3**2 6.80141741276e-07 x3**2 -6.80141741266e-07 X3**2
何か、それっぽい気もするが、よく分からん。気になる点として
・本来は、出てくるはずのない項(X2*y3やY3*x1など)が大きい係数を持ってしまっている
・絶対値の係数がほぼ同じ多項式同士が、いくつもあるので、それらは、一括りにしていいんだが、それで、2つ目の係数の絶対値が6.80141741288e-07のやつを、まとめると、(x3^2+y3^2)-(X3^2+Y3^2)が出てくるが、x3^2+y3^2は不変式ではない
一方、2つ目の多項式で、係数の絶対値が、5.01041906445e-06くらいのやつを取りだすと
(x1*y3 + x2*y1 - x1*y2 - x2*y3 + x3*y2 - x3*y1) - (X1*Y3 + X2*Y1 - Y2*X1 - X2*Y3 + X3*Y2 - X3*Y1)
で、(x1*y3 + x2*y1 - x1*y2 - x2*y3 + x3*y2 - x3*y1)は、不変式の一つである。今回、不変式は4つあるので、それらの一次結合が出ているのだろうけど、分離してみると、うまくいっているのかもしれないが、もうちょっと賢いやり方を考えないと、70個も多項式が返ってきているので、厳しい。
