마할라노비스 제곱거리와 카이제곱분포을 통한 다변량 이상치 탐색

 

 

다변량 자료를 다루기 전, 먼저 알고 넘어갈 사항이 있다. 정규분포와 마찬가지로, 다변량 정규분포는 몇몇 다변량 분석 절차에서 중요한 역할을 한다. 또한, 모집단의 모수에 대한 추론은 주로 모집단에 대해 다변량 정규분포의 가정하에 이루어 진다. 

 

1. 다변량 정규분포

 

확률벡터 $X = X(X_1,X_2, \cdots,X_p)^T$가 평균이 $\mu$이고 공분산행렬이 $\Sigma$인 다변량 정규분포를 따를 때 (기호로, $X$ ~ $N_p(\mu,\Sigma$)), 확률밀도함수는 다음과 같다.

$f(x;\mu,\Sigma) = (2\pi)^{-p/2} det(\Sigma)^{-1/2} exp[-\frac{1}{2} (x-\mu)^T \Sigma^{-1} (x-\mu)] $

 

위의 확률밀도함수는 동일한 값을 가지는 지수부에 대해 함수값이 동일하므로,  $(x-\mu)^T \Sigma^{-1} (x-\mu)=C$ 를 만족하는 $x$들의 자취(타원체)를 spectrum decomposition(스펙트럼 분해)를 이용해 전개시킬 수 있다. 

 

$(X-\mu)^T P\Lambda^{-1}P^T(X-\mu) = C$
$(X-\mu)^T \left( \sum_{i=1}^{2} \lambda^{-1} e_i e_i^T \right) (X-\mu) = C$

$\sum_{i=1}^{2} \lambda^{-1} \left [ (X-\mu)^T e_i \right ]^2 = C$ 

$\frac{\left [(X-\mu)^Te_1\right]^2}{\lambda_1} + \frac{\left [(X-\mu)^Te_2\right]^2}{\lambda_2} = C$

 

이를 통해, 상수 c의 변화에 따른 contour을 나타낼 수 있다. $\mu_1,\mu_2$는 타원의 중심이 되고, $\lambda_1,\lambda_2$ 값으로 타원의 장축, 단축의 길이를 구할 수 있다. 또한, 변수를 eigenvector로 정사영시킨 형태이기 때문에 eigenvector의 방향 $e_1,e_2$는 장축, 단축의 방향을 나타낸다. 상수 C의 값에 따라 각각 다른 크기의 타원이 그려지지만 장축과 단축의 비율은 항상 일정하다. 

 

이변량 정규분포에 대한 상수-밀도 등고선 

 

이렇게 다변량 정규분포에 대해 아주 간단히 알아보았다. 오늘 글의 주제와 다변량 정규분포와 무슨 상관이 있을까?

다변량 정규분포의 지수부 형태를 자세히 보면, 마할라노비스 거리 식과 비슷하다는 것을 알 수 있다. 

 

$[(X-\mu)^T \Sigma^{-1}(X-\mu)] \approx [(X_r - X_s)^TS^{-1}(X_r - X_s)]^{1/2}$ 

 

이는 각 변수의 분산과 공분산(상관성)구조를 함께 고려한 통계적 거리로, 확률등고선을 보고 관측될 가능성이 더 높은 점이 거리가 더 가깝다고 판단하는 것이다.  $X_1,X_2,\cdots,X_n$ ~ $N_p(\mu,\Sigma)$이고, $\left | \Sigma\right | > 0$일 때, $X_i$의 $\bar X$ 로부터의 마할라노비스 제곱거리(Mahalanobis squared distance) 또는 일반화 제곱거리(generalized squared distance) $d_i^2$은 근사적으로 $x^2(p)$를 따른다. 직관적으로 생각해보면 당연한 이야기이다. 확률변수들이 다변량 정규분포를 따를 때, 각각의 확률변수는 일변량 정규분포를 따르게 된다. (하지만, 반대의 상황은 성립하지 않는다는 것을 유의하자.) 따라서, 일변량 정규분포를 따르는 확률변수(i)를 표준화하면 표준정규분포를 따르고 제곱하면 자유도가 1인 카이제곱분포를 따르게 된다. 이를 p개 더 했으니 자유도가 p인 카이제곱분포를 따르겟죠 ?

 

$d_i^2 = [(X-\bar X)^T S^{-1}(X-\bar X)]$ ~ $x^2(p)$

 

여기서, S는 표본공분산행렬이다.

 

 

마할라노비스 제곱거리($d_i^2$)는 독립이거나 정확한 카이제곱분포를 따르지는 않으나, 관측값이 다변량 정규분포를 따르는지를 판단하는 측도로 사용될 수 있다. 카이제곱그림(chi-square plot)은 다음의 두 단계로 그려진다. 

 

(1) 순서화된 마할라노비스 제곱거리를 구한다 : $ d_{(1)}^2 \leq d_{(2)}^2 \leq \cdots d_{(n)}^2$

(2) 점 $\left \{ x_p^2 \left(\frac{i-0.5}{n} \right), d_{(i)}^2 \right \}$을 그린다.  -> 카이제곱 분위 수와 순서화된 마할라노비스 제곱거리를 같이 그리는 것으로 qqplot( 카이제곱, d^2)으로 생각 !

 

위 그림이 직선에 가까우면 관측 자료가 다변량 정규분포를 잘 따른다고 판단한다. 또한 직선으로부터 크게 벗어나는 점을 이상치 자료로 판단한다. 아래 그림을 보면, 마할라노비스 제곱거리가 큰 상위 3개 정도의 자료가 다변량 정규분포로부터 벗어난 점으로 보인다. 

 

qqplot을 통해 마할라노비스 제곱거리과 카이제곱 분포 비교.

 

 

마할라노비스 제곱거리는 카이제곱그림과 무관하게 그 자체로 다변량 자료의 이상치를 탐색하는 측도로도 사용될 수 있다. 마할라노비스 제곱거리에 기초한 이상치 탐색을 위한 시각적 방법으로 종유석(stalactite) 그림을 소개한다. 

종유석그림의 원리는 다음과 같다. 

 

(1) 랜덤하게 p+1(=m)개의 자료(행)를 선택한다. 

(2) 선택된 자료를 이용하여 $\mu$와 $\Sigma$를 추정한 뒤, n개 자료에 대해 마할라노비스 제곱거리($d_i^2$)을 구하고, 기준치를 초과하는 자료를 이상치로 간주한다. 

(3) $d_i^2$을 오름차순으로 배열한 뒤, 처음(m+1)개에 대응하는 자료(행)를 선택한다. 

(4) (2)~(3)의 과정을 m=n일 때까지 반복한다. 

 

stalac(USairpollution) :: MVA package in r

위 그림은 R을 사용하여 그린 종유석 그림이다. 자료수가 늘어남에 따라 이상치로 취급되는 자료가 적어짐을 알 수 있다. 41개의 자료를 모두 사용할 때, 이상치로 판단되는 자료는 Chicago, Phoenix, Procidence이다.

 

그 밖에도, aq.plot{mvoutlier}함수가 제공하는 수정된 분위수그림(adjusted quantile plot)으로 이상치를 탐색할 수 있다.수정된 분위수그림은 마할라노비스 제곱거리의 경험적 분포함수에 대한 관측값의 순서화된 로버스트 마할라노비스 제곱거리를 계산하여 다변량이상치를 식별할 수 있게 해준다. 

 

+) 이변량 자료의 이상치 탐색을 위해서는 이변량 상자그림을 이용할 수 있다.  bv.boxplot{asbio}, bvbox{MVA}

 

 

참고 )

R 응용 다변량분석 (저자 : 나종화)

https://m.blog.naver.com/PostView.nhn?blogId=mj93m&logNo=221097578389&proxyReferer=https:%2F%2Fwww.google.com%2F