tutorial_bayesian_logistic_regression_basic

Bayesian logistic regression 기초¶

베이지안 로지스틱 회귀는 머신 러닝의 일반적인 도구인 로지스틱 회귀에 대한 베이지안 대응 개념입니다. 로지스틱 회귀의 목표는 주어진 훈련 항목에 대해 1 또는 0을 예측하는 것입니다. 예를 들어, 증상과 개인 정보를 바탕으로 누군가가 아픈지 아픈지 예측할 수 있습니다.

로지스틱회기

베이지안 로지스틱 회기 01

Odds(승산)에서 Logistic 유도하기¶

$P$ : 사건이 발생할 확률 (이길 확률)
$1-P$ : 사건이 발생하지 않을 확률 (질 확률)
$0 \le P \le 1$
$Odds(P) = \frac{P}{1-P}$
$O \le Odds \le \infty$
logit은 Odds에 log를 취한 값으로 logit의 범위는 $-\infty \le \ln(Odds) \le \infty$
$\text{logit}(P) = \ln \left( \frac{P}{1-P} \right)$
$z = \text{logit}(P), \;\;-\infty \le z \le \infty $

$$z = \ln \left( \frac{P}{1-P} \right)\\ e^z = \frac{P}{1-P} \\ P = e^z (1 – p)\\ P = e^z – e^z P\\ P(1 + e^z) = e^z\\ P = \frac{e^z}{1 + e^z}\\ P = \frac{1}{1 + e^{-z}} $$

logistic(sigmoid) function : $P = \frac{1}{1 + e^{-z}}, \;\; -\infty \le z \le \infty$

$\text{Probability Density} = \frac{d}{dz}P = \frac{d}{dz}\left(\frac{1}{1 + e^{-z}}\right) = \frac{1}{e^z+2+e^{-z}}$

용어정리¶

회기(Regression): 종속 변수와 하나 이상의 독립 변수간의 관계를 모델링 하는 통계적 방법을 나타낸다. 가장 기본적인 회귀 분석인 선형 회귀는 독립 변수와 종속 변수 간의 선형관계를 가정한다. 여기서 "선형"은 결과(종속 변수)가 독립변수의 선형 조합으로 나타난다는 의미.

선형 회기: 종속변수와 하나 이상의 독립변수간의 관계가 선형(직선)이라고 가정하고 이를 모델링. 선형회기에서는 데이터가 어떤 직선 또는 평면(또는 그 이상의 차원에서는 초평면) 주위에 분포되어 있음을 가정하고 이 직선 또는 평면을 찾는것

$z = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + … + \epsilon$

여기서 $z$는 종속변수, $x_1,x_2,…,$는 독립 변수들, $\beta_0, \beta_1, \beta_2,…,$ 은 회귀계수, 그리고 $\epsilon$은 오차 항을 나타냄

로지스틱 회귀: 분류문제를 해결 하기 위한것으로 독립 변수의 선형 조합을 계산한 다음, 로지스틱(시그모이드) 함수를 통해 이 값을 확률로 변환한다. 따라서, 로지스틱 회귀는 선형 회기의 결과를 시그모이드 함수에 입력으로 사용하여 0과 1사이의 값(확률)을 생성하는 것으로 볼 수 있다.

$$ P =\hat{y}= \frac{1}{1+e^{-z}} = \frac{1}{1+e^{-(\beta_0 + \beta_1 x_1 + \beta_2 x_2 + … + \epsilon)}} $$

로지스틱 회귀(Logistic Regression)의 사용처는?¶

로지스틱(시그모이드) 함수는 0 ~ 1사이의 값을 가지므로 확률로 해석 할 수 있어 확률 추정에 적합하다. 즉, 특정 항목이 어떤 카테고리에 속할 확률이 얼마인지 예측하는데 사용

이진 분류 문제에 적합: 종속 변수가 이진(binary)인 경우 사용되는 회귀 방법. 예를 들어 스팸 메일 여부, 대출 기본 위험 여부등
확률 예측: 각 범주에 속할 확률을 예측하는 것이기 때문에, 단순한 클래스 레이블만 제공하는 것이 아니라 확률적 해석도 가능

로지스틱 회기에서 파라미터 추정시 파라미터는 무엇인가?¶

선형모델을 찾는데 사용되는 파라미터는 $\beta_0, \beta_1, \beta_2,…,$ 로 회귀계수임

$z = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + … $

가장 쉬운 방법은 MLE(Maximum Likelihood Estimation) – 최대우도(또는 최대가능도) – 가 되는 회귀계수를 찾는 것이다. 하지만 멀티모달이 되는 경우 어려울 수 있다.

i : 데이터 index

1차원 선형회귀의 경우 $$ z_i = \beta_0 + \beta_1 x_i $$ 2차원 선형회귀의 경우 $$ z_i = \beta_0 + \beta_{1} x_{i,1} + \beta_{2} x_{i,2}\\ \beta = \begin{bmatrix} \beta_0 \\ \beta_1 \\ \beta_2 \end{bmatrix}\\ x_i = \begin{bmatrix} 1 \\ x_{i,1} \\ x_{i,2} \end{bmatrix}\\ z_i = \beta^T x_i\\ \hat{y}_i = h_\beta(x_i) = \frac{1}{1+e^{-\beta^Tx_i}}\\ \hat{y}_i : \text{i번째 데이터의 예측값 또는 예측확률}\\ $$

5번째 데이터인 경우 $z_5 = \beta_0 + \beta_{1} x_{5,1} + \beta_{2} x_{5,2}$

$$ \hat{y}_5 = \frac{1}{1+e^{-{z_5}}} = \frac{1}{1+e^{-(\beta_0 + \beta_{1} x_{5,1} + \beta_{2} x_{5,2})}}\\ $$

MLE (Maximum Likelihood Estimation – 최대우도추정)¶

Likelihood 는 아래와 같다

$ \mathcal{L}(\beta|y) = P(y|x;\beta) = P(y_1,y_2,…,y_n|x_1,x_2,…,x_n;\beta)\\ \text{x가 iid 라면}\\ \mathcal{L}(\beta|y) = P(y|x;\beta) = \prod_{i=1}^n P(y_i|x_i;\beta)\\ $

$n$ : 데이터 수
$y_i$ : i 번째 데이터의 실제라벨{0,1}
$x_i$ : i 번째 데이터의 특성 벡터 (d 차원 데이터 라면 $x_i = (x_{i,1},x_{i,2},…,x_{i,d})$)
$\beta$ : 로지스틱 회귀계수(파라미터) 벡터 (d 차원 이라면 $\beta = (\beta_1,\beta_1,…,\beta_d)$)
$P(y_i|x_i;\beta)$ : 파라미터 $\beta$를 사용하여 i번째 데이터의 실제 라벨 $y_i$를 예측 할 확률 즉 i번째 데이터의 likelihood.

라벨 $y_i$가 나올 확률은 $\hat{y}_i$, 라벨 $y_i$가 나오지 않을 확률은 $1- \hat{y}_i$ 이다.
이를 하나의 식으로 표시 하면 $P(y_i) = \hat{y}_i^{y_i} \left( 1 – \hat{y}_i \right)^{y_i}$
베르누이 시행과 동일 성공/실패와 같이 라벨이 A 인가($y=1$) 아닌가($y=0$), $\hat{y}$는 라벨이 A 일 확률
예를 들면 동전이 앞면(라벨 A)이 나올 확률이 $\hat{y}$ 이고 앞면이 나오면 1, 뒷면이 나오면 0, 이경우 i번째 시행에서 앞면이 나올 확률을 구하는 것과 동일

\begin{align*} P(y_i|x_i;\beta) &= h_\beta(x_i)^{y_i} \left( 1 – h_\beta(x_i) \right)^{(1-y_i)}\\ &= \hat{y}_i^{y_i} \left( 1 – \hat{y}_i \right)^{y_i}\\ \end{align*}

여기서 실제라벨 $y_i = \{0,1\}$ 이다
※ $P(y_i|x_i;\beta), P(y_i|x_i,\beta)$ 는 동일하나 의미가 약간 다른데 $P(y_i|x_i;\beta)$는 $x_i$ 는 데이터 이고 ";" 뒤에 $\beta$가 모델 파라미터라는 것을 명확하게 나타낸다.

MLE의 목표는 likelihood를 최대화 하는 $\beta$를 찾는것¶

계산의 편의를 위해 로그 likehood 의 최대화를 $beta$를 찾아도 동일하다. \begin{align*} l(\beta|y) &= \log \mathcal{L(\beta|y)}\\ &=\sum_{i=1}^n \log P(y_i|x_i;\beta)\\ &= \sum_{i=1}^n [ y_i\log(h_\beta(x_i)) + (1-y_i) \log( 1 – h_\beta(x_i))]\\ \end{align*} $ \hat{\beta}_{MLE} = \arg \max_\beta l(\beta) $

베이지안 로지스틱 회귀¶

베이지안 로지스틱 회귀는 로지스틱 회기 모델에 베이즈 추론을 적용한 것으로 로지스틱 회기와 마찬가지로 범주형 종속변수(특히 이진형)를 예측하는데 사용되지만, 파라미터 추정 방법에 차이가 있다.

로지스틱 회기와 베이지안 로지스틱 회기의 주요 차이점¶

1. 추정방법:¶

로지스틱 회기: MLE(Maximun Likelihood Estimation)을 사용하여 파라미터를 추정
베이지안 로지스틱 회기: 파라미터의 사전분포와 데이터에 기반한likelihood를 사용하여 파라미터의 사후 분포를 계산

2. 결과¶

로지스틱 회기: 파라미터의 단일 최적값을 제공
베이지안 로지스틱 회기: 파라미터에 대한 사후분포를 제공하므로, 불확실성에 대한 정보도 함께 얻음

3. 규제(Reqularization)¶

머신 러닝과 통계 모델에서 오버피팅(overfitting)을 방지하기 위한 기술로 모델의 복잡성을 제한하여 모델이 학습 데이터에 너무 과도하게 적합도는 것을 방지

로지스틱 회기: L1 또는 L2 규제와 같은 규제를 사용하여 과적합을 방지 할 수 있다.
베이지안 로지스틱 회귀: 사전분포가 규제의 역할을 할 수 있다. 예를 들어, 파라미터에 대한 사전분포로 정규분포를 사용하는 경우 L2 규제와 유사한 효과를 가질 수 있다.

베이지안 로지스틱 회귀의 주요 장점 중 하나는 파라미터의 불확실성을 직접적으로 모델링하고 추정할 수 있다는 것이다. 이 불확실성은 예측 간격이나 신회 간격을 생성하는데 사용될 수 있다.

참고¶

L1 Regularization(Lasso Regression): 모델의 가중치에 절대값을 기반으로하는 패널티를 추가하여 일부 가중치가 0 이 되어 특성 선택의 효과 있음
L2 Regularization(Ridge Regression): 모델의 가중치의 제곱에 기반하는 패널티를 추가 하기 때문에 가중치가 작아지지만 완전히 0이 되지 않음

Reqularization은 일반적으로 모델의 손실함수(loss function)에 패널티 항을 추가하여 구현되며 패널티 항의 크기는 하이퍼파라미터로 제어되며, 이 값이 크면 패널티의 영향이 커져 모델의 복잡성이 감소한다.

베이지안 로지스틱 회기에서 사후분포(posterior)¶

베이지안 로지스틱 회귀에서는 로지스틱 회기의 파라미터($\beta$)에 대한 사전분포를 설정하고, 테이터를 사용하여 이 파라미터에 대한 사후분포를 계산

$$ p(\beta|X, y) \propto p(y|X,\beta) \times p(\beta) $$

여기서

$p(\beta|X, y)$ : 데이터 $X, y$에 대한 $\beta$의 사후분포
$p(y|X, \beta)$ : 데이터에 대한 likelihood 함수로, 로지스틱 함수를 사용하여 $y$가 주어진 $X$와 $\beta$에 대해 어떻게 분포 되는지 나타낸다.
$p(\beta)$ : $\beta$에 대한 사전분포로, 통상적으로 정규분포나 다른 분포를 사용하여 특정 값에 대한 우리의 믿음(또는 무지)을 나타낸다.

베이지안 로지스틱 회귀의 중요한 부분은 사후분포의 정확한 형태를 직접 구할 수 없다는 것이다. 따라서 수티적인 방법(예, MCMC, 변분 베이즈 추론)을 사용하여 이 분포를 근사한다.

베이지안 로지스틱 회기의 파라미터 추정의 간단예¶

1. 모델 설정¶

데이터는 아래와 같다.

$x$ : 설명 변수 (예: 학생의 공부 시간)
$y$ : 반응 변수 (예: 합격 여부, 1 또는 0)

로지스틱 회귀 모델은 아래와 같음

$ p(y = 1| x, \beta_0, \beta_1) = \frac{1}{1+\exp(-(\beta_0 + \beta_1 x))} $ 여기서, $\beta_0$은 절편, $\beta_1$은 기울기

2. 사전분포 설정¶

파라미터 $\beta_0$와 $\beta_1$에 대한 사전 분포로 표준 정규분포로 가정

$\beta_0,\beta_1 \sim \mathcal{N}(0,1)$

3. 데이터 기반의 사후 분포 계산¶

$ p(\beta_0, \beta_1|x,y) \propto p(y|x,\beta_0,\beta_1) \times p(\beta_0, \beta_1) $

여기서 $p(y|x,\beta_0,\beta_1)$는 likelihood 이며 로지스틱 회기모델을 기반으로 한다. $p(\beta_0, \beta_1)$는 사전 분포

4. 실제 데이터의 예¶

학생 A가 3시간 공부했을 때 합격할 확률을 알고 싶다.

$ p(y=1|x=3, \beta_0,\beta_1) = \frac{1}{1+exp(-(\beta_0 + 3\beta_1))} $

하지만, $\beta_0$ 와 $\beta_1$은 확실하지 않기 떄문에, 이들에 대한 사후 분포를 사용하여 여러 가능한 값들에 대해 위의 확률을 계산하고 평균값을 얻는다.

이러한 계산은 보통 MCMC와 같은 방법을 사용하여 수행되며, 결과적으로 학생 A가 3시간 동안 공부했을 때 합격할 확률의 분포를 얻게 된다.

In [1]:

using Turing
using StatsPlots
using Distributions
using KernelDensity
using StatsFuns
using Random

In [2]:

Random.seed!(1234)

Out[2]:

TaskLocalRNG()

In [3]:

chains = let
    @model function logistic_regression_model(x::Vector{Float64},y::Vector{Int64})
        # 사전분포 정의
        β₀ ~ Normal(0,1)
        β₁ ~ Normal(0,1)

        # Likelihood 정의
        for i = 1:length(x)
            # p = 1 / (1 + exp(-(β₀ + x[i]*β₁)))
            p = logistic(β₀ + β₁* x[i])
            y[i] ~ Bernoulli(p)
        end
    end
    x = [0.0,2.0, 3.0, 5.0, 6.0, 7.0]  # 예: 학생들의 공부 시간
    y = [0, 0, 1, 0, 1, 1]  # 예: 합격 여부 (1: 합격, 0: 불합격)
    model = logistic_regression_model(x,y)
    chains = sample(model,NUTS(),1000)
    chains
end

┌ Info: Found initial step size
└   ϵ = 0.8
Sampling: 100%|█████████████████████████████████████████| Time: 0:00:00

Out[3]:

Chains MCMC chain (1000×14×1 Array{Float64, 3}):

Iterations        = 501:1:1500
Number of chains  = 1
Samples per chain = 1000
Wall duration     = 9.92 seconds
Compute duration  = 9.92 seconds
parameters        = β₀, β₁
internals         = lp, n_steps, is_accept, acceptance_rate, log_density, hamiltonian_energy, hamiltonian_energy_error, max_hamiltonian_energy_error, tree_depth, numerical_error, step_size, nom_step_size

Summary Statistics
  parameters      mean       std      mcse   ess_bulk   ess_tail      rhat   e ⋯
      Symbol   Float64   Float64   Float64    Float64    Float64   Float64     ⋯

          β₀   -0.6356    0.8231    0.0391   443.3511   415.5404    1.0034     ⋯
          β₁    0.3136    0.2875    0.0134   456.6196   333.5096    0.9993     ⋯
                                                                1 column omitted

Quantiles
  parameters      2.5%     25.0%     50.0%     75.0%     97.5% 
      Symbol   Float64   Float64   Float64   Float64   Float64 

          β₀   -2.1207   -1.2120   -0.6402   -0.0538    0.9919
          β₁   -0.2284    0.1247    0.3085    0.4840    0.9186

In [4]:

let
    plot(chains) |> display    
end

In [5]:

let
    density(chains[:lp], legend=false) |> display
    idx = argmax(chains[:lp])
    chains[:lp][idx] |> display
end

-5.4359280771505905

베이지안 추론으로 구한 파라미터 $\beta_0, \beta_1$을 사용하여 공부한 시간에 따른 합격 확률 계산¶

MAP(Maximum A Posteriori estimatiom)를 통한 합격 확률 계산
위 코드에서 Turing에서 계산한 chains[:lp] 는 로그 사후 분포로서 이 값이 최대가 되는 값의 $\beta_0,\beta_1$을 사용하여 합격 확률 계산

In [6]:

let
    idx = argmax(chains[:lp])
    β0,β1 = chains["β₀"],chains["β₁"]
    β0_MAP = β0[idx]
    β1_MAP = β1[idx]
    x = 3 # 공부한 시간
    
    p = logistic(β0_MAP + x * β1_MAP)
    
    println("$x 시간 공부하는 경우 합격 확률 = $p")
end

3 시간 공부하는 경우 합격 확률 = 0.5465004041333105

x 시간 공부 했을 때 합격 예측 확률을 모든 파라미터 샘플에 대해 계산¶

In [7]:

let
    β0,β1 = chains["β₀"],chains["β₁"]
    x = 3
    probs = logistic.(β0 + β1*x )
    
    # density(probs, xlabel="Probability", 
    #    ylabel="Probability Density", legend=false)
    
    mean_prob = mean(probs)
    println("$x 시간 공부하는 경우 평균 합격 확률 = $mean_prob")
    
    # 합격확률 의 density 구하기 
    pdf = kde(probs[:,1])
    
    # 가장 가능성이 높은 합격확률 구하기
    max_likelihood_of_probability = pdf.x[argmax(pdf.density)]
    println("$x 시간 공부하는 경우 가장 가능성이 높은 합격 확률 = $max_likelihood_of_probability")
    
    plot(pdf, title="Probability density of \nprobability of passing with $x hours of study",
        titlefontsize=12,
        xlabel="Probability", 
        ylabel="Probability Density of ", legend=false) |> display
end

3 시간 공부하는 경우 평균 합격 확률 = 0.5642519735551992
3 시간 공부하는 경우 가장 가능성이 높은 합격 확률 = 0.6471715599760544

In [ ]:

[Julia, Turing] Tutorial –Bayesian logistic regression 기초