레이블이 KL transform인 게시물을 표시합니다. 모든 게시물 표시
레이블이 KL transform인 게시물을 표시합니다. 모든 게시물 표시

2014년 1월 5일 일요일

Object tracking with sparse representation (L1 tracker)

Sparse 표현은 최근 cv 분야에서 폭 넓게 인용되고 연구되는 주제이다. 물체 추적 문제에서 가림(occlusion) 문제의 해결은 여전히 어려운데, sparse표현이 이를 다루고 있다.

Mei[1]는 L1 tracker를 제안하여 물체 추적 시 부분 가림을 trivial template를 이용하여 해결 하였다.  목표 물체를 표현하는 방법을 몇 가지 살펴 보자. 


(1) PCA에 의한 관찰 물체 표현


그림은 목표 물체의 영상을 표현하기 위해 pca를 적용한 경우를 보여준다.  주축(pca basis) 몇 개와 주축에서의 좌표 값(target coefficient)이 곱해져서 목표 영상을 생성한다: 

y = Uz +e

여기서, y는 관찰 벡터, U는 기저이고, z는 좌표 값, e는 오차이다.  y를 안다면, z는 z=U'y로 계산 가능하다.  잔차(residual)은 ||y-UU'y||^2이다.  


(2) L1 tracker

Mei는 그림에서 보인 바와 같이 sparse표현을 제안 했다. 물체의 가림을 다루기 위해 trivial template를 도입하였다.  물체에서 가림 없이 보이는 부분은 target template가 나타내고, 가려지는 부분은 trivial template에 곱해지는 trivial 계수가 표현한다.  예를 들면, 얼굴의 크기가 32x32라면 trivial 계수는 1024의 요소를 가진 벡터이다.  target 계수는 여러 템플레이트 중에서 목표 물체와 가장 유사한 템플릿의 인덱스를 나타낸다.  

y=Az+e = [A a][z e]' = Bc

여기서, y는 관찰 벡터, A는 템플릿들의 행렬, z는 템플릿을 선택하는 계수, e는 trivial template에 대한 계수 벡터이다. 행렬 a는 32x32 크기의 trivial 템플릿(즉, 1024크기의 벡터)가 1024개 있는 행렬로(즉, 1024x1024이다) 하나의 벡터에는 해당 픽셀에 대응되는 위치에 1이 저장되어 있고 나머지가 모두 0인 벡터이다.   trivial 템플릿의 계수벡터 e는 1024 벡터이다.  계수 벡터의 값이 occlusion을 표현하게 된다. 
   
수식의 풀이는 L1 minimization을 통한다.  


여기서, ||.||2는 l2놈이고, ||.||1은 l1놈이다.  l1놈은 모든 계수의 절대치의 합의 값으로 생각하면 된다.  식의 의미를 살펴보면, Bc는 가능하면 관찰 벡터 y와 일치해야 하고, c벡터는 sparse(l1놈이 작아지려면 벡터 내에 값이 듬성듬성 들어 있어야 한다)해야 한다는 의미이다.  즉, occlusion은 최소화하면서 y와 잘 일치하도록 계수치들을 결정하면 된다.  단점은 연산량이 크다. 


(3) PCA 기저를 이용한 L1 tracker 

Wang[2]은 SKL과 L1 tracker를 결합하였다.  그림 (c)는 위 그림 (b)의 target 템플릿을 pca 기저(basis)로 바꾼 것이다. 여러 템플릿 영상을 pca의 기저로 표현하면 템플릿의 누적, 갱신 등이 효과적이다.  즉, SKL(sequential KL)을 적용할 수 있다. 

y = Uz+e = [U a][z e]'

여기서, y는 관찰 벡터, U는 pca기저, z는 기저에 대한 계수, e는 trivial 템플릿에 대한 계수벡터이다.    이 식을 다시 쓰면


이다.  위 그림에서 (b)와 (c)를 비교해 보면, (b)의 경우는 물체나 trivial 템플릿 계수 벡터나 둘 다 sparse하다.  즉, target template들은 모두 하나의 물체를 나타내는 영상이므로, 이들 사이에는 상관성이 존재한다. (c)의 경우는 pca기저 벡터가 서로 독립(orthonormal)하기 때문에 기저 계수 벡터는 sparse하지 않고, trivial 계수만 sparse하다.  일반적으로 pca기저는 몇개 안되며(wang논문에서는 16개 사용), trivial template수는 아주 많다(wang에서는 1024개).

식에서 y-Uz는 y가 Uz(pca 기저에 의해 표현되는 템플릿 패치)에 가까와야 한다는 것을 표현하고, y-Uz의 값이 큰 부분(즉, 가림이나 background clutter)은 e에 의해 제거(보상)된다. 단, e가 패치 내의 모든 픽셀에 대해 작용하지 않도록 y-Uz값을 분석해서 잔차가 큰 부분에만 선택적으로 e가 작용하게 한다. 



참고 논문
[1] X.Mei, H. Ling, Robust visual tracking using L1 minimization, CVPR 2009.
[2] D.Wang, et.al., Online object tracking with sparse prototypes, IEEE Trans. IP, 2011.






  

2013년 11월 9일 토요일

Sequential KL 변환

영상 내에 있는 특정 물체(예를 들면 얼굴)를 추적 한다고 하자. 이 때, 물체를 표현하기 위해 여러 개의 영상 템플릿을 사용하는 PCA(principal component analysis)을 적용할 수 있다.

물체의 외관은 추적 동안 여러가지 원인으로 크게 변한다. PCA 기법을 통해 물체를 추적할 때, 물체를 표현하는 주축(principal axis)을 효과적으로 갱신하는 것이 필요하다.

연속 영상에서 물체 추적동안 새로 얻어지는 프레임들을 고려하여 주축을 갱신할 필요가 있을 때 sequential KL 알고리즘을 사용할 수 있다. (Sequential PCA 또는 Sequential visual tracking으로 부를 수도 있다.)


skl 변환은 sequential kl 변환이다.  연속 신호에서 데이터가 차례로 도착한 것을 처리하거나, 큰 행렬의 SVD를 계산해야 할 때, 큰 행렬을 작은 행렬로 분할하여 여러번 반복 처리할 때 사용한다.

예를 들면, A'=[A E] 행렬을 생각해 보자. A행렬은 10장의 영상 템플릿으로 이루어진 행렬이다. E행렬은 새로 매칭된 5장의 프레임에서 추출된 5장의 영상 템플릿으로 이루어진 행렬이다.
이 때, 이용할 수 있는 템플릿이 모두 15개가 되었으므로, 이것 모두가 포함된 A' 행렬에 대해 새로 주축을 구하는 것이 필요하다. 그런데, A행렬에 대해서는 특이값 분해를 통해 이미 구해 놓은 주축이 있다고 가정하자.
효율적 계산을 위해서는 A'행렬을 다시 특이값 분해를 하는 일 없이, 새롭게 추가된 E행렬에 대해서만 특이값 분해를 하고, 그 결과를 A행렬에 대한 특이값 분해의 결과와 결합한다. 이렇게 결합된 결과가 A'행렬에 대한 특이값 분해의 값과 동일하도록 만드는 일이 SKL이다.    

R-SVD[2]는 A의 분해에 기반하여 효과적으로 행렬 A'의 분해 A' = [A|E] = U'D'V''를 계산할 수 있다.  [U,D,V]=svd(A)이라고 가정하자.  이 때, A'의 특이값 분해는


과 같이 구해진다.  이 알고리즘을 해석해 보자.





MxL크기의 행렬 B가 있다.  특이값 분해하면

B = UDV'

이다.  U는 왼쪽 특이 벡터들로 구성된 MxL크기의 직교(orthonormal) 행렬이다. D는 LxL 크기의 행렬이다. 대각요소가 B행렬의 특이 값이다.
V는 각 열이 오른쪽 특이 벡터인 LxL의 직교 행렬이다.


이 때, 분해해야 할 더 큰 행렬 B*가

B*=(B|E)

로 주어 졌다.
E는 P개의 추가 (이미지)데이터가 있는 MxP 크기의 행렬이다.  B*의 특이값 분해를 수행하자.

(Step 1) 행렬 (U|E)에 직교화(orthonormalization) 과정(주3)을 수행한다.  직교 행렬 Ua = (U|Et)를 만든다. E에서 Et를 구하기 위해 QR분해[1]를 이용한다.    (참고: 표기에서 Et(E tilt)는 E'(transpose)가 아님. 직교 행렬임)

(Step 2) (L+P)x(L+P)크기의 Va행렬을 만든다.


Ip는 PxP의 단위 행렬이다. Va는 직교 행렬이다.

(Step 3) B*=Ua.Da.Va' 이므로, Da = Ua'.B*.Va 이다. 대입하면


이다.
Et'.B.V = 0(주1)이고, U'.B.V=D이다. 따라서


이다. Da는 D 행렬에 오른쪽 P개의 열이 추가된 새로운 행렬이다.

Da 행렬의 각 요소에 대해 살펴 보자.

.D는 기저 U로 표현 했을 때의 B 행렬에 대한 좌표 값이다.
.U'E는 새 데이터 E가 U기저에 투영된 좌표이다.  즉, U기저에서 E의 표현이다.
.Et'E는 새 데이터가 Et 기저에 투영된 좌표이다. 즉, Et기저에서 E의 표현이다.
즉, 새로 추가된 데이터를 새로 확장된 기저 (U|Et)에 대해 표현한 좌표 값이 (U'E; Et'E)이다.  이를 다시 쓰면 (U'; Et').E이다.

대각요소가 아닌 요소들이 있으므로 Da를 SVD하자.

Da = Ub.Db.Vb'

이고, B*의 SVD는

B* = Ua.Da.Va'= Ua.(Ub.Db.Vb').Va' = (Ua.Ub).Db.(Va.Vb)'

이다.
여기서, Db는 (L+P)x(L+P)의 대각 형렬이다. (Ua.Ub)는 Mx(L+P)의 column-orthonormal 행렬이고, (Va.Vb)는 (L+P)x(L+P)크기의 orthonormal 행렬이다.

상기한 바와 같이 B*의 SVD는 Ua, Da, Va의 계산과 작은 크기인 Da행렬의 SVD만을 요구한다.  Ua, Da, Va의 계산은 B행렬의 SVD에서 얻어진 데이터를 사용한다.



(예1) 간단한 행렬에서 테스트

>> clear; clc;
>> A = rand(10,5); % 행렬 A를 10x5크기로 생성(old 데이터로 가정)
>> B = rand(10,4); % 행렬 B를 생성 (new 데이터로 가정)

>> [u0,d0,v0] = svd([A B]); % 두 행렬의 결합 확장 행렬의 특이값 분해
     % 기존 데이터와 새 데이터의 모든 데이터가 결합된 행렬의 기저벡터를 얻을 수 있음
     % A, B가 합해진 전체 데이터이고, 이에 대한 u0, d0가 필요.
     % 그런데, 비디오처리를 할 경우, 계속 새 영상
     % 데이터가 입력된다고 보면, [A B C ...]로 행렬 계속 커짐
     % 큰 행렬의 svd는 계산량, 메모리 문제가 있음


%*** (1) skl 변환 [2]
%     원래의 고전적 방법
>> [ua,da,va]=svd(A); % A를 특이값 분해 (기존 데이터)
>> [q,r]=qr([ua*da B]); % new 데이터 B의 직교벡터 추출
     % q행렬 내에는 직교벡터 ua 5개와 B의 직교벡터 4개가 얻어짐
     % A대신 직교행렬과 특이행렬의 곱 ua*da 넣음
>> [u1,d1,v1]=svd(r); % q*r = q*u1*d1*v1'이므로
     % *** u0 == q*u1, d0 == d1 ***
     % A에서 ua, da만 사용했으므로 q*r는 [A B]가 아님

>> ua=q*u1; % 다음 단계를 위해 ua 갱신
>> da=d1; % 다음 단계를 위해 da 갱신
     % 즉 기존 데이터의 u,d만 유지하면 새 데이터가 입력될 때 계산량, 메모리 문제 해결



%*** (2) 보다 개선된 방법
%  [ua*da B]를 qr분해하는 것은 비효율적
% (이미 ua는 직교행렬로 분해되어 있기 때문)
% 따라서 추가되는 B만 qr분해 하면 됨
>> o1=B-ua*ua'*B; % B를 ua기저에 사영하고, 기저 ua로 표현되지 않는 B의 요소 추출(주4)
>> [q1,r1]=qr(o1); % 잔차를 가지고 새로운 직교벡터 q1을 얻음
>> rr=[da   ua'*B; zeros(10,5)   q1*o1];
     % ua'*B: B를 ua기저에서 표현, q1*o1: 잔차의 q1기저 표현
>> [u2,d2,v2]=svd(rr); % rr(20x20)
>> uu = [ua q1]*u2; % uu=10x20. 결합 데이터의 기저
>> dd = d2;  % u2=20x20. 결합 데이터의 특이값
>> ua = uu(:,1:10); % u0와 같음
>> da = dd(1:10,1:9); % d0와 같음





위 알고리즘은 새로운 데이터가 도착해서 이용가능한 데이터가 늘어났을 때, 기존 행렬과 새로 추가되는 행렬 모두에 대해 zero mean을 가정한다. 즉, 평균 값의 변화를 고려하지 않았다.   따라서 sample mean의 변화를 고려하면서 주축을 갱신하는 것이 필요하다.

다시 말하면, C=[A B]행렬의 특이 값 분해를 할 때, C행렬에 대해 C의 평균을 제거한 후, 특이값 분해를 수행해야 한다.  A 및 B행렬의 특이 값 분해도 마찬가지이다. 그런데, A행렬, B행렬의 평균은 모두 C행렬의 평균과 다르며 제 각각이다.  즉, 평균은 항상 달라진다. 이를 고려하여 특이값 분해를 할 필요가 있다.


%*** (3) 평균 변화를 고려한 방법 [3]
%   평균 변화를 고려한 방법. 새로운 데이터가 추가되면 (즉, 기존 A
%   행렬에 B행렬이 추가되어 C행렬이 되었다면
%   A평균, B평균, C평균 값이 다 다르다) 평균도 달라지므로 이를
%   고려해야 한다 [3].
>> C = [A B];
>> n = size(A,2); m = size(B,2);
>> ma = mean(A,2); mb=mean(B,2);
>> mc = n*ma/(n+m) + m*mb/(n+m);
>> Ma = repmat(ma,[1,5]);
>> Mb = repmat(mb,[1,4]);
>> Mc = repmat(mc,[1,9]);
>> Am = A-Ma; % A 행렬을 zero-mean으로
>> Bm = B-Mb; % B            "
>> Cm = C-Mc; % C            "
>>
>> [u0,d0,v0]=svd(Cm); % 목표값(true).
                            % 즉, Am, Bm에 skl 적용으로 u0, d0를 구해야 함.


참고문헌 [3]은 평균변화를 고려한 Sequential KL변환에 대해 잘 설명하고 있다.  참고문헌에서 scatter matrix Sc가 나오는데 covariance matrix와 다음의 관계가 있다:

Sc = (n+m).cov(C);
cov(C) = Sc/(n+m);
cov(C) = CC'=(UDV')(VDU')=U(D^2)U';

따라서 C행렬의 scatter 행렬 Sc의 svd를 구하는 것은 어떤 행렬 C의 U와 D^2을 구하는 것과 같다.


그런데, 참고문헌 [3]의 appendix에 따르면 평균이 제거된 C행렬의 scatter행렬은

Sc = Sa+Sb+ sqrt(.)(Ib-Ia)
   = [(A-Ia) (B-Ib) sqrt(.)(Ib-Ia)] [(A-Ia) (B-Ib) sqrt(.)(Ib-Ia)]'
   = Bt.Bt';

임이 증명되어 있다.

즉, C의 제곱(C.C')이 행렬 Bt의 제곱으로 이루어져 있다.  따라서 C의 svd는 Bt행렬의 svd로 구하면 된다.  Bt 행렬은 (A-Ia)부분과 [(B-Ib) sqrt(.)(Ib-Ia)]의 두 부분의 결합으로 이루어져 있으며, skl을 적용할 수 있다. 즉, (A-Ia)의 svd 결과와 추가되는 [B-Ib sqrt(.)(.)]행렬의 svd를 결합한다.


>> [ua,da,va]=svd(Am);
>> v1 = sqrt((n*m)/(n+m))*(mb-ma);
>> Bh = [Bm v1];
>> o1 = Bh-ua*ua'*Bh; % Bh를 기저 ua로 표현 후의 잔차
>> [q1,r1]=qr(o1); % 잔차에서 새로운 기저 추출
>> rr = [da   ua'*Bh; zeros(10,5)   q1*o1];
>> [u2,d2,v2]=svd(rr);
>> dd = d2;
>> uu = [ua q1]*u2;

% C는 열이 9개이나 d2의 의미 있는(0이 아닌) 특이값이 9개가 안되는 경우
% 이 경우는 8개만 나옴
% uu의 1~8열까지만 유의미 함
>> dd = diag(dd);
>> cutoff = sum(dd.^2)*1e-6; % 임계치 정의
>> keep = find(dd.^2>=cutoff);
>> da = dd(:,keep);
>> ua = uu(:,keep);




(주1) Et'BV = Et'(UDV')V = (Et'U)D(V'V) = (0)D=0
(주2) SVD(X,0) 명령은 SVD(X)보다 더 작은 크기의 U, D, V행렬을 만든다.  만일 크기 mxn의 X 행렬에서 m>n이면 U의 첫번째 n열만 계산한다.  S의 크기는 nxn이다.   만일 m<=n이면 SVD(X,0)는 SVD(X)와 같다.
(주3) [U|E]행렬을 직교화해야 한다. U는 이미 직교행렬이다. E는 새로 추가된 행렬이고 직교행렬이 아니다.  직교화를 위해 qr분해를 이용하는데, qr분해는 기존 벡터들에 수직한 새로운 벡터를 추가해 나가는 방식으로 직교 벡터를 만들어 간다.  svd와 qr분해에서 만들어내는 직교 벡터들 사이의 관계를 통해 이해한다 (보충 설명 필요).
(주4) o1=B-ua*ua'*B:  기저 ua로 표현되지 않는 B의 요소 추출한다는 의미를 살펴보자. svd를 통한 행렬의 주축은 데이터의 산포(scatter)가 큰 몇 개의 축을 선택한다는 의미이다. 이 축들의 아닌 방향으로는 데이터의 산포가 없다는 의미이다.  그런데, 새로운 데이터 B는 기존 주축으로 표현이 되는 데이터도 있고 그러지 않은 데이터도 있을 것이다. 따라서, 기저 ua를 따르지 않는 데이터는 이를 표현하기 위한 새로운 주축이 필요하다.  즉, 잔차 o1에 대한 주축을 새로 선정하면 된다.  



[1] QR분해: A행렬을 QR 분해하면 각 열이 서로 수직이고 열벡터의 크기가 1인 직교(orthonormal) Q행렬과 상삼각행렬(upper triangular matrix) R의 곱으로 분해 할 수 있다. 즉, A=Q.R을 얻는다.
[2] A. Levy and M. Lindenbaum, Sequential Karhunen-Loeve Basis Extraction and its application to Images, IEEE Trans. Image Processing 9(8), 2000.
[3] D.A. Ross, J. Lim, et al., Incremental Learning for Robust Visual Tracking

KL 변환

먼저 차원축소 부분을 읽기 바란다.

Karhunen-Loeve (KL) transform으로 불린다. 저 차원의 부분 공간으로 아주 큰 영상이나 벡터의 집합을 근사하는 방법이다.

mxn 크기의 2차원 배열로 주어지는 영상은 M(=mxn) 크기의 벡터 하나로 바꿀 수 있다. 만일 영상이 N개가 있다면, 전체 데이터의 크기는 MxN이다.
MxN 데이터에 대해 기저를 찾아보자.

일반적으로 MxN 벡터 집합에 대한 KL 기저(basis)를 계산하는 것은 O(MN^2)의 연산과 O(MN)의 메모리를 요구한다.  이러한 계산량은 데이터의 실시간 처리를 위해서는 상당히 크다.


KL 변환: N개의 (images) 벡터 (Ai),i=1...N가 주어질 때, K (K<N)개의 직교(orthonomal) 기저벡터의 집합을 찾는 것이다. 단, 기저 Uj에 의해 확장(span)된 공간에 사영(projection)된 원 데이터의 좌표로 재생한 데이터는 원 데이터를 잘 근사해야 한다.

행렬 A는 크기가 MxN (M>>N)이다.
A의 특이값 분해(SVD: Singular Value Decomposition)를 이용하여 기저벡터를 얻는다.  이때, KL 기저는 (특이값이 0이 되지 않는) K개의 왼쪽 특이 벡터(U의 첫번째 K개의 열)이다.

R-SVM[1]은 SVD를 수정한 것으로 M>>N일 때 효과적이다:

(1) A = QR: qr분해를 통해 (열)직교(column-orthonormal) Q(MxN) 행렬, 상 삼각 행렬 R(NxN)을 얻는다. 계산량은 4MN^2이다.
(2) R = Ua.Da.Va': 표준 SVD알고리즘으로 R을 특이값 분해 한다. 계산량은 O(N^3)이다.
(3) U = Q.Ua, 계산량은 2MN^2이다.

정리해서 다시 쓰면,

A = QR = Q(Ua.Da.Va') = (Q.Ua).Da.Va' = UDV'

이다. 따라서, U = Q.Ua이다.

계산량을 살펴보면 (1)과 (3)이 6MN^2이고, (2)가 N^3이다. 즉 (1),(3)번이 좌우하고 (2)는 무시할 만하다.


(예1) 간단하게 테스트해 보자.

>> clear; clc;
>> A = rand(10,9); % 행렬 A를 10x9크기로 생성
>> [u0,d0,v0] = svd(A); % 행렬의 특이값 분해
>>
>> [q,r] = qr(A); % A 행렬의 QR 분해
>> [u1,d1,v1] = svd(r); % r의 특이값 분해
>> q*u1 % u0와 같음, 즉, u0.d0.v0' = (q*u1).d1.v1'이다.
>> d1 % d0와 같음


[1] G. Golub and C. Van Loan, Matrix Computation. Johns Hopkins Univ. Press, 1992.