[data analysis] 베르누이와 이항확률분포

베르누이분포(Bernoulli distribution)와 이항확률 분포

내용

• 베르누이분포(Bernoulli distribution)
• 이항확률분포(Binomial distribution)

베르누이분포(Bernoulli distribution)

한번의 시행에서 성공 또는 실패(1또는 0)의 결과만을 보이는 확률분포를 베르누이분포(Bernoulli distribution)라고 합니다. 확률변수는 두개의 값만을 포함합니다. 확률변수는 두개의 값만을 포함하며 확률질량함수는 식 1과 같이 나타낼 수 있습니다.

$$\begin{array}{rl}f(x)& =P(X=x)\\ & =\{\begin{array}{ll}p& \text{for}x=1\\ 1-p& \text{for}x=0\end{array}\end{array}$$

(식 1)

위 확률질량함수(PMF)는 하나의 매개변수(parameter) 즉, 확률 p에 의해 결정됩니다. 그러므로 이 분포는 식 2와 같이 나타냅니다.

X ~ Bernoulli(P)

(식 2)

import numpy as np
import pandas as pd  
from scipy import stats
from sympy import *
import matplotlib.pyplot as plt

예 1)

1개의 주사위를 시행하는 경우 확률변수는 다음과 같습니다.

표 주사위 시행에서의 확률변수

주사위 눈	x
1 or 3	1
other	0

이 분포의 확률질량함수(PMF)와 E(x)?

확률함수와 기대값은 식 3과 같이 계산할 수 있습니다.

$$\begin{array}{rl}f(x)& ={\left(\frac{1}{3}\right)}^x{\left(\frac{2}{3}\right)}^{1-x}\\ E(x)& =1\cdot \frac{1}{3}+0\cdot \frac{2}{3}\end{array}$$

(식 3)

베르누이 분포의 PMF, CDF, 기대값(평균), 분산등의 통계량은 scipy.stats 모듈의 bernoulli() 클래스를 사용하여 계산할 수 있습니다.

pmf=stats.bernoulli.pmf(1, 1/3)
print("PMF: %.3f"%pmf)

PMF: 0.333

mu, var=stats.bernoulli.stats(1/3, moments="mv")
print("mean: %.3f, variation: %.3f"%(mu,var))

mean: 0.333, variation: 0.222

이항확률분포(Binomial distribution)

이항분포(Binomial distribution)는 위의 베르누이 분포를 여러 번 시행하는 경우의 확률분포를 의미합니다. 예를 들어 동전을 3번 던지는 경우 앞면이 나오는 수를 랜덤변수 X로 한다면 랜덤변수의 범위는 S_x = {0, 1, 2, 3}가 됩니다. 이 경우 앞면이 나오는 확률을 p라고 하면 각 변수 값에 대한 확률은 식 4와 같이 나타낼 수 있습니다.

f(0) = P(X = 0) = P(TTT) = (1 − p)3	(식 4)
f(1) = P(X = 1) = P(TTH or THT or HTT) = 3(1 − p)2p
f(1) = P(X = 1) = P(TTH or THT or HTT) = 3(1 − p)2p
f(3) = P(X = 3) = P(HHH) = p3

위의 과정들을 조합공식을 사용하면 이항분포의 확률질량함수는 식 5과 같이 정의 됩니다.

$$\begin{array}{rl}f(x)& =P(X=x)\\ & =(\genfrac{}{}{0ex}{}{n}{k})p^k(1-p{)}^{n-k}\\ & k:0,1,2,\cdots ,n\end{array}$$

(식 5)

n번 시행할 경우 k번의 성공확률은 위와 같은 식으로 계산됩니다. 다시말하면 이항분포는 시행횟수 n과 확률 p에 의해 특성됩니다. 즉, 이 분포는 모수와 함께 식 6과 같이 나타냅니다.

X ∼ B(n, p)

(식 6)

이항분포의 확률질량함수는 scipy.stats.binom() 클래스의 메소드 pmf를 사용하여 계산할 수 있습니다.

예 2)

B(10, 0.3)과 B(20, 0.6)의 이항분포를 시각화합니다.

pmf1=stats.binom.pmf(range(0, 11), 10, 0.3)
print(np.around(pmf1, 3))

[0.028 0.121 0.233 0.267 0.2   0.103 0.037 0.009 0.001 0.    0.   ]

pmf2=stats.binom.pmf(range(0, 21), 20, 0.6)
print(np.around(pmf2, 3))

[0.    0.    0.    0.    0.    0.001 0.005 0.015 0.035 0.071 0.117 0.16
 0.18  0.166 0.124 0.075 0.035 0.012 0.003 0.    0.   ]

그림 1. 두 이항분포.

fig, ax=plt.subplots(figsize=(4,3))
ax.bar(range(0, 11), pmf1, color="brown", alpha=0.3, label="B(10, 0.3)")
ax.bar(range(0, 21), pmf2, color="g", alpha=0.3, label="B(20, 0.6)")
ax.set_xlabel("x")
ax.set_ylabel("Probability")
ax.legend(loc="best")
plt.show()

그림 1에 의하면 이항분포의 모수 n, p는 분포의 형태의 변화를 보여줍니다. 이항분포는 베르누이 분포를 여러 번 시행하는 것으로 각 베르누이 시행의 합이 이항분포가 됩니다. 즉, 동전을 1회 던지는 시행은 Bernoulli(p)를 확률질량함수로하는 베르누이 확률분포입니다. 이러한 시행을 n번 반복 시행할 경우 각각의 베르누이 시행의 합은 이항분포 즉, B(n, p)를 구성합니다.

[이항 분포]

x₁, x₂, · · ·의 각 사건이 Bernoulli(p)을 따르는 랜덤변수이면 각 변수들의 합 (x = x₁ + x₂ + · · ·)의 분포는 이항분포 B(n, p)가 됩니다.

위 정의에 의하면 동일한 확률을 가진 두 개의 이항분포의 사건들은 그 각 사건들을 더하는 것으로 새로운 확률변수를 생성하는 것이 가능 합니다. 즉, 식 7과 같이 두 이항변수들(X와 Y)의 합은 다음과 같이 새로운 확률변수 Z으로 나타낼 수 있습니다.

X ∼ B(n, p), Y ∼ B(m, p)	(식 7)
Z ∼ B(n + m, p)

예 3)

어느 공장 생산품 중에 불량률이 0.01이라면 생산품 중 랜덤하게 30개를 선택할 경우 적어도 2개의 불량품이 포함될 확률?

이 사건의 표본 공간은 S_x = {0, 1, 2, · · · }이며 이항분포로 식 8와 같이 나타낼 수 있습니다.

X ∼ B(30, 0.01)

(식 8)

이확률분포에서 2개이상의 불량품이 포함될 확률은 식 9와 같이 계산됩니다.

P(X ≥ 2) = 1 − P(X < 2)

(식 9)

Pmore2=1-np.sum(stats.binom.pmf([0, 1], 30, 0.01))
print("불량품이 2개 이상이 될 확률은 %.3f 입니다. "%Pmore2)

불량품이 2개 이상이 될 확률은 0.036 입니다. 

round(1-stats.binom.cdf(1, 30, 0.01), 3)

0.036

예 4)

주사위 던지기 세번 시행에서 특정한 수에 따라 상금이 정해지는 규칙을 가지는 게임을 합니다.

표 주사위 게임 규칙

횟수	0	1	2	3
상금	-1000	1000	2000	3000

이 게임의 기대값을 결정합니다.

$$\begin{array}{rl}X& \sim B(3,\frac{1}{6})\\ f(x)& =(\genfrac{}{}{0ex}{}{3}{x}){\left(\frac{1}{6}\right)}^x{\left(\frac{5}{6}\right)}^{3-x}\end{array}$$

(식 10)

X=np.array([-1000, 1000, 2000, 3000])
pf=stats.binom.pmf(range(4), 3, 1/6)
print(np.around(pf, 3))

[0.579 0.347 0.069 0.005]

E=np.sum(X*pf)
round(E, 3)

-78.704

결과로는 약 79원의 손실이 예상됩니다.

예 5)

다음 조건하에 주식 매매를 합니다.

• 시가(open price) 대비 종가(close price)가 1% 이상의 증가되는 확률: 0.45
• 시가에 매수 후 종가가 1% 증가 시 매도
• 확률변수 x는 매도하는 횟수

다음을 계산합니다.

1. 10번 매수 후 매매가 이루어지는 평균횟수?
2. 10번 매수에서 5번 이상 매매가 이루어질 확률?

10회 횟수동안에 대한 표본공간과 확률분포는 식 11과 같이 나타낼 수 있습니다.

Sx = {0, 1, 2, …, 10},   X ~ B(10, 0.045)

(식 11)

E=stats.binom.stats(10, 0.45, moments="m");E

4.5

위 결과는 평균적으로 10번 중 4번의 거래가 이루어질 것으로 기대됩니다.

10번의 매입에서 5번이상 매매가 이루어질 확률은 식 12와 같이 누적확률로 계산할 수 있습니다.

P(X ≥ 5) = 1 − P(X ≤ 4)

(식 12)

round(1-stats.binom.cdf(4, 10, 0.45), 3)

0.496

이항분포의 기대값과 분산은 모멘트생성함수(M_X(t))의 1차 미분과 2차 미분으로 결정할 수 있습니다. 식 13은 이항분포의 PMF를 적용한 모멘트 생성함수입니다.

$$\begin{array}{rl}{M}_x(t)& =E(e^{tx})\\ & =\sum _{x=0}^n{e}^{tx}(\genfrac{}{}{0ex}{}{n}{x})p^x(1-p{)}^{n-x}\\ & =\sum _{x=0}^n(\genfrac{}{}{0ex}{}{n}{x})(p{e}^t{)}^x(1-p{)}^{n-x}\\ & ={(p{e}^t+1-p)}^n\\ \therefore & (p+q{)}^n=\sum _{x=0}^n(\genfrac{}{}{0ex}{}{n}{x})p^x{q}^{n-x}\end{array}$$

(식 13)

식 13의 미분을 계산하기 위해 파이썬 패키지인 sympy를 사용합니다.

t, p, n = symbols("t p n")
M=(p*exp(t)+1-p)**n 
dM=M.diff(t)
print(dM)

n*p*(p*exp(t) - p + 1)**n*exp(t)/(p*exp(t) - p + 1)

위 결과를 일반적인 방법으로 표현하면 식 14와 같습니다.

$$\frac{np(p{e}^t-p+1{)}^n{e}^t}{p{e}^t-p+1}$$

(식 14)

E=dM.subs(t, 0); E

np

분산은 이항분포 MGF(식 13)의 2차 미분을 고려하여 식 15와 같이 계산됩니다.

$${\sigma }^2=\frac{d^2(M_x(t))}{d{t}^2}(0)-(E(X){)}^2$$

(식 15)

ddM=M.diff(t, 2)
ddM_0=ddM.subs(t, 0); ddM_0

np(np-p+1)

var=ddM_0-E**2
simplify(var)

np(1-p)

위 결과를 정리하면 이항분포의 기대값과 분산은 식 16과 같습니다.

E(X) = np	(식 16)
Var(X) = np(1 − p)

예 6)

어느 학생이 5지선다형 문제 15개에서 랜덤으로 답을 선택하였을 경우 5개 이상 답을 맞힐 확률, 기대값과 분산?

1개의 문제에서 답을 맞힐 확률 $p=\frac{1}{5}$ 를 가지는 이항확률입니다(식 17).

$$\begin{array}{rl}{S}_x& =\{0,1,2,\cdots \}\\ f(x)& =(\genfrac{}{}{0ex}{}{15}{x}){\left(\frac{1}{5}\right)}^x{\left(\frac{4}{5}\right)}^{15-x}\end{array}$$

(식 17)

pMore5=1-stats.binom.cdf(4, 15, 1/5); round(pMore5, 2)

0.16

mu, var=stats.binom.stats(15, 1/5, moments="mv")
print("평균: %.2f, 분산: %.2f"%(mu, var))

평균: 3.00, 분산: 2.40