Loading [MathJax]/jax/output/CommonHTML/jax.js
본문 바로가기

컴퓨터 비전/밑바닥부터 시작하는 딥러닝5

3. 다변량 정규 분포

이전에 배운 정규분포는 하나의 스칼라<키, 무게 등>에 대한 정규분포였고 이번에 배울 것은 여러개의 실수로 이루어진 벡터의 정규분포에 대해 배운다.

3.1 넘파이와 다차원 배열

3.1.1 다차원 배열

다차원 배열은 값<원소> 여러개를 한꺼번에 처리하기 위한 데이터 구조이다. 원소의 배열에는 방향이 있으며 이 방향을 이라고 하고, 축의 개수를 차원이라고 한다. 

 

스칼라: 1

백터 (123)

행렬 (123456)

 

스칼라는 단순 숫자 하나를 뜻하고, 벡터는 하나의 축을 따라 나열되는 형태를 띈다. 행렬은 두개의 축을 따라 나열되는데 가로 방향을 행<row>, 세로 방향을 열<column>이라고 한다. 그리고 코드상에서는 이러한 다차원 배열을 텐서라고 한다. 

벡터는 여기서는 단순 숫자의 나열이지만, 중요한건 축이 2개이기 때문에 열 백터인지, 행 백터인지 구별은 해야한다. 일단 해당 책에서는 모든 벡터는 열 벡터로 취급한다고 한다. 또한 벡터나 행렬은 기호로 표기할 때 수식에서는 볼드체로 표기하여 스칼라와 구분하는 식으로 표현한다. 

차원은 벡터와 행렬에서 다르게 사용된다. 벡터에서는 원소의 갯수가 차원의 수이다. 따라서 벡터에 원소가 3개면 이건 3차원 벡터이다. 하지만 행렬에서는 축의 갯수와 차원의 수가 같다. 따라서 원소가 3개라고 해서 3차원 행렬이 되지는 않는다. 

3.1.2 넘파이 다차원 배열

import numpy as np
x = np.array([1, 2, 3])

print(x.__class__)
print(x.shape)
print(x.ndim)

위 코드를 이용하면 간단하게 넘파이를 사용해서 벡터를 생성하고 이에 대한 정보를 출력할 수 있다.

<class 'numpy.ndarray'> 3, 1

위와 같이 결과가 나온다. 여기서 .shape는  다차원 배열의 형상을, .ndim는 차원의 수를 나타내다. 잘 보면 원소가 3개인 1차원 배열인 것을 알 수 있다. 

import numpy as np
x = np.array([[1, 2, 3], [4, 5, 6]])

print(x.__class__)
print(x.shape)
print(x.ndim)

위와 같이 배열을 선언한다면 다음과 같이 2x3 행렬 즉 2차원 배열임을 알 수 있다. 

<class 'numpy.ndarray'> 2,3 2

3.1.3 원소별 연산

import numpy as np

W = np.array([[1, 2, 3], [4, 5, 6]])
X = np.array([[1, 2, 3], [4, 5, 6]])

print(W+X)
print("---")
print(W-X)

위와 같이 행렬끼리 연산을 한다면 컴공과들이 대부분 배우는 선형대수학에서와 마찬가지로 행렬끼리의 합이 연산되어 원소별 연산이 가능해진다.

[[ 2 4 6] [ 8 10 12]]

[[0 0 0] [0 0 0]]

여기서 차이점은 곱 연산이다. 행렬에서 곱 연산은 dot product, cross product 등 여러 연산이 있는데 아래와 같이 단순 곱하기 연산이라면 아다마르 곱이라고 하여 원소별 곱 연산을 수행한다.

import numpy as np

W = np.array([[1, 2, 3], [4, 5, 6]])
X = np.array([[1, 2, 3], [4, 5, 6]])

print(W*X)

[[ 1 4 9] [16 25 36]]

3.1.4 벡터의 내적과 행렬 곱

벡터의 내적은 각 원소간 곱을 전체 더한 값이다. 따라서 내적을 하기 위해선 두 벡터의 원소의 수가 동일해야 하고 그 결과값은 스칼라 값으로 나온다. 

행렬 곱은 흔히 아는 외적으로 그 연산이 상대적으로 복잡하다. 

위와 같은 방식으로 연산이 진행되고, 따라서 연산 간에 앞 행렬의 행과 뒤 행렬의 열이 같은 개수이어야만 행렬곱이 가능하다. 그리고 일반적인 곱셈과 달리 두 피연산자의 순서에 따라 결과가 달라지므로 연산 시 주의해야 한다.

코드로 하면 다음과 같다. 

import numpy as np

a = np.array([1, 2, 3])
b = np.array([4, 5, 6])
y = np.dot(a, b)
print(y)

A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])
Y = np.dot(A, B)
print(Y)

파이썬에서는 해당 연산을 np.dot 함수로 제공한다. 해당 연산은 인자로 벡터가 들어오면 내적을, 행렬이 들어오면 행렬곱을 수행한다. 따라서 연산 결과는 다음과 같다

32 [[19 22] [43 50]]

연산자를 사용하려면 @ 연산자를 사용하면 된다 

import numpy as np

a = np.array([1, 2, 3])
b = np.array([4, 5, 6])
y = a@b
print(y)

A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])
Y = A@B
print(Y)

3.2 다변량 정규분포

키를 대상으로 정규분포를 생성한다면 키라는 하나의 스칼라 값을 관측값으로 하여 정규분포를 구성하였다. 그 다음으로는 실수 여러개로 구성된 즉 벡터로 구성된 벡터의 정규분포이다. 키와 몸무게로 구성된 2차원 벡터로 가정한다면, 관측값은 <160cm, 48kg>, <175cm, 70kg> 이런 식으로 구성될 것이다. 

2차원 정규분포는 위와 같이 생성된다. 2차원인 이유는 관측값이 <키, 몸무게>와 같이 2차원 벡터 값이기 때문이다. 

3.2.1 다변량 정규분포 공식

다변량 정규분포이기 때문에 여기서 사용할 확률변수는 기본적으로 D차원 벡터를 사용한다. 그리고 다차원 벡터 x에 대한 정규분포식은 다음과 같다.

 

N(x;μ,)=1(2π)D||exp{12(xμ)1(xμ)}

 

위 식에서 나오는 기호에 대해 정리하자면

 

||: 행렬식

(xμ): 전치행렬

σ1: 역행렬

 

다음과 같다. 그리고 위 식에서 사용된 다변량 정규분포의 매개변수인 뮤와 시그마는 각각 평균벡터와 공분산 행렬을 뜻한다. 평균벡터는 확률변수와 동일한 차원의 벡터이고, 공분산 행렬은 DxD 행렬이 된다. 특징으로는 공분산 행렬의 대각선 성분은 각 변수의 분산을, 나머지 성분은 각 변수 사이의 공분산을 나타낸다. 

공분산

먼저 확률변수 x의 분산은 다음 식으로 계산이 가능하다.

$Var[x_{i}]=E[({x_{i}-\mu_{i})^{2}]$

공분산은 분산을 일반화 한 것이다. 이 공분산은 다음과 같은 식으로 계산이 가능하다.

Cov[xi,xj]=E[(xiμi)(xjμj)]

여기서 중요한 점은 i와 j가 같은 경우 즉 공분산 행렬의 대각선 부분은 분산 계산식과 동일하다는 점이다. 이러한 점 때문에 공분산은 분산의 일반화 버전이라고 하는 것이다. 그리고 i와 j가 순서가 바뀌어도 식의 결과값은 동일한 것을 볼 수 있다. 따라서 공분산 행렬은 무조건 대칭행렬이다. 따라서 전치를 하더라도 행렬이 변하지 않는다. 

전치

(xμ)

위 기호는 전치 행렬을 뜻한다. 이 전치는 기존 행렬에서 행과 열을 바꾸면 된다. 1행은 1열이 되고 1열은 1행이 되어 행렬이 새롭게 구성된다. 좀 더 직관적으로 말하면 대각선을 기준으로 선대칭을 한다고 생각하면 쉽다. 

추가적으로 벡터의 내적을 전치를 사용해 표현할 수도 있다. 만약 x, y 두 벡터가 열벡터라면 그 내적은 

xy

위 식으로 표현이 가능하다. numpy에서는 이를 행렬의 속성 중 T로 제공한다. 코드는 다음과 같다. 

import numpy as np

a = np.array([[1, 2, 3], [4, 5, 6]])
print(a)
print("------------------------")
print(a.T)

위 코드를 출력하면 처음은 일반적으로 2x3 행렬이 출력되고, T 속성을 사용하여 출력하면 전치행렬로 되어 3x2 행렬로 출력되는 것을 확인할 수 있다.

[[1 2 3]

[4 5 6]]

------------------------

[[1 4]

[2 5]

[3 6]]

행렬식

||

위 식은 행렬식을 뜻한다. 행렬식이란 정사각 행렬의 특징을 나타내는 하나의 지표로 하나의 스칼라 값으로 계산된다. 계산식 자체는 다음과 같다.

이러다보니 2x2는 사람이 구해도 쉽게 할 수 있는데 3x3부터는 계산이 복잡해지고 4x4는 사실상 사람이 하긴 많이 어렵다<귀찮다...>

다행히 numpy에서 det 함수로 해당 계산식을 제공하기 때문에 굳이 계산을 열심히 할 필요는 없다. 여기서 det은 determinent의 약자이다.

import numpy as np

a = np.array([[1, 2, 3],[4, 5, 6],[7, 8, 9]])
print(np.linalg.det(a))

위처럼 np.linalg <linear algebra, 선형대수>에서 .det 함수를 사용하면 구할 수 있다.

-9.51619735392994e-16

원래 결과는 0인데 아쉽게도 해당 함수는 약간의 오차를 가지고 있다. 따라서 그냥 0에 가깝게 나온 정도로 만족하자...

역행렬

NxN 행렬에서 대각선 성분 값이 전부 1이고 나머지 성분 값이 0인 행렬을 단위 행렬<I>이라고 한다. 그리고 특정 행렬과 행렬곱을 했을 때 그 결과가 단위 행렬이 나올 때 그 두 행렬 관계를 역행렬 관계라고 한다. 수식으로는 다음과 같이 표현한다.

1

구하는 식은 행렬식과 마찬가지로 행렬이 커질수록 어려워지니 그냥 코드로만 알아보자...

import numpy as np

a = np.array([[1, 2, 3],[4, 5, 6],[7, 8, 9]])
print(np.linalg.inv(a))

위처럼 np.linalg에서 inverse의 약자인 inv함수를 사용하면 역함수를 구할 수 있는데, 문제점은 위 식대로 하면 원래는 역함수가 나오면 안된다. 역행렬을 구할 때 행렬식으로 나누는 부분이 있는데 위 행렬은 행렬식이 0이기 때문에 원래는 구할 수 없는 값이다. 그럼 값이 나오면 안되는데, 이것도 결국 부동소수점 계산에 대한 오차가 있어 값이 나와버린다...

[[ 3.15251974e+15 -6.30503948e+15 3.15251974e+15]

[-6.30503948e+15 1.26100790e+16 -6.30503948e+15]

[ 3.15251974e+15 -6.30503948e+15 3.15251974e+15]]

물론 누가봐도 너무 큰 값이기 때문에 잘못나온 값이란 건 쉽게 알 수 있다. 아무튼 제대로 나오는 값으로 계산하면 그 결과는 다음과 같다.

[[-1.79166667 0.91666667 -0.125 ]

[ 1.58333333 -0.83333333 0.25 ]

[-0.125 0.25 -0.125 ]]

이렇게 잘 출력되는 것을 볼 수 있다. 

3.2.2 다변량 정규분포 구현

아까 위에서 봤던 다변량 정규분포 식을 다시 한번 보자

 

N(x;μ,)=12(2π)D||exp12(xμ)1(xμ)

 

그리고 위 수식을 코드화 한 것이 다음 코드이다. 

import numpy as np

def multivariate_normal(x, mu, cov):
    det = np.linalg.det(cov)
    inv = np.linalg.inv(cov)
    D = len(x)
    z = 1/np.sqrt((2*np.pi)**D*det)
    y = z * np.exp((x-mu).T @ inv @ (x-mu)/-2.0)
    return y

먼저 필요한 행렬들을 변수로 미리 계산한 뒤 최종적으로 y값을 반환하는 함수이다. 중요한 점은 여기서 사용된 인자들은 다음 조건을 만족해야 한다.

인자 형상<D는1 이상의 정수>
x <확률변수> 1, D 또는 D, 
mu <평균벡터> 1, D 또는 D, 
cov <공분산 벡터> D, D

x와 mu는 무조건 열벡터이어야 하는데 numpy 특성 상 이번 경우는 D,  형상이어도 계산 자체는 문제가 없다. 이제 위 선언한 함수를 이용하여 2차원 벡터를 이용하여 확률 밀도를 구해보자.

import numpy as np

def multivariate_normal(x, mu, cov):
    det = np.linalg.det(cov)
    inv = np.linalg.inv(cov)
    D = len(x)
    z = 1/np.sqrt((2*np.pi)**D*det)
    y = z * np.exp((x-mu).T @ inv @ (x-mu)/-2.0)
    return y

x = np.array([[0], [0]])
mu = np.array([[1], [2]])
cov = np.array([[1, 0],
                [0, 1]])

y = multivariate_normal(x, mu, cov)
print(y)

 

위 코드를 실행하면 하나의 확률 밀도 값이 떨어진다.

[[0.01306423]]

행렬의 형태를 띄긴 하지만, 사실상 하나의 값 즉 스칼라 값인 것을 확인할 수 있다. 행렬의 형태로 나오는건 피연산자들이 np.array 타입이기 때문이다.

3.3 2차원 정규분포 시각화

3.3.1 3D 그래프 그리기

matplotlib 라이브러리를 이용하면 3D로 그래프를 그릴 수 있다. 

import numpy as np
import matplotlib.pyplot as plt

X = np.array([[-2, -1, 0, 1, 2],
              [-2, -1, 0, 1, 2],
              [-2, -1, 0, 1, 2],
              [-2, -1, 0, 1, 2],
              [-2, -1, 0, 1, 2]])

Y = np.array([[-2, -2, -2, -2, -2],
              [-1, -1, -1, -1, -1],
              [0, 0, 0, 0, 0],
              [1, 1, 1, 1, 1],
              [2, 2, 2, 2, 2]])

Z = X ** 2 + Y ** 2

ax = plt.axes(projection='3d')
ax.plot_surface(X, Y, Z, cmap='viridis')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
plt.show()

위 코드를 사용하면 다음 결과를 확인할 수 있다. 

5x5인 2차원 행렬을 3차원의 그래프로 표현하면 다음과 같다. 3차원으로 되는건 2차원의 행렬이 값을 가지기 때문에 이를 표현하려면 3차원이 필요하기 때문이다. 이건 값이 이산적이기 때문에 각진 그래프 모양을 가지지만, np.meshgrid 함수를 사용하면 부드러운 곡선을 표현할 수 있다. 

import numpy as np
import matplotlib.pyplot as plt

xs = np.arange(-2, 2, 0.1)
ys = np.arange(-2, 2, 0.1)

X, Y = np.meshgrid(xs, ys)
Z = X ** 2 + Y ** 2

ax = plt.axes(projection='3d')
ax.plot_surface(X, Y, Z, cmap="viridis")
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
plt.show()

위 코드는 행렬을 선언하는 부분만 5x5를 정적으로 선언하는 것을 meshgrid 함수로 변경한 코드이다. 그 결과는 다음과 같다.

해당 함수는 1차원 배열에 대해 즉 xs, ys 배열에 대해 3차원 그래프에 대응하는 격자점용 2차원 배열이 생성된다. 

[[-2. -1.9 -1.8 ... 1.7 1.8 1.9]

[-2. -1.9 -1.8 ... 1.7 1.8 1.9]

[-2. -1.9 -1.8 ... 1.7 1.8 1.9]

...

[-2. -1.9 -1.8 ... 1.7 1.8 1.9]

[-2. -1.9 -1.8 ... 1.7 1.8 1.9]

[-2. -1.9 -1.8 ... 1.7 1.8 1.9]]

 

[[-2. -2. -2. ... -2. -2. -2. ]

[-1.9 -1.9 -1.9 ... -1.9 -1.9 -1.9]

[-1.8 -1.8 -1.8 ... -1.8 -1.8 -1.8]

...

[ 1.7 1.7 1.7 ... 1.7 1.7 1.7]

[ 1.8 1.8 1.8 ... 1.8 1.8 1.8]

[ 1.9 1.9 1.9 ... 1.9 1.9 1.9]]

이런 식으로 생성된다.

3.3.2 등고선 그리기

앞에서 plot_surface 함수를 사용했는데 이건 곡면을 그리는 함수이다. 이와 달리 등고선을 그리려면 contour 함수를 사용하면 된다.

import numpy as np
import matplotlib.pyplot as plt

x = np.arange(-2, 2, 0.1)
y = np.arange(-2, 2, 0.1)
X, Y = np.meshgrid(x, y)
Z = X ** 2 + Y ** 2

ax = plt.axes()
ax.contour(X, Y, Z)
ax.set_xlabel('x')
ax.set_ylabel('y')
plt.show()

3.3.3 2차원 정규분포 그래프

먼저 다변량 정규분포 식을 다시 한번 확인해보자

 

N(x;μ,)=1(2π)D||exp(12(xμ)1(xμ))

 

이걸 코드화 하면 다음과 같다.

import numpy as np
import matplotlib.pyplot as plt

def multivariate(x, mu, cov):
    det = np.linalg.det(cov)
    inv = np.linalg.inv(cov)
    D = len(x)
    
    z = 1/np.sqrt((2 * np.pi) ** D * det)
    y = z * np.exp((x - mu).T @ inv @ (x - mu) / -2.0)
    return y

mu = np.array([0.5, -0.2])
cov = np.array([[2.0, 0.3], [0.3, 0.5]])

xs = ys = np.arange(-5, 5, 0.1)
X, Y = np.meshgrid(xs, ys)
Z = np.zeros_like(X)

for i in range(X.shape[0]):
    for j in range(X.shape[1]):
        x = np.array([X[i, j], Y[i, j]])
        Z[i, j] = multivariate(x, mu, cov)

fig = plt.figure()

ax1 = fig.add_subplot(1, 2, 1, projection="3d")
ax1.set_xlabel('x')
ax1.set_ylabel('y')
ax1.set_zlabel('z')
ax1.plot_surface(X, Y, Z, cmap = 'viridis')

ax2 = fig.add_subplot(1, 2, 2)
ax2.set_xlabel('x')
ax2.set_ylabel('y')
ax2.contour(X, Y, Z)

plt.show()

결과는 다음과 같다. 간단하게 코드를 보자면 fig,subplot 함수의 경우 하나의 결과 즉 캔버스에 여러개의 그래프를 그릴 수 있도록 해주는 함수로 인자 순서는 행, 열, 그래프 영역번호 순서이다. 따라서 1, 2, 1 순서로 전달한다면 이것은 1행 2열의 왼쪽 영역에 그래프를 그리라는 의미이다. 

mu의 경우 평균벡터로 등고선의 중앙, 산봉우리의 정상이 위치하는 자리이다. 그리고 cov의 경우 공분산으로 산의 넓이 또는 퍼짐 정도를 결정하는 요소이다. 행렬의 원소 중 대각선 행렬 [0, 1], [1, 0]는 공분산 즉 x, y 축 두 방향 모두에 대한 퍼짐 정도를 나타내고 [0, 0]은 x축, [1, 1]은 y축에 대해 퍼짐 정도를 나타낸다. 

3.4 다변량 정규분포의 최대 가능도 추정 

3.4.1 최대 가능도 추정하기

일단 먼저 다변량 정규분포의 정규분포식부터 다시 확인하자. 

 

N(x;μ,)=1(2π)D||exp(xμ)1(xμ)

 

다음으로 실제 데이터 샘플D를 얻은 경우에 대한 확률밀도는 각 샘플에 대한 확률을 전부 곱하면 되기 때문에 그 식은 다음과 같다.

 

p(D,μ,)=N(x(1);μ,)N(x(2);μ,)N(x(3);μ,)N(x(N);μ,)=Nn=1N(x(n);μ,)

 

위 식은 매개변수를 뮤와 시그마로 했을 때 샘플 D를 얻을 수 있는 확률을 의미한다. 그리고 위 식의 결과를 가능도라고 한다. 우리가 하고자 하는건 최대 가능도 추정으로 이 말은 가능도를 최대로 해주는 매개변수 즉 뮤와 시그마 값을 찾는 과정이다. 여기서 가능도 계산은 이전 장에서 일반 정규분포에서와 마찬가지로 로그를 취해서 계산을 진행한다. 이제 이전장과 마찬가지로 미분을 진행하는데 

 

Lμ=0,L=0

 

이렇게 구할 것이다. 근데 문제는 이전에 했던 일반 정규분포와 달리 이번에는 L은 스칼라지만, 뮤와 시그마는 벡터란 점이다. 벡터, 행렬에서의 미분은 각 원소들을 미분하면 된다. 

대충 이런식인데 구체적인 계산 과정은 걍 생략하고 결과만 보자면 다음 식이 떨어진다.

 

ˆμ=1NNn=1x(n),ˆ=1NNn=1(x(n)ˆμ)(x(n)ˆμ)

 

이렇게 식이 나온다. 이제 저 두 값이 0이 되는 지점을 찾으면 된다.위 식에서 나오는 벡터의 형상을 간단하게 확인하자면 다음 표와 같다.

변수 형상
x(n) Dx1
x(n)μ Dx1
(x(n)μ) 1xD
(x(n)ˆμ)(x(n)ˆμ) DxD

3.4.2 최대 가능도 추정 구현

위 식을 코드로 구현하면 다음과 같다.

import numpy as np

np.random.seed(0)

N = 10000
D = 2
xs = np.random.rand(N, D)

mu = np.sum(xs, axis=0)/N
cov = 0
for n in range(N):
    x = xs[n]
    z = x - mu
    z = z[:, np.newaxis]
    cov += z @ z.T
    
cov /= N
print(mu)
print(cov)

간단히 설명하면 np.sum에서 axis=0의 의미는 0번째 축을 따라 sum연산을 하라는 의미이고 z.[:, np.newaxis]는 새로운 축을 만들라는 의미로 기존 z의 형상은 <D, >이지만 <D, 1>로 즉 2차원 열벡터가 추가된다. 마지막으로 cov += z @ z.T에서 z 행렬을 전치한 이유는 위에서 공분산 식은 

Cov[xi,xj]=E[(xiμi)(xjμj)]

다음과 같은데 잘 보면 E에서 앞의 행렬은 i를 기준으로 하고 뒤 행렬은 j를 기준으로 한다. 따라서 전치를 취한 뒤 반복문으로 행렬곱을 수행하면 쉽게 계산이 가능하다. 그리고 그 결과는 다음과 같다

[0.49443495 0.49726356]

[[ 0.08476319 -0.00023128]

[-0.00023128 0.08394656]]

근데 매번 수식을 쓰는건 너무 귀찮기도 하고 사실 numpy에서 함수를 제공하고 있기 때문에 그냥 함수를 가져다 쓰면 된다.

import numpy as np

np.random.seed(0)

N = 10000
D = 2
xs = np.random.rand(N, D)

mu = np.mean(xs, axis=0)
cov = np.cov(xs, rowvar=False)

print(mu)
print(cov)

여기서 cov 함수에서 인자 중 rowvar는 각 을 하나의 데이터로 취급하여 공분산 행렬을 진행한다는 의미이다. 기본값은 True로 이렇게 되면 각 을 하나의 데이터로 취급하여 연산을 진행한다. 

3.4.3 실제 데이터 사용

2챕터에서 사용한 키 데이터를 사용할 건데 이번에는 키와 몸무게가 쌍으로 이루어진 데이터를 사용한다.

import numpy as np
import os
import matplotlib.pyplot as plt

path = os.path.join(os.getcwd(), "height_weight_chap3.txt")
xs = np.loadtxt(path)
small_xs = xs[:500]
plt.scatter(small_xs[:, 0], small_xs[:, 1])
plt.xlabel("Height(cm)")
plt.ylabel("Weight(kg)")
plt.show()

일단 위 코드를 통해 2차원 평면에서 해당 데이터의 분포를 시각화한다.

이제 이 데이터를 위에서 봤던 3D그래프와 등고선으로 표현하면 다음과 같은 코드를 이용하면 된다. 해당 코드는 위 코드를 거의 그대로 쓰고 조금만 변형한 것이다.

import numpy as np
import matplotlib.pyplot as plt

def multivariate(x, mu, cov):
    det = np.linalg.det(cov)
    inv = np.linalg.inv(cov)
    D = len(x)
    
    z = 1/np.sqrt((2 * np.pi) ** D * det)
    y = z * np.exp((x - mu).T @ inv @ (x - mu) / -2.0)
    return y

mu = np.array([0.5, -0.2])
cov = np.array([[2.0, 0.3], [0.3, 0.5]])

xs = ys = np.arange(-5, 5, 0.1)
X, Y = np.meshgrid(xs, ys)
Z = np.zeros_like(X)

for i in range(X.shape[0]):
    for j in range(X.shape[1]):
        x = np.array([X[i, j], Y[i, j]])
        Z[i, j] = multivariate(x, mu, cov)

fig = plt.figure()

ax1 = fig.add_subplot(1, 2, 1, projection="3d")
ax1.set_xlabel('x')
ax1.set_ylabel('y')
ax1.set_zlabel('z')
ax1.plot_surface(X, Y, Z, cmap = 'viridis')

ax2 = fig.add_subplot(1, 2, 2)
ax2.set_xlabel('x')
ax2.set_ylabel('y')
ax2.contour(X, Y, Z)

plt.show()

달라진 부분은 기본적으로 랜덤 데이터가 아닌 실제 데이터를 사용하기에 생겨난 xs 선언 부분 그리고 실제 데이터에 맞게 meshgrid를 데이터 범위에 맞게 적절히 선언하는 부분, 마지막으로 등고선 뿐만 아니라 등고선에 실제 데이터를 뿌려줄 수 있도록 ax2.scatter 함수를 사용한 부분 이 3가지 이외에는 전부 동일한 코드이다. 결과는 다음과 같다.

결론

이번 챕터에서는 기존 정규분포 모델링에서 더 나아가 다차원 데이터에 대한 다변량 정규분포에 대한 모델링까지 다루어 보았다.

'컴퓨터 비전 > 밑바닥부터 시작하는 딥러닝5' 카테고리의 다른 글

5. EM 알고리즘  1 2025.02.21
4. 가우스 혼합 모델  0 2025.02.08
2. 최대 가능도 추정  0 2025.01.06
1장. 확률분포  0 2025.01.03