9 Stiefel manifolds ぱらぱらめくる『Nonparametric Inference on Manifolds: With Applications to Shape Spaces』

• この章がここに入っているのはなぜなんだろう？？？
• (Direct) Similarity Shape Spaces と、その変形であるRelfection Shape Spacesは行列であって、その成分の２乗和が１という制約を基本としている
  • 
  • これに回転同一視をして、商にしたり、線型独立性制約をつけたり、するバリエーションがある
• 他方、Stiefel manifoldは以下に示すように、行列の形は同じだが、条件は違う
  • テキストとの取り方が異なるが、以下のように書くことにする
  • 
• 本来は次元空間にある個の点が、次元座標を持っているとき、Stiefel manifoldの行列を使うと、次元空間座標に等長変換できる
  • 以下に示す方法で変換行列と変換後座標を算出できる
    • (QR分解)。
  • のときは、以下の方法でもできる
    • ,

library(GPArotation)
library(Matrix)

m <- 3
k <- 5

my.rstiefel <- function(m,k){
	tmp <- Random.Start(k)
	return(tmp[,1:m])
}

x <- my.rstiefel(m,k)

round(t(x) %*% x, 5)

# ランクをm未満にする
real.dim <- m-1

# real.dimランクのm次元ベクトルをk個作る
mu <- matrix(rnorm(real.dim*k),ncol=real.dim)
mu <- cbind(mu,matrix(0,k,m-real.dim))
tmp.rot <- Random.Start(m)
mu <- mu %*% tmp.rot


rankMatrix(mu)

U <- mu %*% solve(sqrtm((t(mu) %*% mu)))
P <- sqrtm(t(mu) %*% mu)

qr.out <- qr(mu)
U. <- qr.Q(qr.out)
P. <- qr.R(qr.out)

# 分解はちゃんと出来ている
round(U %*% P - mu,5)
round(U. %*% P. - mu,5)

# 両者は違う
U - U.
P - P.
# 等長変換になっていない
round(t(P) %*% P - t(mu) %*% mu,5)
round(t(P.) %*% P. - t(mu) %*% mu,5)
# m'=real.dim < mなので、UはStiefel manifoldに乗っていない

round(t(U) %*% U,5)
round(t(U.) %*% U.,5)


# ランクをmにする
real.dim <- m

# real.dimランクのm次元ベクトルをk個作る
mu <- matrix(rnorm(real.dim*k),ncol=real.dim)
mu <- cbind(mu,matrix(0,k,m-real.dim))
tmp.rot <- Random.Start(m)
mu <- mu %*% tmp.rot


rankMatrix(mu)

U <- mu %*% solve(sqrtm((t(mu) %*% mu)))
P <- sqrtm(t(mu) %*% mu)

qr.out <- qr(mu)
U. <- qr.Q(qr.out)
P. <- qr.R(qr.out)

# 分解はちゃんと出来ている
round(U %*% P - mu,5)
round(U. %*% P. - mu,5)
# 両者は違う

U - U.
P - P.
# 等長変換になっている
round(t(P) %*% P - t(mu) %*% mu,5)
round(t(P.) %*% P. - t(mu) %*% mu,5)
# m'=real.dim == mなので、UはStiefel manifoldに乗る

round(t(U) %*% U,5)
round(t(U.) %*% U.,5)