포아송 분포

만약 우리가 금은방을 운영하고 있다고 가정해보겠습니다. 퇴근 후 한 시간에 도둑이이 10명 올 확률은 어떻게 될까요? 100명이 올 확률은 어떻게 될까요?(도둑이 오면 안되겠지만요...) 이와 같이 일정 기간 동안에 확률이 낮은 특정 사건이 일어날 확률을 나타내기 위해 활용하는 것이 포아송 분포입니다. 

다시 위의 예를 좀 더 깊게 들여보다면 저 확률을 이항분포로 나타낼 수 있지 않을까? 하는 생각이 드실 수도 있습니다. 다시 말하면 1분에 도둑이 올 확률이 0.01 오지 않을 확률이 0.99라면 이는 결과가 두 가지 뿐인 베르누이 시행으로 간주할 수 있습니다. 따라서 한 시간(60분)은 베르누이 시행을 60번 시행했다고 볼 수 있을 것입니다.

하지만 도둑이 1분에 한명만 오는 것이라고 한정할 수 없습니다. 즉 1초에 한명만 오는 시행으로 볼 수도 있고 0.01초에 한명만 오는 시행으로 볼 수도 있을 것입니다. 즉 베르누이 시행의 횟수를 무한대로 나타낼 수 있다는 것입니다.

따라서 포아송 분포는 베르누이 시행을 무한대로 시행한 경우로 나타내며, 아래와 같이 정의 할 수 있습니다.

이 때 람다(λ)는 확률변수 X(도둑의 수)의 평균값이며, 도둑의 예를 들면 한시간동안 오는 평균적인 도둑 수입니다.

(** 만약 성공할 확률 p가 0에 가깝지 않고, n이 무한대는 아니나 충분히 큰 수이면, 베르누이 시행을 n번 시행했을 때 성공 횟수는 정규분포에 가깝게 근사가 됩니다. 포아송 분포는 베르누이 시행해서 p가 아주 작고 n이 무한대로 표현될 수 있을때,  근사된 분포입니다.)

다시 정리하면, 포아송 분포는 매 순간 사건 발생이 가능하나, 순간의 발생 확률이 아주 작은 경우에 활용합니다. 따라서 포아송 분포는 연속적인 시간에서 매 순간의 발생확률 대신, 특정 기간에 발생하리라 생각되는 평균 발생 횟수를 이용하여 실제로 발생하는 사건의 횟수의 관한 문제를 다루는 모형입니다.

실제로 포아송 분포를 적용하기 위해서는 다음과 같은 세 가지의 가정을 만족해야 합니다.

1) 사건의 평균 발생횟수는 구간의 길이에만 영향을 받는다. 2) 한 순간에 2회 이상의 사건이 발생할 확률은 거의 0에 가깝다. 3) 한 구간에서 발생한 사건은 다른 구간에서 발생한 사건에 영향을 주지 않는다.

Python으로 포아송 분포를 시뮬레이션 해보겠습니다.

## 필요 Module Import
from scipy.stats import poisson
import matplotlib.pyplot as plt
import numpy as np
from collections import Counter
import pandas as pd

### 평균 3가정 하여 분포 추정
poisson.stats(mu, moments ="mvsk") ##mean, var, skewness, kurtosis 확인 

## 추정한 분포 시각화
fig, ax = plt.subplots(1, 1)
x = np.arange(poisson.ppf(0.01, mu), poisson.ppf(0.99, mu))## ppf -> 누적 분포 함수 확률이 0.01~0.99사이인 x 값들 반환
ax.plot(x, poisson.pmf(x, mu), 'bo', ms=7,alpha=0.5,label='poisson pmf') ## pmf - > 확률 질량 함수
ax.vlines(x, 0, poisson.pmf(x, mu), colors='g', lw=5, alpha=0.5) ## 수직선
ax.legend()
plt.show()

## 시뮬레이션 수행
poisson_ = poisson(mu) # 분포 추정
data_poisson=poisson_.rvs(1000) #추정한 분포 기반 Sample 1000개 생성
##(data_poisson=np.random.poisson(3, 1000) ## numpy 기반 Sample 1000개 생성)

## 성공횟수(x)별 횟수 도출 및 Dataframe 저장
count = Counter(data_poisson)
count_data = pd.DataFrame.from_dict(count, orient="index").reset_index()
count_data =count_data.rename(columns={"index":"number", 0:"count"})

## 성공횟수(x)별 정렬
count_data=count_data.sort_values(by=["number"], ascending=True)

## 시각화
plt.plot(count_data["number"], count_data["count"],'bo', ms=7,alpha=0.5)
plt.vlines(count_data["number"],0, count_data["count"],colors='r', lw=5, alpha=0.5)
plt.xlabel("number of success")
plt.ylabel("count")
plt.title("Simulation Result")
plt.xlim(-0.5, 7.5)
plt.show()

참조:

https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.poisson.html