Zguu's Data Science Blog

[베이지안 with Python] 'The Price is Right' Inference

 
  • Dec 09, 2019

프로그래머를 위한 베이지안 with 파이썬 서적의 Loss Function 챕터의 추론을 정리했다.

The Price is Right

이 TV쇼의 규칙

$\ $’The Price is Right’ 쇼케이스에서 가격을 최적화하는 베이지안 추론을 진행해보자. 규칙은 다음과 같다.

  1. 두 참가자가 쇼케이스에서 겨루게 된다.
  2. 각 참가자는 각각 다른 구성의 상품을 보게 된다.
  3. 상품을 관찰한 뒤 참가자는 자신들 앞에 놓인 상품의 가격을 써낸다.
  4. 써낸 가격이 실제 가격보다 높다면 그 참가자는 탈락한다.
  5. 써낸 가격이 실제 가격보다 $250 이내라면 참가자는 승리하고 상품을 경품으로 받는다.

사전확률분포

지금까지의 방송들을 종합해서, 실제 가격이 정규분포를 따른다고 가정한다.

$$True Price \thicksim Normal(\mu_p, \sigma_p)$$
$$\mu_p = 35000, \sigma_p = 7500$$

각가의 상품에 대한 가격들 또한 정규분포를 따른다.

$$Prize_i \thicksim Normal(\mu_i,\sigma_i),\ i = 1,2$$

각각의 상품이 \(\mu_i\)를 가질 것이라고 추측하지만, 이에 대한 불확실성은 \(\sigma_i\) 파라미터를 통해 표현된다.
실제 가격 True Price 는 다음과 같이 표현될 것이다. \(Price_1 + Price_2 + \epsilon\). 각 상품은 캐나다 토론토로 가는 여행권, 그리고 눈 제설기이다.
각 상품들의 가격들을 다음과 같이 정규분포로 표현하자.

$$Snow\ blower \thicksim Normal(3000,500)$$
$$Toronto \thicksim Normal(12000,3000)$$

만약 우리가 68.2%의 확률로 토론토로 가는 여행권이 1 표준편차만큼 떨어져있다고 믿는다면, 68.2%의 신뢰도로 이 여행권은 [9000, 15000] 의 범위에 존재할 것이다.

%matplotlib inline
import scipy.stats as stats
from IPython.core.pylabtools import figsize
import numpy as np
import matplotlib.pyplot as plt

figsize(12.5, 9)

norm_pdf = stats.norm.pdf

plt.subplot(311)
x = np.linspace(0, 60000, 200)
sp1 = plt.fill_between(x , 0, norm_pdf(x, 35000, 7500),
                color = "#348ABD", lw = 3, alpha = 0.6,
                label = "historical total prices")
p1 = plt.Rectangle((0, 0), 1, 1, fc=sp1.get_facecolor()[0])
plt.legend([p1], [sp1.get_label()])

plt.subplot(312)
x = np.linspace(0, 10000, 200)
sp2 = plt.fill_between(x , 0, norm_pdf(x, 3000, 500),
                 color = "#A60628", lw = 3, alpha = 0.6,
                 label="snowblower price guess")

p2 = plt.Rectangle((0, 0), 1, 1, fc=sp2.get_facecolor()[0])
plt.legend([p2], [sp2.get_label()])

plt.subplot(313)
x = np.linspace(0, 25000, 200)
sp3 = plt.fill_between(x , 0, norm_pdf(x, 12000, 3000),
                 color = "#7A68A6", lw = 3, alpha = 0.6,
                 label = "Trip price guess")
plt.autoscale(tight=True)
p3 = plt.Rectangle((0, 0), 1, 1, fc=sp3.get_facecolor()[0])
plt.legend([p3], [sp3.get_label()]);


$\ $위의 세 그래프는 각각 전체가격과, 제설기 가격, 토론토 여행권에 대한 사전확률분포를 보여준다. 이 사전확률분포들을 활용한 사후확률분포를 추정해보자.

사후확률분포

import pymc3 as pm

data_mu = [3e3, 12e3]

data_std =  [5e2, 3e3]

mu_prior = 35e3
std_prior =  75e2
with pm.Model() as model:
    true_price = pm.Normal("true_price", mu=mu_prior, sd=std_prior)

    prize_1 = pm.Normal("first_prize", mu=data_mu[0], sd=data_std[0])
    prize_2 = pm.Normal("second_prize", mu=data_mu[1], sd=data_std[1])
    price_estimate = prize_1 + prize_2

    logp = pm.Normal.dist(mu=price_estimate, sd=(3e3)).logp(true_price)
    error = pm.Potential("error", logp)


    trace = pm.sample(50000, step=pm.Metropolis())
    burned_trace = trace[10000:]

price_trace = burned_trace["true_price"]

figsize(12.5, 4)

import scipy.stats as stats

x = np.linspace(5000, 40000)
plt.plot(x, stats.norm.pdf(x, 35000, 7500), c = "k", lw = 2,
         label = "prior dist. of suite price")

_hist = plt.hist(price_trace, bins = 35, normed= True, histtype= "stepfilled")
plt.title("Posterior of the true price estimate")
plt.vlines(mu_prior, 0, 1.1*np.max(_hist[0]), label = "prior's mean",
           linestyles="--")
plt.vlines(price_trace.mean(), 0, 1.1*np.max(_hist[0]), \
           label = "posterior's mean", linestyles="-.")
plt.legend(loc = "upper left");


$\ $우리가 원래 전체가격에 대해서 믿고있던 평균 값은 35000이었다. 하지만 각 상품들의 사전확률분포를 포함하여 전체가격에 대한 사후확률을 추론한 결과, 바로 위의 그래프처럼 약 20000의 값으로 확률분포의 평균 값을 추정하였다. 이는 원래 가격이 35000에 비해 15000이나 낮아진 값이며, 각 상품이 갖고 있던 불확실성을 반영한 결과이다. 일반적인 빈도주의자 관점에서는 불확실성을 포함시키지 않기 때문에, 전체 가격이 단순히 \(\mu_1 + \mu_2 = 35000\)이라고 말할 수 있지만, naive Bayesian 관점에서는, 사후분포의 평균 값을 새로운 추정 값으로 결정하며 이 과정에서 자연스럽게 불확실성(표준편차들에 의한)이 반영됐다. 하지만 활용 가능한 여분의 정보가 더 있다. 이를 포함시키기 위해 loss function을 생각해보자.

손실함수(Loss function)

def showcase_loss(guess, true_price, risk = 80000):
    if true_price < guess:
        return risk
    elif abs(true_price - guess) <= 250:
        return -2*np.abs(true_price)
    else:
        return np.abs(true_price - guess - 250)

위와 같은 loss function이 제시됐을 때, 이 손실함수가 내포하고 있는 논리는 무엇일까? 처음에 이 쇼케이스의 진행 방식에 생략된 것들이 있는데 그것들은 다음과 같다.

  1. 우리가 실제가격(true price)보다 낮은 가격을 제시(guess)하고, 제시한 가격과 실제 가격의 차이가 250 이내라면, 실제 가격의 2배 를 보상으로 받는다.
  2. 사람마다 이 쇼케이스를 대하는 마음가짐이 다를텐데, 그 default risk를 risk로 표현한다. 누군가는 백만장자라서 해당 쇼케이스에서 아무것도 얻지 못해도 상관이 없을 것이다. 이런 경우의 risk는 0으로 볼 수 있고, 완전 반대의 경우는 그 risk가 매우 클 것이다. 이 경우에는 risk를 80000으로 설정했다. risk 변화에 따른 loss 최소화는 아래에서 더 얘기한다.
  3. 마지막 else 문에서는, 우리가 제시한 가격이 실제가격보다 낮고, 실제가격과 차이가 250 초과로 차이가 날 때인데, 왜 저렇게 loss를 잡는지는… 아직 잘 이해가 안된다. 해당 손실함수를 활용하여, risk에 따른 손실함수의 추이를 확인해보자.
def showdown_loss(guess, true_price, risk = 80000):
        loss = np.zeros_like(true_price)
        ix = true_price < guess
        loss[~ix] = np.abs(guess - true_price[~ix])
        close_mask = [abs(true_price - guess) <= 250]
        loss[close_mask] = -2*true_price[close_mask]
        loss[ix] = risk
        return loss

guesses = np.linspace(5000, 50000, 70)
risks = np.linspace(30000, 150000, 6)
expected_loss = lambda guess, risk: \
    showdown_loss(guess, price_trace, risk).mean()

for _p in risks:
    results = [expected_loss(_g, _p) for _g in guesses]
    plt.plot(guesses, results, label = "%d"%_p)

plt.title("Expected loss of different guesses, \nvarious risk-levels of \
overestimating")
plt.legend(loc="upper left", title="Risk parameter")
plt.xlabel("price bid")
plt.ylabel("expected loss")
plt.xlim(5000, 30000);

Optimization

기대손실을 최소화하기 위해 scipy.optimize 모듈에서 fmin 함수를 사용한다.

import scipy.optimize as sop

ax = plt.subplot(111)

for _p in risks:
    _color = next(ax._get_lines.prop_cycler)
    _min_results = sop.fmin(expected_loss, 15000, args=(_p,),disp = False)
    _results = [expected_loss(_g, _p) for _g in guesses]
    plt.plot(guesses, _results , color = _color['color'])
    plt.scatter(_min_results, 0, s = 60, \
                color= _color['color'], label = "%d"%_p)
    plt.vlines(_min_results, 0, 120000, color = _color['color'], linestyles="--")
    print("minimum at risk %d: %.2f" % (_p, _min_results))

plt.title("Expected loss & Bayes actions of different guesses, \n \
various risk-levels of overestimating")
plt.legend(loc="upper left", scatterpoints = 1, title = "Bayes action at risk:")
plt.xlabel("price guess")
plt.ylabel("expected loss")
plt.xlim(7000, 30000)
plt.ylim(-1000, 80000);

각 리스크 수준에서 얻게되는 최소손실 값은 다음과 같다.

minimum at risk 30000: 14723.45
minimum at risk 54000: 13500.92
minimum at risk 78000: 11900.78
minimum at risk 102000: 11649.08
minimum at risk 126000: 11649.08
minimum at risk 150000: 11329.30
© Zguu's Data Science Blog 2020
Search