移動とマッチングに関するメモ

  • 空間に点集合が2つ(V_1,V_2)ある。一つ目の集合(V_1)から二つ目の集合(V_2)への移動であると考えて、どの点がどの点に移動したとみなすことにするかを考える
  • 同数の点に関するマッチングの決め方であるとみなせる
  • マッチングにはハンガリアン・アルゴリズムが使える場合がある
    • その条件は、n点からn点へのマッチングであるとして、nxn行列の非負正方行列(K)が定義できて、第1集合(V_1)の第i点(x_{1,i}を第2集合(V_2)の第j点(x_{2,j})にマッチングさせたときにスコアk_{i,j}が与えられて、(1,2,...,n) \to (s_1,s_2,...,s_n)なるマッチングにおいて\sum_{i=1}^n k_{i,s_i}の最大化もしくは最小化問題にできる場合である
  • k_{i,j}の決め方は一通りではない
    • k_{i,j}=\sqrt{\sum_{u=1}(x_{1,i}(u)-x_{2,j}(u))^2}のように、2点のユークリッド距離にしてもよいし
    • k_{i,j}=\sum_{u=1}|x_{1,i}(u)-x_{2,j}(u)|のように、2点のマンハッタン距離にしてもよい
    • k_{i,j}=-\sum_{u=1}x_{1,i}(u)\times x_{2,j}(u)(内積)のようなものを考えてもよい。この場合は、式変形に注意すれば[tex:k_{i,j} = \sum_{u=1} *1 - (\frac{\sum_{i=1}^n x_{1,i}(u)}{n} - \frac{\sum_{i=1}^n x_{2,i}(u)}{n}))^2]と大小関係が保たれていることがわかる
      • この例は、すべての2点ペアを結ぶベクトルの多次元分散の大小に相当している
    • さらに、2点ペアのベクトルが似通っていることを求めるのであれば、すべてのペアのベクトルが似ていることを求めるかわりに、近い関係にあるそれにのみ求めるという方法も考えられ(それは非線形なアプローチなのだが)るが、それをnxn行列にi,jのペアとして実現することはできない(もしくは困難)なので、それはハンガリアン・アルゴリズムでの処理には合わないもののようだ

# ノード対応づけをするために試行錯誤した関数
# マッチングでできたペアのユークリッド距離
library(clue)
my.move <- function(x,y,w.x=NULL,w.y=NULL){
	if(!is.matrix(x)){
		x <- matrix(x,ncol=1)
	}
	if(!is.matrix(y)){
		y < matrix(y,ncol=1)
	}
	if(is.null(w.x)){
		w.x <- rep(1,length(x[,1]))
	}
	if(is.null(w.y)){
		w.y <- rep(1,length(y[,1]))
	}
	D <- x %*% t(y)
	D <- D-min(D)
	xy <- rbind(x,y)
	m.xy <- apply(xy,2,mean)
	x.d <- t(t(x)-m.xy)
	y.d <- t(t(y)-m.xy)
	W <- outer(w.x,w.y,"*")
	D <- D * W
	
	as.vector(solve_LSAP(D,maximum=TRUE))
}
# ペアのベクトルの分散
my.move.D <- function(x,y,w.x=NULL,w.y=NULL){
	if(!is.matrix(x)){
		x <- matrix(x,ncol=1)
	}
	if(!is.matrix(y)){
		y < matrix(y,ncol=1)
	}
	if(is.null(w.x)){
		w.x <- rep(1,length(x[,1]))
	}
	if(is.null(w.y)){
		w.y <- rep(1,length(y[,1]))
	}
	D <- (as.matrix(dist(rbind(x,y))))[1:length(x[,1]),(1+length(x[,1])):(2*length(x[,1]))]
	#D <- D-min(D)
	xy <- rbind(x,y)
	m.xy <- apply(xy,2,mean)
	W <- outer(w.x,w.y,"*")
	D <- D * W
	
	as.vector(solve_LSAP(D,maximum=FALSE))
}
# 相互のベクトルの近さで工夫する??
my.move.L <- function(x,y,w.x=NULL,w.y=NULL){
	if(!is.matrix(x)){
		x <- matrix(x,ncol=1)
	}
	if(!is.matrix(y)){
		y < matrix(y,ncol=1)
	}
	if(is.null(w.x)){
		w.x <- rep(1,length(x[,1]))
	}
	if(is.null(w.y)){
		w.y <- rep(1,length(y[,1]))
	}
	D <- x %*% t(y)
	D <- D-min(D)
	xy <- rbind(x,y)
	m.xy <- apply(xy,2,mean)
	x.d <- t(t(x)-m.xy)
	y.d <- t(t(y)-m.xy)
	L.x <- sqrt(apply(x.d^2,1,sum))
	L.y <- sqrt(apply(y.d^2,1,sum))
	LL <- outer(L.x,L.x,"*")
	L.q <- quantile(LL,0.05)
	LL.q <- sign(LL < L.q)
	W <- outer(w.x,w.y,"*")
	D <- D * W * LL.q
	
	as.vector(solve_LSAP(D,maximum=TRUE))
}
v1 <- runif(2)
v2 <- runif(2)
t <- seq(from=0,to=1,length=100)
x <- matrix(0,length(t),2)
for(i in 1:length(t)){
	x[i,] <- t[i] * v1 + (1-t[i]) * v2
}
x[,1] <- x[,1]-mean(x[,1])
x[,2] <- x[,2]-mean(x[,2])
theta <- pi/3
M <- matrix(c(cos(theta),sin(theta),-sin(theta),cos(theta)),2,2)
x <- rbind(x,x%*%M)
theta <- pi/6
M <- matrix(c(cos(theta),sin(theta),-sin(theta),cos(theta)),2,2)

y <- x %*% M
y <- y + 0.3
y <- (abs(y))^1.3
out <- my.move(x,y)
out.D <- my.move.D(x,y)
out.L <- my.move.L(x,y)
par(mfcol=c(2,3))
plot(rbind(x,y),pch=20,col=rep(1:2,each=length(x[,1])),main="内積")
for(i in 1:length(out)){
	segments(x[i,1],x[i,2],y[out[i],1],y[out[i],2],col=3)
}
plot(out,main="内積")
plot(rbind(x,y),pch=20,col=rep(1:2,each=length(x[,1])),main="ユークリッド距離")
for(i in 1:length(out)){
	segments(x[i,1],x[i,2],y[out.D[i],1],y[out.D[i],2],col=3)
}
plot(out.D,main="ユークリッド距離")
plot(rbind(x,y),pch=20,col=rep(1:2,each=length(x[,1])))
for(i in 1:length(out)){
	segments(x[i,1],x[i,2],y[out.L[i],1],y[out.L[i],2],col=3,main="近傍重視・内積")
}
plot(out.L,main="近傍重視・内積")
par(mfcol=c(1,1))

*1:x_{1,i}(u)-x_{2,j}(u