공분산과 공분산 행렬

어떤 현상이나 타깃을 관찰하다 보면 여러 특징(features)들이 있습니다. 예를 들어

• 사람을 관찰할 경우 키, 몸무게, 나이 등등
• 학생 학습 성취도를 관찰할 경우: 국어 점수, 수학 점수 등

그런데 경험적으로 키가 큰 사람이 몸무게도 많이 나가고,  수학 점수가 뛰어난 학생이 과학 점수도 좋은, 이른바 양의 상관관계를 띄는 경우도 있고,  같은 거리를 가는데 느리게 가면 시간이 많이 걸리는 등 속도와 시간의 음의 상관관계도 있고, 또 아예 관계가 없어 보이는 경우도 있습니다. 

 

우린 누구누구랑 상관관계가 아주 쩔어주는 사람들~

 

두 특징(또는 변량, 변수)의 관계를 따지는 유명한 방법이 공분산을 산출해 보는 것입니다.

 

두 변수 $X, Y$의 공분산은 아래와 같이 정의합니다.

$X,Y$의 공분산(covariance)

$X, Y$의 표본(sample, 즉 관찰 데이터) $n$개가 각각
$x_1, x_2, \cdots , x_n$, $y_1, y_2, \cdots, y_n$ 일 때, 공분산은
$$ \mathbb{Cov}(X,Y) = \frac 1n \sum_{i=1}^n (x_i-\mu_x)(y_i-\mu_y)$$

여기서
$\mu_x, \mu_y$ 는 각각 $X, Y$의 평균입니다. 즉,
$$\mu_x = \frac 1n\sum_{i=1}^n x_i ~,~ \mu_y = \frac1n \sum_{i=1}^n y_i $$
입니다.

기댓값 기호 $\mathbb{E}(\cdot)$를 사용해서 표시해보면 

$$\mathbb{Cov}(X,Y) = \mathbb{E}((X-\mu_X ) (Y- \mu_Y)) =\mathbb{E}(XY)-\mu_x \mu_y =\mathbb{E}(XY)-\mathbb{E}(X)\mathbb{E}(Y)$$

라고도 쓸 수 있습니다.

 

그럼 $X$와 $X$ 의 공분산은 뭘까요? 

$$\mathbb{Cov}(X,X) = \frac 1n \sum_{i=1}^n (x_i-\mu_x)^2 = \mathbb{V}(X)$$

바로 $X$의 분산이 됩니다.

 

---

공분산 행렬(covariance matrix)은 공분산을 원소로 하는 행렬입니다. 

features( 즉, 통계 변량)이 $p$개 있고,  $X_1, X_2, \cdots, X_p$ 라 합시다. 관찰 데이터는 $n$개의 sample이 있다 하고요. 이때 공분산 행렬 $\mathbf{C}$는 $p\times p$ 행렬로써 $(i,j)$ 원소가

 $$\mathbf{C}_{ij} = \mathbb{E}(X_i X_j) - \mathbb{E}(X_i)\mathbb{E}(X_j) $$

인 행렬입니다.

만일 

$X_i$ 의 샘플이 $x_{1i}, x_{2i}, \cdots, x_{ni}$ 이고 $X_j$의 샘플이 $x_{1j}, x_{2j},\cdots, x_{nj}$ 라 하면

$$\mathbf{C}_{ij} = \frac 1n \sum_{r=1}^n x_{ri}x_{rj} - \frac1{n^2} \sum_{r=1}^n x_{ri} \sum_{r=1}^n x_{rj} $$

또는

 

$$\mathbf{C}_{ij} = \frac 1n \sum_{r=1}^n (x_{ri}-\mu_i)\cdot (x_{rj} - \mu_j)$$

입니다. 첨자 때문에 힘들어 보여도 별 거 covariance의 공식을 떠올리면 별 거 아닙니다.

 

---

 

그렇다며 공분산 행렬을 쉽게 구하는 방법은 뭘까요?

 

$$ X = \pmatrix{
           \mid & \vert &  & \vert \\
           X_1-\mu_1 & X_2 -\mu_2 & \cdots & X_p-\mu_p \\
           \vert & \vert &  & \vert  } $$

라 정의합시다.  $\mu_i$는 $i$ 번째 변수의 평균(기댓값)입니다.

 

즉 일례로, $X$의 첫 번째 열은

$$X_1 -\mu_{1} = \pmatrix{ x_{11}-\mu_1 \\ x_{21}-\mu_1\\ \vdots \\ x_{n1}-\mu_1 }$$

입니다. 이때,

 

$$ \frac1n X^T \cdot X = \mathbf{C}\tag{1} $$

가 됩니다. 일견 당연해 보이지만 살짝 산수를 해 본다면,

$$
\begin{align}
\frac1n(\mathbf{X}^T \mathbf{X})_{ij}     & = \frac1n \sum_{r=1}^n (\mathbf{X}^T)_{ir} \mathbf{X}_{rj} \\
             & =\frac1n \sum_{r=1}^n \mathbf{X}_{ri} \mathbf{X}_{rj} \\
           & =\frac1n \sum_{r=1}^n (x_{ri}-\mu_i)(x_{rj}-\mu_j) = \mathbf{C}_{ij}             
\end{align}
$$

이 성립하죠.

 

 

---

python으로 공분산 행렬을 구해볼까요?

2022.05.13 - [엑셀] - 트렌드 직선의 비밀(선형회귀) 엑셀실습 #1에서 다루었던 예제를 가지고 분석해 보겠습니다.

어느 배달업체의 빅??데이터

예전 선형회귀 관련 글에서 다루었던 배달업체의 배달 거리와 배달시간 조사내역입니다. 배달 15건을 대상으로 관찰한 데이터입니다.

 

두 가지 방법으로 구해 서로 비교해 보겠습니다. 두 가지 방법은

1.  수식(1)처럼 공분산의 정의대로 구해보기
2. python numpy 안의 cov 함수 사용

입니다.

 

code를 보시죠.

 

import numpy as np

def calculate_cov_matrix():
    data = [[2, 20],
            [3, 24],
            [6, 36],
            [8, 47],
            [2.6, 22],
            [4.8, 32],
            [7, 47],
            [4, 42],
            [2, 21],
            [2.2, 21],
            [3.8, 30],
            [2.4, 25],
            [2.6, 18],
            [5.4, 38],
            [5.1, 28]
            ]
    nData = len(data)

    x = np.array(data)
    x_mean = x
    x_mean[:, 0] = x[:, 0] - np.mean(x[:, 0])
    x_mean[:, 1] = x[:, 1] - np.mean(x[:, 1])

    print('mean adjust of data')
    print(x_mean)

    res = x_mean.T @ x_mean / nData

    print('\n')
    print('covariance matrix(from definition):')
    print(res)

    c = np.cov(x.T, bias=True)
    print('\n')
    print('covariance matrix(python function):')
    print(c)

if __name__ == '__main__':
    calculate_cov_matrix()

 

code를 간략히 보시죠.

 

 

 

    nData = len(data)

○ data의 개수 즉, 15를 뜻합니다.

 

 

 

    x = np.array(data)
    x_mean = x
    x_mean[:, 0] = x[:, 0] - np.mean(x[:, 0])
    x_mean[:, 1] = x[:, 1] - np.mean(x[:, 1])

○ 기본적인 numpy 함수들을 쓰기 위해 data를 numpy의 array 형으로 변환합니다. (np.array)

○ x의 각 열에서 각 열의 평균을 빼 준 행렬을 구하기 위해 x와 똑같은 복사본 x_mean을 준비하고

○ 0 열과 1열 각각 의 평균을 각 열에서 빼줍니다.

 

 

 

    res = x_mean.T @ x_mean / nData

○ X의 transpose와 X의 곱을 구해서 data 개수로 나눠줍니다. 이것이 수식(1)입니다!

 

 

 

    c = np.cov(x.T, bias=True)

○ numpy에서 제공하는 cov 함수를 사용합니다. 그런데 주의할 점은 

 

numpy의 cov의 input은 행이 통계 변수(features)이고 열이 관찰 데이터가 되게끔 넣어야 한다는 것,

(우리의 $X$와 행/열이 거꾸로)

 

입니다. 따라서 우리가 위의 code에서 만든 x를 넣으면 안 되고, x.T 를 넣어줘야 한다는 것!!

 

○ bias를 허용합니다. 즉 표본 공분산은 "관찰 데이터의 개수 - 1 = 14 "로 나눠줘야 하나, 관찰데이터의 개수인 15로 나눠주라는 의미입니다. 2022.05.09 - [엑셀] - STDEV 함수에서 다룬 내용과 일맥상통하는 내용입니다.

 

○ 좀 더 자세한 사항은 numpy의 tutorial를 참고하시기 바랍니다.

 

 

 

이제 결과를 한번 보겠습니다.

 

 

 

mean adjust of data
[[ -2.06       -10.06666667]
 [ -1.06        -6.06666667]
 [  1.94         5.93333333]
 [  3.94        16.93333333]
 [ -1.46        -8.06666667]
 [  0.74         1.93333333]
 [  2.94        16.93333333]
 [ -0.06        11.93333333]
 [ -2.06        -9.06666667]
 [ -1.86        -9.06666667]
 [ -0.26        -0.06666667]
 [ -1.66        -5.06666667]
 [ -1.46       -12.06666667]
 [  1.34         7.93333333]
 [  1.04        -2.06666667]]


covariance matrix(from definition):
[[ 3.4344     15.84933333]
 [15.84933333 90.72888889]]


covariance matrix(python function):
[[ 3.4344     15.84933333]
 [15.84933333 90.72888889]]

Process finished with exit code 0

 

직접 절차대로 구해 본 것과 numpy의 cov 함수를 이용한 것이 똑같죠. 

참고로, 만일 unbiased 된 공분산을 구하고 싶을 땐,

 

    c = np.cov(x.T, bias=False)

로 해주면 됩니다. 그러면 결과는

 

 

covariance matrix(python function):
[[ 3.67971429 16.98142857]
 [16.98142857 97.20952381]]

으로  기존과 다르게 나옵니다.

 

이런 분석은 엑셀을 통해서도 많이 이루어지는데요. 엑셀 시트에서는 어떠한 방법으로 공분산 행렬을 구하는지 다음 글에서 알아보도록 하겠습니다.