원문 : http://jakevdp.github.io/blog/2014/06/06/frequentism-and-bayesianism-2-when-results-differ/
이전 포스팅에서, 빈도주의과 베이지안주의가 과학적 데이터 분석에 연관되어 있기 때문에 이에 대한 간단한 실전예제를 보였다. 여기에서, 빈도주의와 베이즈론주의간의 근본적 철학의 차이를 다뤘고, 두 접근은 몇몇 간단한 문제에 대해서는 같은 결과를 낸다는 것을 보였다. 이렇게 간단한 문제에 대해서 종종 같은 결과를 낸다는 것을 보여주는것은 쉬웠으나, 복잡한 문제들에 대해서는 상당히 다른 결과를 낼 수 있는것 또한 사실이다. 이 차이는 다음의 두가지 상황에서 가장 명백히 드러난다.
1. 성가신(Nuisance) 가진 파라미터들을 다룰때.
2. 빈도주의론자들의 신뢰구간(confidence intervals과 베이지안들의 신뢰구간(credible intervals) 사이의 미묘한(때로는 간과되는) 차이.
두번째 점은 약간 더 철학적이고 깊이 있는 이야기다. 이 두번째 점에 대해서는 이후의 포스팅에서 이야기 하고, 여기서는 첫번째, Nuisance 파라미터들을 처리하는 것에서 빈도주의와 베이지안 사이의 차이에 대해서 이야기 하겠다. 이전 포스팅까지는 공정한 입장(빈도주의도, 베이지안도 아닌)에서 포스팅을 하였으나, 여기서는 베이지안의 입장에서 다루고 있다. 이 포스팅은 또한 다소 논란이 될 2번째 상황에 대해 간단히 다루기도 한다.
What is a Nuisance Parameter?
nuisance 파라미터는 분석의 목적에 상관이 없으나, 일정정도의 관심을 결정하는데 필요한 파라미터를 말한다. 예를 들어, 이전 포스팅에서 다뤘던, 두번째 (Varying Photon Counts : The Bayesian Approach)적용에서, 우리는 관찰된 photons에서 평균
이것을 설명하기위해, nuisance 파라미터가 쓰이는 2가지 예제를 들겠다. 이 예제에서 빈도주의론자와 베이지안들 둘다의 입장에서 이 문제를 다룰 것이다.
Example #1 : The Bayesian Billiard Game
nuisance 파라미터의 한 예제로 시작하겠다. 여차저차해서 토마스 베이즈가 쓴 1763년의 사후 논문으로 돌아간다. 이 문제의 2004년의 Eddy로부터 가져왔다.
배경은 Alice와 Bob이 바로 알수 없는 한 프로세스의 결과에 대해 내기를 하는 게임이다.
Alice와 Bob은 한 방에 들어갔다. 커튼뒤에는 당구대가 있었는데, 이 둘은 볼 수 없었지만, 그들의 친구인 Carol은 볼 수 있었다. 캐롤은 당구대에 공 하나를 굴리고, 이 공이 착지한 자리에 표시를 했다. 이 표시를 유지한채, 캐롤은 다음 공을 테이블에 굴렸다. 만약 이 공이 마크한 곳의 왼쪽에 착지한다면, 앨리스가 한점을 얻고, 만약 오른쪽에 착지를 한다면, 밥이 한점을 얻는다. 우리는 예제를 위해 캐롤이 던지는게 불편향되었다고 가정한다. 즉 공들은 테이블상에서 어디에 착지하든 동일한 확률을 갖는다. 6점에 먼저 도달한 사람이 이기는 게임이다.
첫번째 착지점에 의해 표시된 곳은 바로 nuisance 파라미터가 된다. 즉 이것은 현재 모르는 것이고, 우리가 실제 구하려는 것은 아니나, 목표(다음 착지점, 점수를 내기위한 착지점)를 위해서는 반드시 필요한 것이다. 만약 처음 던진 공이 아주 오른쪽에 착지했다면, 마크를 중심으로 왼쪽 공간이 더 많으므로, 다음 던진 공은 앨리스에게 우리 할 것이다. 대신 만약 처음 던진 공이 아주 왼쪽에 착지 했다면, 밥에게 유리해질 것이다.
이러한 배경을 가정하고, 의문 하나가 든다.
하나의 특정한 게임에서, 8번 굴린후에, 앨리스는 5점 그리고 밥은 3점을 얻었다. 밥이 이 게임을 이길 확률은 어떻게 되는가?
직관적으로, 여러분은 아마도 8점중에서 5점을 얻었기 때문에, 마크의 위치가 그녀에게 더 유리하게 놓여졌을것이라고 깨달을 것이다. 이것을 가정하면, 다음 시행 역시, 앨리스에게 유리할 것이다. 그리고 그녀는 밥이 이기기 전에, 그녀에게 유리한 던지기를 할 기회가 3번이나 있다. 따라서 그녀가 거의 이긴것으로 보인다. 그러나 정량(계량)적으로, 밥이 짜릿한 승리를 얻을 확률이 어떻게 되는가?
A Naïve Frequentist Approach
전통적인 빈도주의론자의 접근을 따르는 누군가는 다음과 같은 이유를 들것이다.
결과를 내기 위해, 우리는 마커가 표시된 곳의 중간 추정이 필요하다. 우리는 이 표시의 위치를 어느 시행이 주어지든 간에 엘리스진영에 착지할 확률 P로 정량화할 것이다. 우리는 바로 P의 최대우도추정이 다음과 같다는것을 볼 수 있다.
즉 , 그는 연속으로 3번을 이겨야 한다. 그러므로 우리는 다음 확률추정을 찾아야 한다.
p_hat = 5. / 8.
freq_prob = (1 - p_hat) ** 3
print("Naïve Frequentist Probability of Bob Winning: {0:.2f}".format(freq_prob))
다르게 말해, 우리는 밥에게 다음의 이길 확률, 혹은 승산비(odds, not probability)를 준다.
print("Odds against Bob winning: {0:.0f} to 1".format((1. - freq_prob) / freq_prob))
따라서, 우리는 빈도주의자들의 생각을 사용해 밥이 한번 이길때마다, 엘리스는 18번 이길것임을 추정했다. 다음은 베이지안적으로 접근해보자.
Bayesian Approach
우리는 또한 이 문제에 대해 베이지안의 입장에서 접근 할 수 있다. 이 접근은 우리가 처음에 몇몇 정의를 더 해야 한다.
우리는 다음의 확률변수(random variables)를 고려할 것이다.
- B = Bob wins
- D = observed data, i.e. D=
= (5, 3)
- p = 현재 공이 엘리스 쪽으로 착지할 알수없는 확률
우리는 
nuisance 파라미터를 다루는데 일반적인 베이지안 방법은 marginalization이다. 이것은 전체 nuisance 파라미터의 범위에 대한 결합확률(joint probability)을 합치는 것이다. 이 경우, marginalization은 우리가 일단 결합분포 
이 identity는 확률의 법칙과 조건부확률의 정의를 따른다. 다시말해 이것은 확률법칙의 기본적인 결과이고, 항상 참일것이다. 심지어는 빈도주의론자도 이것을 인정할 것이다. 그들은 단지 
Building our Bayesian Expression
이 결과를 계산하기 위해, 계산할 수 있는 다른 양에 대해서 이것을 표현할 수 있을 때까지
다음은 베이즈법칙을 사용해 
마지막으로 우리가 시작했던 같은 확률성질을 사용해, 분모에 있는
이제 우리가 원하는 확률이 우리가 계산할 수 있는 3가지 계량들로 표현될 수 있다. 이제 차례대로 각각의 계량들을 보자.
: 이 term은 우리가 위에서 사용한 빈도주의자들의 우도(likelihood)이다. 예제를 통해 설명하자면, 마커의 위치 p와 앨리스가 5번을 이기고 밥이 3번을 이긴다는 것이 주어지면, 밥이 연속으로 3번 이겨 6번의 승리를 할 확률이 어떻게 되는가? 이 확률을 나타내면
이 된다.
: 이것 또한 쉽게 구할수 있다. 예제를 통해 설명하면, 확률p가 주어지면 8번의 시행중 5번의 긍정결과(성공)에 대한 우도가 어떠한가? 답은 잘알려진 이항분포로 부터 얻을 수 있다. 식으로 나타내면
:이 확률 p는 우리의 사전확률이다. 문제 정의에 의해, p가 0~1사이의 확률로 동등하게 뽑힌 것을 가정할 수 있다. 즉
이고 적분의 범위는 0~1이다.
두 적분은 0~1에서 계산된다.
이 적분은 다소 어려워 보인다. 사실 이 식은 베터함수의 특별한 경우이다. 그러나 간단하게 우리는 Scipy의 베터함수를 사용해 바로 이 적분을 계산할 것이다.
from scipy.special import beta
bayes_prob = beta(6 + 1, 5 + 1) / beta(3 + 1, 5 + 1)
print("P(B|D) = {0:.2f}".format(bayes_prob))
따라서 odds는 다음과 같다.
In [4]:
print("Bayesian odds against Bob winning: {0:.0f} to 1".format((1. - bayes_prob) / bayes_prob))
따라서 베이지안결과는 우리에게 10 대 1의 승산비를 주었다. 이는 빈도주의 접근으로 찾은 18 대 1의 승산비와는 상당히 다르다. 그렇다면 어느것이 맞는 것일까?
A Brute Force/Monte Carlo Approach
잘 정의되고 간단한 환경의 이러한 예제는 올바른 답을 얻기 위해, 사실 상대적으로 몬테카를로 시뮬레이션을 사용하기 쉽다. 이것은 기본적으로 가능한 결과 들의 강제적인 도표작성이다.( 쉽게 말해 인위적으로 모든 상황들을 실행시켜, 각 상황의 후보가 문제를 해결하는지를 확인하는것, Brute-force search) 이 예제를 들어 설명하면, 우리는 많은 수의 임의적 게임들을 생성한 후 밥이 우승할 (3번 연달이 쭉 이기는) 관련있는 게임들만 센다. 현재 문제는 특히나 더 간단하다. 왜냐하면, 관련된 많은 확률변수들이 균등분포(uniform distribution)하기 때문이다. 다음과 같이 numpy 패키지를 사용할 수 있다.
import numpy as np
np.random.seed(0)
# play 100000 games with randomly-drawn p, between 0 and 1
p = np.random.random(100000)
# each game needs at most 11 rolls for one player to reach 6 wins
rolls = np.random.random((11, len(p)))
# count the cumulative wins for Alice and Bob at each roll
Alice_count = np.cumsum(rolls < p, 0)
Bob_count = np.cumsum(rolls >= p, 0)
# sanity check: total number of wins should equal number of rolls
total_wins = Alice_count + Bob_count
assert np.all(total_wins.T == np.arange(1, 12))
print("(Sanity check passed)")
# determine number of games which meet our criterion of (A wins, B wins)=(5, 3)
# this means Bob's win count at eight rolls must equal 3
good_games = Bob_count[7] == 3
print("Number of suitable games: {0}".format(good_games.sum()))
# truncate our results to consider only these games
Alice_count = Alice_count[:, good_games]
Bob_count = Bob_count[:, good_games]
# determine which of these games Bob won.
# to win, he must reach six wins after 11 rolls.
bob_won = np.sum(Bob_count[10] == 6)
print("Number of these games Bob won: {0}".format(bob_won.sum()))
# compute the probability
mc_prob = bob_won.sum() * 1. / good_games.sum()
print("Monte Carlo Probability of Bob winning: {0:.2f}".format(mc_prob))
print("MC Odds against Bob winning: {0:.0f} to 1".format((1. - mc_prob) / mc_prob))
몬테카를로 접근은 밥에게 10대 1 의 승산비를 주었다. 이것은 베이지안 접근과 동일한 결과이다. 명백히 Naïve Frequentist 접근은 틀렸음을 알수 있다.
Discussion
이 예제는 nuisance 파라미터 p를 다루는데 몇몇 다른 접근법들을 보였다. 몬테카를로 시뮬레이션은 우리에게 베이지안 접근법과 동일한 참 확률에 대한 brute-force에 가까운 추정을 주었다. (우리가 한 가정이 유효하다는 가정하에) Naïve Frequentist 접근은 nuisance 파라미터 p 에 대한 최대우도법을 사용하여 잘못된 결과를 내었다.
나의 이 말이 큰 문제를 일으키기 전에, 나는 빈도주의론 자체가 틀렸다는 것을 말하는게 아니라는 것을 강조한다. 위의 잘못된 결과는 나이브(naive)한 접근이라는것에 더 초점이 맞춰져야 한다. nuisance 파라미터 같은 파라미터를 다루는 것에 대한 빈도주의론자들의 방법은 확실히 존재한다. 예를 들어 이것은 p 에 대한 종속을 고립시키기 위해, 이론적으로 데이터의 변환과 제약을 적용할 수 있다. 그러나 이 특정한 예제에 대해 베이지안의 p에 대한 marginalization을 사용하지 않는 어떤 다른 접근법도 찾을 수 없었다.
다른 잠재적인 논쟁은 전통적인 빈도주의 접근법에 불리한 방향으로 나온 질문 그 자체에 있다. 빈도주의는 영가설 검증(null test) 또는 신뢰구간에 대해 답을 말하기를 바랬을수도 있다. 다시말해 빈도주의론자들은 올바른 답이 몇몇 값은 p 확률(예를 들어 0.05, 위의 예제에서 든 p 와는 다름을 주의해라!)에 있고 유사한 시행들의 
이러한 두가지 빈도주의론자들의 반응들에는 한가지 명확한 공통점이 있다. 둘다 다소 노력과 특별한 전문지식이 필요하다는 것이다. 아마도 알맞은 빈도론자의 접근은 통계학 박사를 가진 사람에게는 명확해 보일것이다. 그러나 명확하게 바로 이 문제에 대해 답을 하려 노력하는 통계학 비전문가에게는 명확해 보이지 않다. 이러한 측면에서는, 베이즈론이 이러한 종류의 문제에 대해 더 나은 접근법을 제공한다고 생각한다. 베이즈론안에 있는 확률의 잘 알려진 몇몇 법칙들을 수학적으로 간단히 고치면, 우리는 다른 전문적인 지식없이도, 올바른 정답에 바로 도달할 수 있다.
다음에 nuisance 파라미터를 다룬 더 데이터에 기반한 예제들을 다룰 것이다.
Example #2 : Linear Fit with Outliers
nuisance 파라미터의 개념이 도움이 될 수 있는 한가지 상황은 데이터의 이상치(outlier)들을 설명한다는 것이다. 관찰된 변수들 
x = np.array([ 0, 3, 9, 14, 15, 19, 20, 21, 30, 35,
40, 41, 42, 43, 54, 56, 67, 69, 72, 88])
y = np.array([33, 68, 34, 34, 37, 71, 37, 44, 48, 49,
53, 49, 50, 48, 56, 60, 61, 63, 44, 71])
e = np.array([ 3.6, 3.9, 2.6, 3.4, 3.8, 3.8, 2.2, 2.1, 2.3, 3.8,
2.2, 2.8, 3.9, 3.1, 3.4, 2.6, 3.4, 3.7, 2.0, 3.5])이 데이터를 시각화했다.
%matplotlib inline
import matplotlib.pyplot as plt
plt.errorbar(x, y, e, fmt='.k', ecolor='gray');
우리의 목표는 데이터에 가장 잘 맞는 하나의 선을 찾는것이다. 보다싶이 몇개의 이상치들이 있음을 볼 수 있다. 그러나 non-robust 최대우도법으로 시작해 보자. 이전 포스트에서 본것과 같이, 다음 간단한 최대우도법의 결과는 빈도주의자 혹은 베이즈론자(균등 사전확률)에 따라 다른 결과가 나온다. 이런 간단한 문제에서는, 모든 접근법이 같은 결과를 낸다.
우리는 기울기와 절편이 벡터
이 모델을 고려하여, 우리는 각 점에 대해 가우시안 우도를 계산할 수 있다.
전체 우도는 전체 개별 우도의 곱이다. 이것을 계산하고 로그를 취하면 다음과 같이 나온다.
이것은 이전 포스팅을 읽었다면 상당히 익숙하게 보일 것이다. 마지막 식은 이 모형이
동일하게, 우리는 summation term을 최소화 할 수 있다. 이는 손실함수로 알려져 있다.
이 손실함수는 손실제곱(합)(squared loss, Sum of Squared error, SSE))으로도 알려졌다. 우리는 손실제곱이 가우시안 로그 우도로부터 파생되었음을 간단히 볼 수 있다.
Standard Likelihood Approach
이전 포스팅의 논리에 따르면, 우리는
간편함을 위해, 우리는 Scipy의 optimize 패키지를 사용해 손실함수를 최소화 할 것이다. ( 손실제곱의 경우, 행렬 방법을 사용하여, 더 계산을 효율적으로 할 수 있다. 그러나 우리는 간편성을 위해 그냥 numerical 최소화를 사용할 것이다.)
from scipy import optimize
def squared_loss(theta, x=x, y=y, e=e):
dy = y - theta[0] - theta[1] * x
return np.sum(0.5 * (dy / e) ** 2)
theta1 = optimize.fmin(squared_loss, [0, 0], disp=False)
xfit = np.linspace(0, 100)
plt.errorbar(x, y, e, fmt='.k', ecolor='gray')
plt.plot(xfit, theta1[0] + theta1[1] * xfit, '-k')
plt.title('Maximum Likelihood fit: Squared Loss');
이상치들이 회귀선에 불균형한 영향을 미치고 있음을 위 실험에서 알 수 있다. 이것은 손실제곱함수의 성질 때문이다. 만약 단 하나의 이상치, 예를들어 10 표준편차만큼 회귀선에서 떨어져있다고 한다면, 이것의 손실에 대한 영향은 2 표준편차만큼 떨어진 25개의 점들의 영향들 보다 클 것이다.
명백히, 손실제곱합(SSE)은 과도하게 이상치들에 민감하다. 이것은 우리의 회귀선에 문제를 일으킨다. 이 문제를 빈도주의입장에서 해결하는 방법은 손실함수를 더 강건(robust)하게 조정하는것이다.
Frequentist Corrections for Outliers : Huber Loss
다양한 가능한 손실함수들의 종류는 정말 무한대만큼 많다. 그러나 상대적으로 이 상황에 더 잘 맞는 손실함수는 Huber loss 이다. Huber loss는 2차차선(곡선)에서 1차선(직선)으로 변환을 하는 하나의 중요한 값(Critical value)을 정의한다. Huber loss 와 표준 손실제곱을 몇몇 c(critical value, 중요한 값들)에 대해 비교하는 그래프를 그려보자.
t = np.linspace(-20, 20)
def huber_loss(t, c=3):
return ((abs(t) < c) * 0.5 * t ** 2
+ (abs(t) >= c) * -c * (0.5 * c - abs(t)))
plt.plot(t, 0.5 * t ** 2, label="squared loss", lw=2)
for c in (10, 5, 3):
plt.plot(t, huber_loss(t, c), label="Huber loss, c={0}".format(c), lw=2)
plt.ylabel('loss')
plt.xlabel('standard deviations')
plt.legend(loc='best', frameon=False);
모델에 잘 맞는 점들에 대해서는 Huber loss는 squared loss가 같다. 그러나 이상치들의 영향은 줄었다. 예를 들어 회귀선에서 20 표준편차만큼 떨어진 한 점은200의 손실제곱합(SSE)을 갖는다. c=3인 Huber loss는 손실제곱합이 55이다. Huber loss를 사용하여 가장 적합된 회귀선을 찾아보자. 비교를 위해 밝은 회색으로 Squared loss의 결과를 그렸다.
def total_huber_loss(theta, x=x, y=y, e=e, c=3):
return huber_loss((y - theta[0] - theta[1] * x) / e, c).sum()
theta2 = optimize.fmin(total_huber_loss, [0, 0], disp=False)
plt.errorbar(x, y, e, fmt='.k', ecolor='gray')
plt.plot(xfit, theta1[0] + theta1[1] * xfit, color='lightgray')
plt.plot(xfit, theta2[0] + theta2[1] * xfit, color='black')
plt.title('Maximum Likelihood fit: Huber loss');
육안으로, 우리가 바란대로 잘 적합되었음을 알 수 있다.
그러나 한 베이지안은 이 새로운 손실함수의 도입의도가 약간 의심스럽다는 것을 지적할 것이다. 우리가 보았듯, squared-loss는 가우시안 우도로부터 직접 도출될 수 있다. Hubor loss 다소 약간 임의적이다. 다시말해 이것은 어느 (모집단)분포에서 오는가? 어떻게 c값을 정해야만 하는가? 이상치들에 대해 선형 손실함수를 사용하는 것에 대한 어느 좋은 이유라도 있는가? 또는 대신 이상치들을 그냥 단순이 제거하여야만 하는가? 어떻게 이러한 선택들이 우리의 결과모델에 영향을 미칠것인가?
A Bayesian Approach to Outliers : Nuisance Parameters
이상치를 설명하는것에 대한 베이즈주의 접근은 일반적으로 모델을 수정하여 이상치가 설명되게끔 조정하는 과정이 개입된다. 이 예제 데이터에서, 아주 분명하게 단순직선은 우리의 데이터에 좋은 적합이 아니라는 것이 명확하다. 따라서 이상치들을 설명하는 유연하고 더 복잡한 모델을 만들어보자. 한가지 선택은 signal(내 생각으로는 정규분포하는 회귀선이 설명하는 부분)과 background(내 생각으로는 정규분포하는 오차)를 혼합한 모형을 고르는 것이다.
우리가 여태 했던것에 약간의 nuisance 파라미터를 더 도입해 약간 더 확장된 것이다. 자세히 말해 
우리의 회귀모형은 현재 훨씬 더 복잡하다. 이것은 2개가 아닌 22개의 파라미터들이 있다. 그러나 대부분의 파라미터들은 위에서 말한 marginalization을 이용해 구할 수 있는 nuisance 파라미터이다. 위에서 든 당구게임의 예에서 모든 p의 경우를 통합했다. 이 우도를 사용하는 한 함수를 만들어보자. 이전 포스팅에서와 같이, 우리는 파라미터공간을 탐색하기 위해 emcee 패키지를 사용한다.
사실 이것을 계산하기 위해, 우리는 사전확률, 우도함수, 사후확률을 설명하는 함수들을 정의하는 것으로 시작할 것이다.
# theta will be an array of length 2 + N, where N is the number of points
# theta[0] is the intercept, theta[1] is the slope,
# and theta[2 + i] is the weight g_i
def log_prior(theta):
#g_i needs to be between 0 and 1
if (all(theta[2:] > 0) and all(theta[2:] < 1)):
return 0
else:
return -np.inf # recall log(0) = -inf
def log_likelihood(theta, x, y, e, sigma_B):
dy = y - theta[0] - theta[1] * x
g = np.clip(theta[2:], 0, 1) # g<0 or g>1 leads to NaNs in logarithm
logL1 = np.log(g) - 0.5 * np.log(2 * np.pi * e ** 2) - 0.5 * (dy / e) ** 2
logL2 = np.log(1 - g) - 0.5 * np.log(2 * np.pi * sigma_B ** 2) - 0.5 * (dy / sigma_B) ** 2
return np.sum(np.logaddexp(logL1, logL2))
def log_posterior(theta, x, y, e, sigma_B):
return log_prior(theta) + log_likelihood(theta, x, y, e, sigma_B)이제는, 우리는 파라미터 공간을 탐색하기 위해 MCMC 샘플링을 한다.
# Note that this step will take a few minutes to run!
ndim = 2 + len(x) # number of parameters in the model
nwalkers = 50 # number of MCMC walkers
nburn = 10000 # "burn-in" period to let chains stabilize
nsteps = 15000 # number of MCMC steps to take
# set theta near the maximum likelihood, with
np.random.seed(0)
starting_guesses = np.zeros((nwalkers, ndim))
starting_guesses[:, :2] = np.random.normal(theta1, 1, (nwalkers, 2))
starting_guesses[:, 2:] = np.random.normal(0.5, 0.1, (nwalkers, ndim - 2))
import emcee
sampler = emcee.EnsembleSampler(nwalkers, ndim, log_posterior, args=[x, y, e, 50])
sampler.run_mcmc(starting_guesses, nsteps)
sample = sampler.chain # shape = (nwalkers, nsteps, ndim)
sample = sampler.chain[:, nburn:, :].reshape(-1, ndim)
일단, 우리는 이 표본들을 가졌다. 우리는 마르코프 체인의 아주 좋은 특성을 사용할 수 있다. 그들의 분포는 사후확률을 만들기 때문에, 우리는 nuisance 파라미터를 간단히 무시하는 것(독립을 가정한다는건가?)으로 nuisance 파라미터들을 marginalize(결합)한다.
우리는 표본의 첫 2 column들을 갖고 실험을 하여 기울기와 절편의 marginalize된 분포를 볼 수 있다.
plt.plot(sample[:, 0], sample[:, 1], ',k', alpha=0.1)
plt.xlabel('intercept')
plt.ylabel('slope');
점들의 분포의 기울기가 약0.45정도 그리고 절편은 약 0.31 정도임을 볼 수 있다. 그러나 먼저 이 흔적으로부터 알수 있는 다른 정보들이 무엇인지를 보자.
MCMC 표본들을 분석하는 것 중 한가지 좋은 점은 nuisance 파라미터의 선택이 완전히 대칭적이라는 것이다. 우리가 
plt.plot(sample[:, 2], sample[:, 3], ',k', alpha=0.1)
plt.xlabel('$g_1$')
plt.ylabel('$g_2$')
print("g1 mean: {0:.2f}".format(sample[:, 2].mean()))
print("g2 mean: {0:.2f}".format(sample[:, 3].mean()))
이 둘중 어느것도 그렇게 강한 제약은 없었으나, 
이제 모든 정보들을 사용하여, 원래 데이터에 대한 최선의 marginalize된 회귀선을 그려보자. 덧붙여 모델이 이상치로 감지하는 점들을 빨간 원으로 그렸다.
theta3 = np.mean(sample[:, :2], 0)
g = np.mean(sample[:, 2:], 0)
outliers = (g < 0.5)
plt.errorbar(x, y, e, fmt='.k', ecolor='gray')
plt.plot(xfit, theta1[0] + theta1[1] * xfit, color='lightgray')
plt.plot(xfit, theta2[0] + theta2[1] * xfit, color='lightgray')
plt.plot(xfit, theta3[0] + theta3[1] * xfit, color='black')
plt.plot(x[outliers], y[outliers], 'ro', ms=20, mfc='none', mec='red')
plt.title('Maximum Likelihood fit: Bayesian Marginalization');
검은색 선이 우리의 직관과 딱 들어맞는다. 더욱이 이상치로 밝혀진 점들은, 실제 우리가 판단한 이상치이다. 비교를 위해, 회색선은 이전의 2가지 접근법(최대우도법과 huber loss에 기반한 빈도주의적 접근법)으로 그려진 회귀선을 그렸다.
Disscussion
여기서 우리는 선형회귀를 이상치의 존재로 나누었다. 전형적인 가우시안 최대우도법은 이상치를 설명하는데 실패했다. 그러나 빈도주의론적 접근법에서는 loss function을 수정(Huber loss)하여 해결하였고, 베이즈론의 접근법에서는 많은 수의 nuisance 파라미터를 복합모델(mixture model)에 써서 이러한 문제를 해결하였다.
두 접근법 모두 장단점이 있다. 여기서 빈도주의적 접근은 상대적으로 간단하고 계산효율적이나, 이론적인 근거가 잘 세워지지 않은(not particularly well-motivated) loss function을 사용하는것에 기반한다는 점이다. 베이즈론의 접근은 매우 근거가 잘 정립되어 있고 아주 좋은 결과들을 내지만, 다소 주관적인 사전확률선정을 필요하고 코딩하는 시간과 계산비용이 훨씬 더 많이 든다.
Conclusion... for now
이전의 포스팅에서, 빈도주의와 베이지안주의 사이의 철학적인 차이에 대해 논의를 했고 그들의 이러한 차이에도 불구하고, 그들은 간단한 문제에 대해서는 같은 결과를 준다는것을 보였다. 여기서는 이 두 접근의 결과가 달라지는 곳에서의 문제(nuisance 파라미터에 대한 설명)를 다루었다.
베이즈론 입장의 당구게임 예제에서, 우리는 naive 빈도주의 접근이 잘못된 결과로 이끈다는것을 보았다. 반면에 naive 베이즈 접근은 올바른 정답을 이끌었다. 이것은 빈도주의 자체가 틀렸다는 말이 아니고 우리가 이 빈도주의 접근을 적용할때 매우 주의해야한다는 것을 의미한다.
선형회귀예제에서, 우리는 우리의 데이터에 있는 이상치 설명에 대한 빈도주의와 베이즈론 각각으로부터의 한가지 가능한 접근을 봤다. 강건한 빈도주의 cost function (Huber loss)을 사용하여 상대적으로 빠르고 더 정확해졌다. 그러나 도입과정에서 근거가 부족하고 해석하기 어려운 결과들을 냈다. 베이지안의 혼합모델은 높은 계산비용과 더 수고스러운 과정들이 있으나, 많은 문제들이 한번에 해결되는 결과를 냈다. 이 경우, marginalizing은 최고의 회귀선을 찾는 한가지 방법과 동시에 데이터에서 이상치를 알아내는 방법이다.
그렇다면 어느 방법이 더 좋은가?, 빈도주의론 혹은 베이즈론?
답은 문제의 크기와 가능한 계산자원뿐만 아니라 아마 빈도주의론와 베이즈론의 방법을 분석가가 얼마나 잘 아느냐에 달렸다. 물리학의 배경을 가진사람은 베이즈론 방법이 좀 더 익숙한 경향이 있다. 왜냐하면 그들은 근본적인 원칙으로 부터 어떤것이든 도출이 가능하기 때문이다. 베이지안주의로 몇몇 확률법칙의 수학적인 수정에 기반해, 다양한 종류의 문제를 다루는 아주 유연한 방법들을 사용할 수 있다. 베이즈룰과 베이지안 확률 해석은 도박의 승산비(odds)를 계산하는 것에서부터 노이즈가 섞인 photometric data에서 태양계밖의 변화를 감지하는 것까지 가상의 어느 통계적 문제든지 풀 수 있게 적용할 수 있다. 베이지안의 패러다임에서, 복잡하고 중요한 분석을 하기위해, 전문용어(p-values, 영가설검정, 신뢰구간, breakdown points 등)와 빈도주의의 애매한 기술들을 이해하고 기억하는데 수년을 소비할 필요가 없다. 사람들이 Bayesian analysis가 어렵다고 하는 오해가 있다. 그러나 대조적으로 많은 과학자들은 이것이 쉽기 때문에 이 분석방법을 쓰고 있다.
쉽고 아름답다라는 말은 제쳐두고서라도, 이렇게 베이즈론을 강력하게 추천하는 중요한 이유가 이것말고 더 있다. 이것은 과학적 데이터의 문맥안에서 빈도주의의 신뢰구간(confidence intervals)와 베이지안의 신뢰구간(credibility regions)의 해석과 함께 설명해야 한다. 이에 대해서는 다음 포스팅에서 자세히 말하겠다.
'ML & DM > Bayesian ML' 카테고리의 다른 글
| Frequentism and Bayesianism: A Practical Introduction (0) | 2016.04.21 |
|---|---|
| Bayesian machine learning (2) | 2016.03.30 |
