唯物是真 @Scaled_Wurm

プログラミング(主にPython2.7)とか機械学習とか

Orthogonal Procrustes problem

前に読んだ論文で出てきたのでメモ

同じ形の行列\(A, B\)が与えられた時に\(||AX - B||_F\)を最小化する行列\(X\)を求める(\(A\)を\(B\)にできるだけ近く変換する行列を求める
ただし\(X^TX=I\)とする(直交行列

\(A^TB\)の特異値分解が\(U\Sigma V^\ast\)のとき、求めるべき行列は\(X=U V^\ast\)で計算できるらしいです
たぶん導出は以下みたいな感じ?(最後の方はよくわからない
\(||AX - B||_F^2=\mathrm{trace}( (AX - B)^T (AX - B))=\mathrm{trace}(A^TAX^TX + B^TB - 2A^TX^TB)\)
\(X^TX=I\)より
\(\mathrm{trace}(A^TAX^TX + B^TB - 2A^TX^TB)=\mathrm{trace}(A^TA + B^TB - 2A^TX^TB)\)
よって\(\mathrm{trace}(A^TX^TB)\)を最大化するような\(X\)を見つければよい
\(A^TB=U\Sigma V^\ast\)より
\(\mathrm{trace}(A^TX^TB) = \mathrm{trace}(U\Sigma V^\ast X^T)\)
したがって\(X=U V^\ast\)の時に\(\mathrm{trace}(U\Sigma V^\ast X^T) = \mathrm{trace}(\Sigma)\)となり最大となる(?)

以下のソースコードでは適当な行列\(A, X\)を考えて、それらの積である\(B\)と、\(A\)を使って\(X\)を逆に推定しています

# -*- coding: utf-8 -*-
import numpy as np
from scipy.linalg import svd

A = np.arange(9)
A.shape = (3, 3)
print u'行列A:'
print A

X_known = np.array([[1, 0, 0], [0, 0, 1], [0, 1, 0]])
print u'推定したい行列X:'
print X_known

B = np.dot(A, X_known)
print u'行列B:'
print B

temp = svd(np.dot(A.T, B))

X_estimated = np.dot(temp[0], temp[2])
print u'推定された行列X:'
print X_estimated

次のような出力が得られて、ちゃんと\(X\)を推定できていることがわかります

行列A:
[[0 1 2]
 [3 4 5]
 [6 7 8]]
推定したい行列X:
[[1 0 0]
 [0 0 1]
 [0 1 0]]
行列B:
[[0 2 1]
 [3 5 4]
 [6 8 7]]
推定された行列X:
[[  1.00000000e+00  -6.10622664e-16   8.32667268e-16]
 [ -9.99200722e-16   7.21644966e-16   1.00000000e+00]
 [  5.82867088e-16   1.00000000e+00  -6.10622664e-16]]