[기초통계 정리 2] - Probability, Likelihood, MLE와 Python 구현

Deep Learning까지 가지 않더라도, 대부분의 머신 러닝의 개념에서는 확률과 우도의 개념이 빈번하게 등장한다. 그래서 기초 통계를 공부하던 도중, Probability와 Likelihood의 개념을 조금 더 자세하고 천천히 살펴보았다. 

본 포스팅은 확률과 우도에 대한 기초적인 개념을 공부하고, 이를 간단한 Python 코드로 짜본 것이다.

---

1. Random Variable : 확률 변수

- 확률 변수의 정의

일반적으로, 변수라고 함은 특정한 하나의 숫자를 대표하는 것이다. 반면 확률 변수는 나올 수 있는 값들이 확률적 분포를 가지는 것이다. 머신 러닝에서의 확률 변수는 특정한 Vector, 혹은 숫자를 생성하는 기계에 비유할 수 있다. 확률은 0에서 1사이의 값으로 표현하지만, 확률 변수는 실수의 범위 내에서 모든 값으로 표현이 가능하다. 또한 확률 변수를 사용하면 모든 표본은 실수의 집합 즉 수직선으로 표시할 수 있다. 표본의 집합인 사건(event)은 이 수직선 상의 숫자의 집합인 구간(interval)으로 표시된다.

- 이산 확률 변수 (Discrete Random Variable)

만약 확률 변수의 값이 학점과 같이 (A,B,C,D,E,F) 의 형태(discrete)로 구성되어 있다면, 이산 확률 변수라고 한다. 착각하기 쉬운 개념은, 이산 확률 변수도 무한할 수 있다는 것이다(변수간의 유한성이 존재할 뿐 변수 자체는 무한성이 가능).

- 연속 확률 변수 (Continuous Random Variable)

이산 확률 변수와 달리, 확률 변수 사이간에 나올 수 있는 변수가 무한한 경우를 말한다. 연속 확률 변수에서는 사건을 interval 단위로 사용한다. 예를 들어 {a < X < b} 라는 식으로 사건을 나타낼 수 있다.

그리고 사건 A가 표본 k의 집합이라고 할 때, 위와 같이 표현이 가능하다.

---

2. Likelihood : 가능도, 우도

0부터 100까지를 수직선의 구간으로 가지는 연속 확률 변수에서, 특정 숫자 50을 뽑을 확률은 100분의 1이 아니다. 정확히는 1/∞ 이므로, 0이다. 0부터 100사이에 있는 어떠한 실수도 0의 확률을 가진다. 그래서 특정 구간(interval)에 속할 확률로써 확률 변수를 표현한다고 했었다. 

0부터 50까지의 숫자를 뽑을 확률을 생각해보자. 이때의 확률은 1/2이다. 이처럼 특정 구간의 확률을 구하는 함수는 PDF(Probability Density Function)이다. 요약하면 PDF는 특정 구간에 속할 확률을 구해주는 함수를 의미한다.

위와 같은 그림에서의 함수가 바로 PDF이다. 그렇다면 여기에서 의문이 생긴다. 연속 사건의 경우에는 구간에 대한 추정은 가능하지만, 특정 사건에 대한 추정은 불가능한가? 사건에 대한 확률을 구하려면 구간으로 때려맞히는 것 밖에는 없는가?

여기서 생겨난 개념이 바로 Likelihood이다. PDF에서의 y값을 사건이 일어날 가능도로 보겠다는 것이다. 이산 확률 변수에서는 특정 사건이 일어날 확률 = 가능도가 되며, 연속 확률 변수에서는 Likelihood가 확률이 아니라, PDF의 y값이 Likelihood가 된다. 수식적으로, 가능도라는 것은 y=f(x)로 표현할 수 있는 수식이다.

정규분포를 따르는 PDF를 생각해보자. 0부터 100까지의 수직선에 있는 변수 중 50을 뽑을 확률은 원래는 0이지만 Likelihood에서는 0.4이다. 여러 가지 사건을 뽑을 경우도 생각해보자. 세 번을 뽑을 때 50, 60, 70이 나올 확률은 0 X 0 X 0 으로 0이다. 하지만 Likelihood 에서는 0.4 X 0.3 X 0.05 = 0.006이다. 즉 50, 60, 70이 나올 가능도는 0.006으로 계산한다는 것이다.

---

3. MLE (Maximum Likelihood Estimator) : 최대 가능도 추정

MLE는 측정값 혹은 관찰값들에 대한 최대 가능도를 추정하는 것으로, 통계적인 직관을 어려운 말로 풀어서 쓴 것에 불과하다. 간단하게는 가능도를 최대로 하는 확률 추정이라고 할 수 있다. 예시로써 살펴보자. 몸무게의 정규분포를 나타내는 PDF의 y값이 다음과 같다고 하자.

<그림 3>

그리고 어떤 사람의 몸무게를 5번 측정한다(몸무게 역시 연속 확률 변수이다. 지금까지의 내용을 연관하여 잘 생각해보자). 각각 64, 64.5, 65, 65.5, 66의 측정값이 나왔다고 할 때, 이러한 측정값이 나올 MLE를 구하는 예시이다. 직관적으로 예측해본다면, 이 사람의 몸무게의 MLE는 65일 것 같다. 

이제부터 위에서 본 Likelihood에 대한 예시대로 계산해보면, 각각이 (64, 64.5, 65, 65.5, 66)으로 나올 Likelihood는 각각의 Likelihood를 모두 곱한 것이다. 그리고 그 결과값 L을 최대로 하는 μ의 MLE를 계산하는 것이다. 과정은 아래와 같다.

- 몸무게의 모집단은 참값은 평균 65(m), 분산 25(a^2)를 따르는 정규분포라고 가정한다.

- 측정값이 x일때의 가능도, 즉 PDF의 y는 위의 <그림 3>의 수식과 같다.

- (64, 64.5, 65, 65.5, 66) 을 각각 <그림 3>의 수식에 대입하고, 그 결과를 모두 곱한 것이 Likelihood(L)이다.

- L이 최대가 되는 estimator(μ)를 구하는 것이 MLE이다. 풀어서 설명하면, 사건의 가능도가 최대가 되는 추정 값을 구하는 것.

- 이를 Naive하게 Gradient Descent 방식으로 구한다.

- 결과적으로 estimator는 65로 수렴한다.

여기서 잘 생각해보면, μ의 MLE를 구하는 과정은 회귀분석에서 Gradient Descent로 파라미터를 추정하는 것과 동일하다는 것을 알 수 있다. 그래서 딥러닝을 비롯한 각종 ML, 통계학에서 MLE의 개념이 중요하다는 것을 알 수 있다. 이 과정을 Python으로 구현해보자.

---

4. Python 코드로 MLE 구현하기

위의 MLE계산 과정을 그대로 파이썬 코드로 구현하였다.

Github 링크

import numpy as np
import math


class MLE():
    def __init__(self, samples, m, std, learning_rate, epochs, verbose=False):
        """
        :param samples: samples for get MLE
        :param learning_rate: alpha on weight update
        :param epochs: training epochs
        :param verbose: print status
        """
        self._samples = samples
        self._m = m
        self._std = std
        self._learning_rate = learning_rate
        self._epochs = epochs
        self._verbose = verbose


    def likelihood(self, x, M):
        """
        Probability Density Function is Normal distribution
        PDF's y is same as likelihood

        :param x:
        :return: likelihood of input x (likelihood of input x is same as y of pdf)
        """
        return (1 / math.sqrt(2*math.pi) * math.pow(self._std, 2)) * np.exp(-(np.power(x - M, 2) / (2*math.pow(self._std, 2))))


    def fit(self):
        """
        training estimator
        M, which minimizes Likelihood, is obtained by the gradient descent method.
        M is the MLE of the samples
        """

        # init M
        self._estimator = np.random.normal(self._m, self._std, 1)

        # train while epochs
        self._training_process = []
        for epoch in range(self._epochs):
            likelihood = np.prod(self.likelihood(self._samples, self._m))
            prediction = np.prod(self.likelihood(self._samples, self._estimator))
            cost = self.cost(likelihood, prediction)
            self._training_process.append((epoch, cost))
            self.update(self._samples, self._estimator)

            # print status
            if self._verbose == True and ((epoch + 1) % 10 == 0):
                print("Iteration: %d ; cost = %.4f" % (epoch + 1, cost))


    def cost(self, likelihood, prediction):
        """
        cost function
        :param likelihood: likelihood of population
        :param prediction: likelihood in samples
        :return: the cost of optimizing the parameters
        """
        return math.sqrt(likelihood - prediction)


    def update(self, x, M):
        """
        update in gradient descent
        gradient is approximated
        :param x: samples
        :param M: estimator
        """
        gradient = np.sum(np.exp(-(np.power(x - M, 2) / (2*math.pow(self._std, 2)))))
        if self._m > self._estimator:
            self._estimator += self._learning_rate * gradient
        else:
            self._estimator -= self._learning_rate * gradient


    def get_mle(self):
        """
        parameter getter
        :return: estimator of MLE
        """
        return self._estimator


# run example
if __name__ == "__main__":

    # samples for MLE
    samples = np.array([64, 64.5, 65, 65.5, 66])

    # assumptions about the population
    mean = np.array([65.0])
    std = 5

    # get MLE
    estimator = MLE(samples, mean, std, learning_rate=0.1, epochs=30, verbose=True)
    estimator.fit()
    result = estimator.get_mle()
    print(result)

---

참고 자료 :

1. http://rstudio-pubs-static.s3.amazonaws.com/204928_c2d6c62565b74a4987e935f756badfba.html

2. https://datascienceschool.net/view-notebook/4bcfe70a64de40ec945639236b0e911d/