확률공간과 확률변수

확률공간의 직관

"이 센서가 다음 측정에서 3.5보다 큰 값을 낼 확률"을 말하려면 세 가지를 정해야 한다. 일어날 수 있는 결과가 무엇인지, 확률을 매길 대상이 무엇인지, 그리고 그 대상에 어떤 숫자를 줄지다. 이 셋을 묶은 것이 확률공간(probability space) $(\Omega, \mathcal{F}, P)$다.

• 표본공간(sample space) $\Omega$: 실험의 가능한 결과 $\omega$ 전체의 집합. 주사위 하나면 $\{1, \dots, 6\}$, 부품 수명이면 $[0, \infty)$다.
• 사건(event): 확률을 매길 수 있는 $\Omega$의 부분집합. "짝수가 나온다" $= \{2, 4, 6\}$. 사건들의 모음을 $\mathcal{F}$라 쓴다.
• 확률 측도(probability measure) $P$: 사건마다 $[0, 1]$의 숫자를 주는 함수.

$P$는 세 가지 공리(Kolmogorov 공리)만 만족하면 된다.

첫째는 확률이 음수가 아니라는 것, 둘째는 "무언가는 반드시 일어난다"는 것, 셋째는 겹치지 않는 사건들의 확률은 더해진다는 것이다. 여사건의 확률 $P(A^c) = 1 - P(A)$, 합사건의 확률 $P(A \cup B) = P(A) + P(B) - P(A \cap B)$ 같은 익숙한 규칙은 모두 이 공리에서 나온다.

확률변수

확률공간은 결과 $\omega$를 그대로 다루지만, 연구에서 관심 있는 것은 결과에서 계산한 숫자다. 주사위 두 개를 던진 결과 $\omega = (3, 5)$에서 관심 있는 것은 합 8이고, 대기 관측 하루에서 관심 있는 것은 강수량 mm다. 결과를 숫자로 바꾸는 규칙이 확률변수다.

정의. 확률변수(random variable) $X$는 표본공간에서 실수로 가는 함수 $X: \Omega \to \mathbb{R}$이다. 무작위인 것은 함수가 아니라 어떤 $\omega$가 실현되는가다. $\omega$가 정해지면 $X(\omega)$는 하나로 정해진다. 그래서 같은 실험 위에 주사위 합, 최댓값, 첫 주사위 값처럼 여러 확률변수를 함께 정의할 수 있고, 이들 사이의 관계(11.3장의 결합분포)를 말할 수 있다.

$X$에 대한 질문은 $\Omega$의 사건으로 바뀐다. $\{X \le x\}$는 $\{\omega : X(\omega) \le x\}$의 줄임말이고, 그 확률을 $P(X \le x)$라 쓴다. 대문자 $X$는 확률변수, 소문자 $x$는 그것이 취한 값이다.

확률변수는 값의 모양에 따라 셋으로 나뉜다.

• 이산(discrete): 셀 수 있는 값만 갖는다. 검출기에 들어온 광자 수, 불량품 수.
• 연속(continuous): 구간의 값을 갖고, 한 점의 확률은 모두 0이다. 부품 수명, 기온.
• 혼합(mixed): 둘이 섞인 것. 하루 강수량은 비가 오지 않은 날 정확히 0이고(이산 부분), 비가 온 날 연속적인 양수다(연속 부분).

혼합형은 예외처럼 보이지만 과학 데이터에서 흔하다. 비가 오는 날이 전체의 40%이고, 비가 온 날의 강수량은 평균 5 mm인 지수분포를 따른다고 두고 시뮬레이션한다.

import numpy as np

rng = np.random.default_rng(0)
n_days = 100_000
wet = rng.random(n_days) < 0.4                                          # (N,) 비가 온 날
rain_mm = np.where(wet, rng.exponential(scale=5.0, size=n_days), 0.0)   # (N,) 강수량

print(f"P(R = 0)    ≈ {np.mean(rain_mm == 0):.3f}   (이론 0.600)")
print(f"P(R = 3 mm) ≈ {np.mean(rain_mm == 3.0):.3f}   (연속 부분의 한 점)")
for r in (0.0, 1.0, 5.0, 20.0):
    F_theory = 0.6 + 0.4 * (1 - np.exp(-r / 5.0))
    print(f"F({r:4.1f} mm) ≈ {np.mean(rain_mm <= r):.3f}   이론 {F_theory:.3f}")

P(R = 0)    ≈ 0.601   (이론 0.600)
P(R = 3 mm) ≈ 0.000   (연속 부분의 한 점)
F( 0.0 mm) ≈ 0.601   이론 0.600
F( 1.0 mm) ≈ 0.673   이론 0.673
F( 5.0 mm) ≈ 0.854   이론 0.853
F(20.0 mm) ≈ 0.993   이론 0.993

값 0에는 확률 0.6이 실려 있지만, 3 mm 같은 연속 부분의 한 점에는 확률이 없다. 이런 변수는 다음 절의 PMF나 PDF 하나로는 표현되지 않는다. 그러나 누적확률 $P(R \le r)$은 이론값($0.6 + 0.4(1 - e^{-r/5})$)과 잘 맞는다. 다음 절에서 볼 CDF가 모든 확률변수에 공통인 표현인 이유다.

분포의 표현

확률변수 $X$의 분포(distribution)는 "$X$가 어느 구간에 있을 확률"을 모두 모은 것이다. 이를 적는 방법이 네 가지 있고, 각각 쓰임이 다르다.

PMF. 이산확률변수의 확률질량함수(probability mass function)는 각 값의 확률이다.

PDF. 연속확률변수의 확률밀도함수(probability density function) $p(x)$는 구간의 확률을 적분으로 준다.

PDF 값 $p(x)$는 확률이 아니라 단위 길이당 확률이다. 짧은 구간 $[x, x + h]$의 확률이 약 $p(x)\,h$다. 그러므로 PDF 값은 1보다 클 수 있다. 폭 0.5인 구간에 확률 1을 고르게 펴면 밀도는 2가 된다.

from scipy import stats

X = stats.uniform(loc=0.0, scale=0.5)          # 0~0.5 균등분포
print("pdf(0.2) =", X.pdf(0.2))
print("P(0.1 < X <= 0.3) = F(0.3) - F(0.1) =", round(X.cdf(0.3) - X.cdf(0.1), 3))

x = X.rvs(size=200_000, random_state=np.random.default_rng(1))   # (N,)
print("표본에서 센 비율                    =", round(np.mean((x > 0.1) & (x <= 0.3)), 3))
for h in (0.1, 0.01, 0.001):
    frac = np.mean(np.abs(x - 0.2) < h / 2)                       # 폭 h 구간에 든 비율
    print(f"h = {h:<5}  P(|X - 0.2| < h/2) ≈ {frac:.4f}   비율 / h = {frac / h:.2f}")

pdf(0.2) = 2.0
P(0.1 < X <= 0.3) = F(0.3) - F(0.1) = 0.4
표본에서 센 비율                    = 0.402
h = 0.1    P(|X - 0.2| < h/2) ≈ 0.2005   비율 / h = 2.01
h = 0.01   P(|X - 0.2| < h/2) ≈ 0.0200   비율 / h = 2.00
h = 0.001  P(|X - 0.2| < h/2) ≈ 0.0021   비율 / h = 2.05

구간이 좁아질수록 그 구간에 든 확률은 0으로 가지만, 확률을 구간 폭으로 나눈 값은 PDF 값 2 근처에 머문다. 이것이 밀도의 뜻이다. $h = 0.001$에서 2.05로 조금 벗어난 것은 구간에 든 표본이 400개 정도로 적어 생긴 표본 잡음이다. 코드에서 scipy.stats의 분포 객체는 pdf(또는 이산 분포의 pmf), cdf, ppf(quantile function), rvs(표본 추출)를 같은 이름으로 제공한다. 이 책은 이 인터페이스를 계속 쓴다.

CDF. 누적분포함수(cumulative distribution function)는

다. 이산·연속·혼합 어느 경우에도 정의된다는 것이 CDF의 장점이다. 성질은 세 가지다. 감소하지 않고, $x \to -\infty$에서 0, $x \to \infty$에서 1로 가며, 오른쪽 연속이다. 구간 확률은 $P(a < X \le b) = F(b) - F(a)$이고, 한 점의 확률 $P(X = x)$는 그 점에서 CDF가 뛰는 높이다. 연속이면 $p(x) = F'(x)$, 이산이면 CDF가 PMF의 값만큼 계단처럼 뛴다.

Quantile function. CDF의 역함수 $F^{-1}(p)$는 "누적확률이 $p$가 되는 값"이다. $F^{-1}(0.5)$가 중앙값, $F^{-1}(0.9)$가 90% 분위수다. CDF가 계단이거나 평평한 구간이 있으면 보통의 역함수가 없으므로 일반화 역함수를 쓴다.

"누적확률이 처음으로 $p$에 닿는 $x$"다. 이 정의로 $F(F^{-1}(p)) \ge p$가 항상 성립하고, 연속이며 순증가하는 CDF에서는 등호가 된다.

T = stats.expon(scale=2000.0)                  # 부품 수명(시간), 평균 2000 h
for p in (0.1, 0.5, 0.9):
    t_h = T.ppf(p)                             # quantile function F^{-1}(p)
    print(f"F^-1({p}) = {t_h:7.1f} h   F(F^-1(p)) = {T.cdf(t_h):.3f}")

K = stats.poisson(mu=3.0)                      # 이산 분포: CDF가 계단
print("Poisson(3)  F(2) =", round(K.cdf(2), 3), "  F(3) =", round(K.cdf(3), 3))
print("F^-1(0.5) =", K.ppf(0.5), "  (F(k) >= 0.5인 가장 작은 k)")

u = np.random.default_rng(2).random(100_000)   # (N,) Unif(0, 1)
t_sim = T.ppf(u)                               # 역변환 샘플링
print(f"F^-1(U)의 표본 중앙값 {np.median(t_sim):.1f} h   이론 중앙값 {T.median():.1f} h")

F^-1(0.1) =   210.7 h   F(F^-1(p)) = 0.100
F^-1(0.5) =  1386.3 h   F(F^-1(p)) = 0.500
F^-1(0.9) =  4605.2 h   F(F^-1(p)) = 0.900
Poisson(3)  F(2) = 0.423   F(3) = 0.647
F^-1(0.5) = 3.0   (F(k) >= 0.5인 가장 작은 k)
F^-1(U)의 표본 중앙값 1390.7 h   이론 중앙값 1386.3 h

세 가지를 읽는다. 첫째, 평균 수명이 2000 h인 지수분포에서 중앙값은 1386 h로 평균보다 짧다. 오른쪽 꼬리가 길어 평균을 끌어올리기 때문이다. "부품 절반이 1386 h 안에 고장 난다"와 "평균 수명 2000 h"는 같은 분포에 대한 다른 진술이다. 둘째, Poisson(3)은 $F(2) = 0.423 < 0.5 \le F(3) = 0.647$이므로 중앙값이 3이다. 계단 CDF에서는 일반화 역함수가 정확히 이렇게 동작한다. 셋째, 균등난수 $U$를 quantile function에 넣으면 그 분포의 표본이 나온다. $P(F^{-1}(U) \le x) = P(U \le F(x)) = F(x)$이기 때문이다. 이 역변환 샘플링(inverse transform sampling)은 11.6장에서 다른 샘플링 방법과 함께 다룬다.

아래 위젯에서 분포와 모수를 바꾸며 위쪽의 PDF(이산 분포는 PMF)와 아래쪽의 CDF를 함께 본다. $x$를 움직이면 $x$ 왼쪽의 넓이(이산 분포는 막대의 합)가 아래 그림의 $F(x)$ 높이와 같다는 것을 확인할 수 있다.

네 표현의 쓰임. PMF·PDF는 우도(12.2장)와 모델의 손실을 쓸 때, CDF는 구간 확률·초과확률·혼합형 변수를 다룰 때, quantile function은 예측 구간·위험 한계(예: "99% 확률로 이 하중을 넘지 않는다")·샘플링에 쓴다. 같은 분포를 어느 표현으로 보느냐에 따라 눈에 띄는 성질이 다르다. 밀도의 봉우리는 PDF에서, 꼬리의 무게는 CDF의 $1 - F(x)$나 quantile에서 잘 보인다.

변수변환

확률변수를 함수에 통과시키는 일은 흔하다. 측정한 전압을 온도로 환산하고, 양수 파라미터를 $\log$ 공간에서 다루고, 신경망 출력에 softplus를 씌운다. $X$의 분포를 알 때 $Y = g(X)$의 분포는 무엇일까.

CDF로 구하기. 가장 일반적인 방법은 CDF의 정의로 돌아가는 것이다. $g$가 증가함수이면

이다. 감소함수이면 부등호가 뒤집혀 $F_Y(y) = 1 - F_X(g^{-1}(y))$가 된다(연속일 때).

밀도 공식. 연속이고 $g$가 순단조이며 미분 가능하면, 위 식을 $y$로 미분해

를 얻는다. 첫 인자는 "$y$로 가는 $x$의 밀도", 둘째 인자는 $g$가 구간 길이를 얼마나 늘리거나 줄이는지를 보정하는 Jacobian이다. $g$가 구간을 두 배로 늘리면 같은 확률이 두 배 길이에 펴지므로 밀도는 절반이 된다. 다변수에서는 이 인자가 Jacobian 행렬식 $|\det \mathbf{J}_{g^{-1}}|$이 되고, normalizing flow(19.4장)가 바로 이 공식으로 밀도를 계산한다. $g$가 단조가 아니면(예: $Y = X^2$) $y$로 가는 $x$가 여러 개이므로 각 가지의 기여를 더한다.

$X \sim \mathcal{N}(0, 0.5^2)$, $Y = e^X$로 확인한다. $g^{-1}(y) = \log y$, $|dg^{-1}/dy| = 1/y$이므로 $p_Y(y) = p_X(\log y)/y$이고, 이것이 로그정규(log-normal) 분포다.

X = stats.norm(0.0, 0.5)
Y = stats.lognorm(s=0.5)                       # exp(X)의 분포 (scipy의 로그정규)
y_grid = np.array([0.5, 1.0, 1.5, 2.0])
p_formula = X.pdf(np.log(y_grid)) / y_grid     # p_X(g^{-1}(y)) |d g^{-1}/dy|
print("공식   p_Y(y):", np.round(p_formula, 4))
print("scipy  p_Y(y):", np.round(Y.pdf(y_grid), 4))

y = np.exp(X.rvs(size=200_000, random_state=np.random.default_rng(3)))   # (N,)
print(f"중앙값  exp(X의 중앙값) = {np.exp(X.median()):.3f}   Y 표본 {np.median(y):.3f}")
print(f"평균    exp(X의 평균)   = {np.exp(X.mean()):.3f}   Y 표본 {y.mean():.3f}")
print(f"최빈값  exp(X의 최빈값) = {np.exp(0.0):.3f}   Y 이론 {np.exp(-0.5**2):.3f}")

공식   p_Y(y): [0.6105 0.7979 0.3829 0.1526]
scipy  p_Y(y): [0.6105 0.7979 0.3829 0.1526]
중앙값  exp(X의 중앙값) = 1.000   Y 표본 1.000
평균    exp(X의 평균)   = 1.000   Y 표본 1.133
최빈값  exp(X의 최빈값) = 1.000   Y 이론 0.779

공식이 scipy의 로그정규 밀도와 일치한다. 아래 세 줄이 이 절의 핵심이다. 단조 증가 변환은 분위수를 보존한다. $P(Y \le e^m) = P(X \le m)$이므로 중앙값은 그대로 $e^0 = 1$로 간다. 그러나 평균은 보존되지 않는다. $e^x$가 볼록하므로 $\mathbb{E}[e^X] > e^{\mathbb{E}[X]}$이고(Jensen 부등식, 11.5장), 이론값은 $e^{0.5^2/2} \approx 1.133$이다. 최빈값도 보존되지 않는다. 밀도에 Jacobian $1/y$가 곱해져 봉우리가 $e^{-0.25} \approx 0.779$로 옮겨 간다. 밀도의 봉우리는 좌표를 바꾸면 움직이는 양이다.

연구에서의 의미. 로그 공간에서 학습한 모델의 예측을 $\exp$로 되돌리면 그것은 원래 공간의 평균이 아니라 (로그 공간 예측이 중앙값이었다면) 중앙값이다. 반응 속도, 농도, 입자 크기처럼 양수이고 오른쪽으로 긴 양을 로그 변환해 회귀하는 일이 흔한데, 이때 되돌린 예측값이 평균보다 체계적으로 작게 나온다. 또 MAP 추정(13.1장)은 사후밀도의 봉우리이므로 파라미터를 어떤 좌표로 적느냐에 따라 바뀐다. 분위수와 CDF는 좌표에 따라 바뀌지 않는다.