Ch.16 § 16.1~16.3 심화 — Standard GLM Likelihoods·Working with GLMs·Weakly Informative Priors

1 개요

Ch.16 심화 시리즈의 첫 번째 편. § 16.1~16.3은 GLM의 이론적 뼈대 + 계산 엔진 + prior 철학을 다룬다.

• § 16.1 — Likelihood 층: Poisson·binomial·probit의 정확한 형태와 과분산 처리.
• § 16.2 — 계산 층: canonical link·offset·잠재 변수·정규 근사 (IWLS).
• § 16.3 — Prior 층: 분리(separation)라는 실무 문제와 Cauchy(0, 2.5) weakly informative prior의 근거.

직관: 세 절의 역할

Gelman의 구성 의도:

1. § 16.1 에서 “어떤 likelihood를 쓸 것인가” 를 정의 — 데이터 타입에 맞는 분포 선택.
2. § 16.2 에서 “어떻게 계산할 것인가” — conjugate 아닌 likelihood를 선형 근사해 Ch.14 기계로 돌림.
3. § 16.3 에서 “prior는 어떻게 둘 것인가” — 분리 문제를 피할 weakly informative 기본값 제안.

세 절이 “likelihood → 계산 → prior” 의 통합 workflow를 형성한다. 이 편이 이해되면 Ch.16 나머지 (16.4 NYC 검문·16.5 MRP·16.6 multinomial) 는 이 프레임의 특수 사례로 자연스럽게 읽힌다.

2 § 16.1 Standard GLM Likelihoods

2.1 GLM 3단 구조 재복습

\[ \eta = X\beta \quad \text{(linear predictor)} \]

\[ \mu = g^{-1}(\eta) \quad \text{(link)} \]

\[ y | X \sim p(y | \mu, \phi), \quad \mathbb{E}(y) = \mu \quad \text{(random)} \]

일반 likelihood (식 (16.1)):

\[ p(y | X, \beta, \phi) = \prod_{i=1}^n p(y_i | X_i \beta, \phi) \]

2.2 Continuous Data

Normal identity (Ch.14 선형 회귀):

\[ y_i \sim N(X_i \beta, \sigma^2), \quad g(\mu) = \mu \]

Normal log — 양수 \(y\) 의 곱셈적 효과:

\[ \log y_i \sim N(X_i \beta, \sigma^2) \]

→ \(y_i \sim \text{LogNormal}(X_i \beta, \sigma^2)\). 경제·소득·면적 데이터의 표준.

Gamma regression: 양수, right-skewed.

\[ y_i \sim \text{Gamma}(\nu, \nu / \mu_i), \quad \mathbb{E}(y_i) = \mu_i, \quad \mathrm{Var}(y_i) = \mu_i^2 / \nu \]

Canonical link: \(g(\mu) = 1/\mu\) (inverse). 실무에서는 log link이 더 흔함.

Weibull: 생존 분석 (Ch.22 survival models 전단계).

2.3 Poisson Regression — 식 (16.2) 유도

Setup: \(y_i \in \{0, 1, 2, \dots\}\) 계수 데이터, Poisson likelihood.

\[ y_i | \mu_i \sim \text{Poisson}(\mu_i) \Leftrightarrow p(y_i | \mu_i) = \frac{\mu_i^{y_i} e^{-\mu_i}}{y_i!} \]

Canonical link: \(g(\mu) = \log \mu\).

\[ \log \mu_i = X_i \beta \Leftrightarrow \mu_i = \exp(\eta_i), \quad \eta_i = X_i \beta \]

전체 likelihood:

\[ p(y | \beta) = \prod_{i=1}^n \frac{(\exp \eta_i)^{y_i} e^{-\exp \eta_i}}{y_i!} = \prod_{i=1}^n \frac{1}{y_i!} e^{-\exp(\eta_i)} \exp(y_i \eta_i) \quad \text{(16.2)} \]

Log-likelihood:

\[ L(\beta | y) = \sum_i \left[ y_i \eta_i - e^{\eta_i} \right] + \text{const} \]

(\(y_i!\) 는 \(\beta\) 무관 상수로 흡수.)

직관: 왜 Poisson regression에서 log link인가

수학적 편의: log link면 \(\mu_i > 0\) 자동 보장 (exp의 image). 다른 link 쓰면 \(\mu_i\) 가 음수 될 위험.

해석적 편의: \(\beta_j\) 해석이 “\(x_j\) 1 단위 증가 → \(y\) 의 \(e^{\beta_j}\) 배” (rate ratio). 예: 주 인구 → 범죄 건수 모델링 시 경제 지표 \(\beta = 0.2\) → “지표 1 단위 증가가 범죄율 \(e^{0.2} \approx 22\%\) 증가”.

이론적 근거: Log이 Poisson의 canonical link — exponential family의 natural parameter. Sufficient statistic 계산·MLE 수렴 이론이 깔끔하다.

2.4 Binomial — Logistic & Probit

Setup: \(y_i | n_i, \mu_i \sim \text{Bin}(n_i, \mu_i)\), \(n_i\) 알려짐.

Logistic (logit link):

\[ g(\mu_i) = \log \frac{\mu_i}{1 - \mu_i} = \eta_i \Leftrightarrow \mu_i = \frac{e^{\eta_i}}{1 + e^{\eta_i}} \]

Likelihood:

\[ p(y | \beta) = \prod_i \binom{n_i}{y_i} \left( \frac{e^{\eta_i}}{1 + e^{\eta_i}} \right)^{y_i} \left( \frac{1}{1 + e^{\eta_i}} \right)^{n_i - y_i} \]

Log-likelihood (binomial coefficient 흡수):

\[ L(\beta | y) = \sum_i \left[ y_i \eta_i - n_i \log(1 + e^{\eta_i}) \right] \]

Probit (probit link):

\[ g(\mu_i) = \Phi^{-1}(\mu_i) = \eta_i \Leftrightarrow \mu_i = \Phi(\eta_i) \]

Likelihood:

\[ p(y | \beta) = \prod_i \binom{n_i}{y_i} \Phi(\eta_i)^{y_i} (1 - \Phi(\eta_i))^{n_i - y_i} \]

2.5 Logit vs Probit — 실무 차이

측면	Logit	Probit
해석	Odds ratio	잠재 정규 임계값
Tail	더 두꺼움	더 얇음
계수 스케일	\(\beta_{\text{logit}} \approx 1.6 \beta_{\text{probit}}\)	—
Gibbs sampler	복잡 (Polya-Gamma 필요)	간단 (잠재 정규)
범주형 이론	—	잠재 연속 변수 해석 자연

경험 법칙: Logit이 해석·경제학 표준, probit이 잠재 변수·심리측정 표준.

2.6 Complementary Log-Log

\[ g(\mu) = \log(-\log(1 - \mu)) \]

비대칭: \(g(\mu) \neq -g(1-\mu)\). 즉 “성공” 과 “실패” 를 비대칭적으로 취급.

응용: Interval-censored 생존 데이터 (한 구간 내 사건 발생 확률).

2.7 Overdispersion — 과분산 확장

Poisson 기본: \(\mathrm{Var}(y) = \mu\). 실제 데이터는 종종 더 큰 분산 — “extra-Poisson variation”.

원인:

• 모델에 없는 변수의 영향.
• Cluster (같은 그룹 관측 간 상관).
• Sampling 외 randomness 원천.

해법 1 — Negative Binomial:

\[ y_i | \mu_i, \phi \sim \text{NegBin}(\mu_i, \phi), \quad \mathrm{Var}(y_i) = \mu_i + \mu_i^2 / \phi \]

\(\phi\) 크면 Poisson 수렴, 작으면 과분산 강함.

해법 2 — 계층 random effects:

\[ y_i | \beta, \epsilon_i \sim \text{Poisson}(\exp(X_i \beta + \epsilon_i)), \quad \epsilon_i \sim N(0, \sigma_\epsilon^2) \]

§ 16.4 NYC 경찰 검문 예제의 핵심. Gelman 선호.

Binomial 과분산 — Beta-binomial:

\[ y_i | n_i, \mu_i \sim \text{BetaBin}(n_i, \mu_i, \phi) \]

또는 \(y_i | \beta, \epsilon_i \sim \text{Bin}(n_i, \text{logit}^{-1}(X_i \beta + \epsilon_i))\), \(\epsilon_i \sim N(0, \sigma^2)\).

3 § 16.2 Working with GLMs

3.1 Canonical Link — Exponential Family 유도

Exponential family (Ch.2 재방문):

\[ p(y | \theta) = h(y) \exp(\theta \cdot T(y) - A(\theta)) \]

\(\theta\) 가 natural parameter, \(T(y)\) sufficient statistic, \(A(\theta)\) log-partition.

Canonical link 는 \(g\) s.t.

\[ g(\mu) = \theta \]

즉 평균 \(\mu\) 를 natural parameter \(\theta\) 로 직접 매핑.

분포	\(\theta\)	\(g(\mu)\)
Normal	\(\mu / \sigma^2\) (변수 분리)	\(\mu\) (identity)
Poisson	\(\log \mu\)	\(\log \mu\)
Binomial	\(\log \frac{\mu}{1-\mu}\)	logit
Gamma	\(-1/\mu\)	\(1/\mu\) (inverse)

이론적 이점:

• MLE의 충분통계량이 \(X^T y\) 로 단순.
• Fisher information이 구조적.
• IWLS 수렴 이론 표준.

실무적 유연성: Probit·complementary log-log 같은 non-canonical link도 많이 씀. Prior·해석 편의에 따라 선택.

3.2 Offset — 알려진 계수 강제

상황: Poisson에서 \(\mu_i = \lambda_i T_i\) 형태, \(T_i\) 는 관측의 노출 (exposure) — 시간, 인구, 면적 등.

\[ \log \mu_i = \log T_i + X_i \beta \]

\(\log T_i\) 가 계수 1 로 고정된 “fake predictor”. 이것이 offset.

구현: \(X\) 에 \(\log T_i\) column 추가하되 해당 계수를 1로 강제. statsmodels.GLM(y, X, offset=np.log(T)).

실무 예 (§ 16.4):

\[ \text{stops}_i \sim \text{Poisson}(\exp(\log \text{crimes}_i + \beta_0 + \beta_{\text{race}} x_i^{\text{race}})) \]

“이전 범죄 대비 검문 rate” 로 해석.

3.3 Interpreting Coefficients — Link Scale vs Original Scale

GLM 계수 \(\beta_j\) 는 link 스케일에서 “1 단위 변화당 \(\eta\) 의 1 단위 변화”를 의미.

원 스케일 효과: 비선형이라 기준점 \(x_0\) 에 의존.

\[ \Delta \mu = g^{-1}(X_0 \beta + (\Delta x) \beta) - g^{-1}(X_0 \beta) \]

예시: 로지스틱에서 \(\beta = 0.5\):

• \(p_0 = 0.5\) 일 때: \(p_0 \to 0.62\) (12%p 상승).
• \(p_0 = 0.05\) 일 때: \(p_0 \to 0.08\) (3%p 상승).
• \(p_0 = 0.95\) 일 때: \(p_0 \to 0.97\) (2%p 상승).

“같은 계수, 다른 효과” — 기준 확률에 따라 효과 크기가 크게 달라진다.

3.4 잠재 연속 변수 해석 — 식 (16.3)

Probit 모형:

\[ \Pr(y_i = 1) = \Phi(X_i \beta) \]

이는 다음과 수학적으로 동등 (식 16.3):

\[ u_i \sim N(X_i \beta, 1), \quad y_i = \mathbb{1}[u_i > 0] \]

증명:

\[ \Pr(y_i = 1) = \Pr(u_i > 0) = \Pr(X_i \beta + \epsilon > 0), \quad \epsilon \sim N(0, 1) \]

\[ = \Pr(\epsilon > -X_i \beta) = 1 - \Phi(-X_i \beta) = \Phi(X_i \beta) \]

3.5 잠재 변수의 실무 가치

직관: “관측된 binary 뒤의 연속 성향”

정치 설문에서 \(y_i\) = “공화당 지지 여부”. 내부적으로는 연속 성향 \(u_i\):

• \(u_i \gg 0\): 압도적 공화당 지지.
• \(u_i \ll 0\): 압도적 민주당 지지.
• \(u_i \approx 0\): 중도, 투표는 불안정.

관측은 \(\mathbb{1}[u_i > 0]\) 로 축약되지만, 모델은 연속 \(u_i\) 위에서 전개 → prior·계층 구조 설계가 자연스러움.

Gibbs sampler 이점:

1. \(u_i | \beta, y_i\): 절단 정규 (truncated normal).
   • \(y_i = 1\): \(u_i | \cdot \sim N(X_i\beta, 1), u_i > 0\).
   • \(y_i = 0\): \(u_i | \cdot \sim N(X_i\beta, 1), u_i < 0\).
2. \(\beta | u\): 단순 선형 회귀 (\(u\) 가 연속이므로 Ch.14).

이 “probit Gibbs sampler” (Albert-Chib 1993) 가 베이즈 discrete data 분석의 고전.

로지스틱에는 Polya-Gamma augmentation (Polson-Scott-Windle 2013) 이 유사한 역할.

3.6 Ordered Multinomial — Cut-points

\(y_i \in \{0, 1, \dots, K\}\). 잠재 변수 + 임계값:

\[ u_i \sim N(X_i \beta, 1), \quad y_i = k \text{ if } c_{k-1} < u_i < c_k \]

\(c_0 = -\infty, c_K = \infty\).

Identifiability: \(c_1 = 0\) 고정 (shift invariance 제거).

3.7 Prior 선택 4가지

1. Noninformative on \(\beta\)

\(p(\beta) \propto 1\). 빈도주의 MLE = posterior mode. 편리하지만 분리 문제 발생 (§ 16.3).

2. Conjugate — Pseudo-data

Ch.14 § 14.8 의 prior-as-data 트릭. 가상 관측 \((y_0, X_0, n_0)\) 을 원 데이터에 append → 가중 회귀로 계산. GLM에서는 정확히 conjugate는 아니지만 정규 근사 프레임에서 유효.

3. Nonconjugate — Normal on \(\beta\)

\[ \beta \sim N(\beta_0, \Sigma_\beta) \]

명시적 사전 지식 인코딩. 계산은 MCMC 또는 정규 근사.

4. Hierarchical

\[ \beta_j \sim N(\alpha, \tau^2), \quad \alpha, \tau \text{ 도 hyperparameter} \]

Ch.15 와 동일 구조. 여러 계수 자동 shrinkage.

3.8 정규 근사와 IWLS — 식 (16.4) 완전 유도

목표: GLM likelihood \(L(\beta | y)\) 를 선형 회귀 likelihood로 근사하여 Ch.14 기계 재사용.

전략: 각 관측의 log-likelihood \(L(y_i | \eta_i, \phi)\) 를 현재 추정 \(\hat\eta_i = X_i \hat\beta\) 주변에서 2차 Taylor 전개.

\[ L(y_i | \eta_i, \phi) \approx L(y_i | \hat\eta_i, \phi) + L'(y_i | \hat\eta_i) (\eta_i - \hat\eta_i) + \frac{1}{2} L''(y_i | \hat\eta_i) (\eta_i - \hat\eta_i)^2 \]

\(L''\) 은 음수 (log-likelihood가 오목) 이므로 \(-L''\) 정의하고 완성 제곱:

\[ L(y_i | \eta_i, \phi) \approx -\frac{1}{2\sigma_i^2}(z_i - \eta_i)^2 + \text{const} \]

여기서 pseudodata \(z_i\) 와 pseudovariance \(\sigma_i^2\):

\[ \boxed{ z_i = \hat\eta_i - \frac{L'(y_i | \hat\eta_i, \hat\phi)}{L''(y_i | \hat\eta_i, \hat\phi)}, \quad \sigma_i^2 = -\frac{1}{L''(y_i | \hat\eta_i, \hat\phi)} } \quad \text{(16.4)} \]

해석: 원 GLM이 근사적으로 가중 선형 회귀

\[ z \approx X\beta + \epsilon, \quad \epsilon_i \sim N(0, \sigma_i^2) \]

Ch.14 weighted regression (§ 14.7) 기계를 그대로 적용.

3.9 Binomial-Logistic 완전 유도

Log-likelihood:

\[ L(y_i | \eta_i) = y_i \eta_i - n_i \log(1 + e^{\eta_i}) \]

1차 미분:

\[ L'(y_i | \eta_i) = y_i - n_i \frac{e^{\eta_i}}{1 + e^{\eta_i}} = y_i - n_i \mu_i \]

(\(\mu_i = e^{\eta_i} / (1 + e^{\eta_i})\).)

2차 미분:

\[ L''(y_i | \eta_i) = -n_i \frac{e^{\eta_i}}{(1 + e^{\eta_i})^2} = -n_i \mu_i (1 - \mu_i) \]

(로지스틱 분산 \(n_i \mu_i (1 - \mu_i)\) 에 음의 부호.)

Pseudodata:

\[ z_i = \hat\eta_i - \frac{y_i - n_i \hat\mu_i}{-n_i \hat\mu_i (1 - \hat\mu_i)} = \hat\eta_i + \frac{y_i - n_i \hat\mu_i}{n_i \hat\mu_i (1 - \hat\mu_i)} \]

\(y_i / n_i \equiv p_i\) (observed proportion) 로 바꾸면:

\[ z_i = \hat\eta_i + \frac{p_i - \hat\mu_i}{\hat\mu_i (1 - \hat\mu_i)} \]

Pseudovariance:

\[ \sigma_i^2 = \frac{1}{n_i \hat\mu_i (1 - \hat\mu_i)} \]

직관: IWLS의 내부 메커니즘

\(z_i \approx \hat\eta_i + (\text{관측 비율과 예측 비율의 차이}) / \text{로지스틱 분산}\).

\(\sigma_i^2 \approx 1 / (n_i \cdot \mu(1-\mu))\) — \(\mu\) 가 0.5 근처면 작음 (정보 많음), 0 또는 1 근처면 큼 (정보 적음).

이것이 logistic regression의 weighting:

• 중간 확률 관측이 회귀에 큰 영향.
• 극단 확률 관측이 작은 영향.

왜 그런가? \(p \approx 0\) 또는 \(1\) 관측이면 “이미 거의 확실” → \(\beta\) 변화에 덜 민감 → 분산이 큼. 이 정보량의 차이를 가중치로 자동 반영.

3.10 IWLS 알고리즘

1. 초기값 \(\hat\beta^{(0)}\) (보통 상수항만).
2. 반복 (\(t = 0, 1, 2, \dots\)):
   1. \(\hat\eta^{(t)} = X\hat\beta^{(t)}\), \(\hat\mu^{(t)} = g^{-1}(\hat\eta^{(t)})\).
   2. Pseudodata \(z^{(t)}\), weights \(W^{(t)} = \text{diag}(1/\sigma_i^{2(t)})\).
   3. 가중 회귀: \(\hat\beta^{(t+1)} = (X^T W^{(t)} X)^{-1} X^T W^{(t)} z^{(t)}\).
3. 수렴까지 반복.

결과: 수렴 후 \(\hat\beta\) 는 MLE, \((X^T W X)^{-1}\) 는 asymptotic covariance.

Newton-Raphson과의 동등성: 이 알고리즘이 바로 \(dL/d\beta = 0\) 의 뉴턴-랩슨 해법과 동일. GLM IWLS가 “hidden Newton”.

3.11 정규 근사 사후 분포

수렴 후:

\[ \beta | y \approx N(\hat\beta, V_\beta), \quad V_\beta = (X^T W X)^{-1} \]

대규모 데이터에서는 이 근사로 충분. 소규모·계층 모형·nonconjugate prior에서는 정규 근사를 proposal로 쓴 MCMC 또는 VI/EP 필요 (Ch.13 § 13.7-13.8).

4 § 16.3 Weakly Informative Priors for Logistic Regression

4.1 분리(Separation)의 정의

완전 분리 (complete separation): \(X\) 의 어떤 선형 결합이 \(y\) 를 완벽 분류.

\[ \exists \beta^* : \; X_i \beta^* > 0 \iff y_i = 1 \; \forall i \]

이런 \(\beta^*\) 가 있으면 \(c \cdot \beta^*\) (\(c \to \infty\)) 도 모든 관측을 완벽 설명 → MLE 존재하지 않음 (\(\beta \to \infty\)).

유사분리 (quasi-separation): 일부 관측은 정확히 분리, 나머지만 혼합. MLE도 여전히 발산.

4.2 실제 사례 — Table 16.1 (1964 선거)

Gelman의 1960~1972 대통령 선거 설문 연구. 공화당 투표 여부를 sex·ethnicity·income 으로 예측.

1964년 데이터:

• 흑인 응답자 \(n = 87\).
• 그 중 공화당 지지자 \(= 0\) (Goldwater의 civil rights 반대로 100% 민주당).

결과: “black” 계수 추정치 \(= -16.83\), SE \(= 420.40\) (R glm). 이 값은 iterative fitting의 중단 시점에 의존하는 의미 없는 숫자. MLE는 실제로 \(-\infty\).

Table 16.1 전체:

연도	\(n\)	black 계수	SE
1960	875	-1.03	0.36
1964	1059	-16.83	420.40
1968	850	-3.64	0.59
1972	1518	-2.58	0.26

다른 연도는 합리적 추정 (\(\beta \in [-4, -1]\)), 1964만 극단값.

4.3 Figure 16.1 — Profile Likelihood

1964 데이터의 “black” 계수에 대한 profile likelihood (다른 모수 MLE 고정):

• \(\beta \to -\infty\) 로 단조 증가.
• Maximum 이 경계.
• 통상적인 MLE·SE 해석 불가.

왜 “분리” 가 실무에서 흔한가

단순 데이터 수 부족으로 보이지만, 실은 매우 흔한 문제:

1. 희귀 카테고리: 소수민족 × 특정 정치성 → 100% 한쪽으로 쏠림.
2. 극단적 지표: 고소득층이 모두 공화당 지지 등.
3. 많은 binary predictors: 변수 많을수록 perfect 분리 일어날 확률 증가.
4. 작은 표본: 관측 수 < 변수 수 × 카테고리 수.

전통적 “해법” (피해야 함):

• 문제 변수 제거 → 가장 강한 예측변수를 버림. 역설적으로 가장 유용한 변수가 분리 일으킴.
• 모델을 “알아서” 빼기 (stepwise). 결정이 불투명.

Gelman의 권장: Prior로 안정화. 아래 weakly informative 기법.

4.4 Weakly Informative Prior — Cauchy(0, 2.5)

Gelman-Jakulin-Pittau-Su (2008) 의 권장:

\[ \beta_j \sim \text{Cauchy}(0, 2.5), \quad j = 1, \dots, k \]

\[ \beta_0 \sim \text{Cauchy}(0, 10) \quad \text{(상수항)} \]

4.5 왜 Cauchy?

\(t\)-family 일반화:

\[ \beta_j \sim t_\nu(0, s) \]

• \(\nu = 1\): Cauchy (heavy tail).
• \(\nu = 7\): \(t_7\) (moderate tail).
• \(\nu = \infty\): Normal.

Heavy tail의 장점: “진짜 강한 효과” 가 있을 때 계수가 prior scale \(s\) 를 넘어 자유롭게 움직일 수 있음. Cauchy가 특히 관대.

0 중심의 장점: “대부분 효과는 작다” 라는 현실적 사전 지식. 분리 상황에서 0 근처로 shrinkage.

4.6 왜 Scale 2.5?

로지스틱 스케일의 의미:

• \(\beta = 1\): logit 1 단위 → odds ratio \(e^1 \approx 2.72\).
• \(\beta = 2.5\): odds ratio \(\approx 12.2\), \(p_0 = 0.5 \to p_0 = 0.92\).
• \(\beta = 5\): odds ratio \(\approx 148\), \(p_0 = 0.5 \to 0.993\).

\(|\beta| > 5\) 는 “거의 모든 것을 완벽히 분류” 에 해당. 실무에서 이런 효과는 매우 드물다.

Scale 2.5는 “한 표준편차 안에서 거의 완벽 분류까지 허용” 의 중간값. 강한 효과를 방해하지 않으면서 분리를 막을 만큼 약한 정보.

4.7 Figure 16.2 — 세 분포 비교

Gelman의 Figure 16.2는 다음 3개 곡선을 겹친다:

1. Cauchy(0, 2.5).
2. \(t_7(0, 2.5)\).
3. “반 개 성공 + 반 개 실패” 단일 binomial 관측의 likelihood for \(\beta\).

관찰: 세 곡선 모두 \(|\beta| < 5\) 을 선호. Cauchy가 tail이 가장 두꺼워 “큰 효과 허용” 에서 가장 관대. \(t_7\) 은 중간, 단일 관측 likelihood는 경계.

Gelman의 선택: Cauchy를 default. Tail 관대함 > 일관성.

4.8 권장 Workflow

1. Binary 예측 변수: 0/1 유지 (표준화 X).
2. 연속 예측 변수: 평균 0, 표준편차 0.5로 표준화.
3. Prior:

\[ \beta_0 \sim \text{Cauchy}(0, 10), \quad \beta_j \sim \text{Cauchy}(0, 2.5) \; (j \geq 1) \]

4. 계산: arm::bayesglm (R) 또는 PyMC/NumPyro.

직관: 왜 연속 변수를 SD=0.5로 표준화하나?

일반적 표준화는 SD = 1. 왜 0.5?

이유: binary 예측 변수 (0/1) 의 SD는 대략 0.5 (\(\sqrt{p(1-p)}\), \(p \approx 0.5\)). 연속 변수도 SD 0.5로 맞추면 binary와 연속이 같은 “단위 효과” 를 갖는다.

즉 Cauchy(0, 2.5) prior가 binary/continuous 구분 없이 동일한 scale에서 작용.

bayesglm 이 이 표준화를 자동 수행하는 이유.

4.9 이 Prior의 이론적 정당화

수학적 효과:

• MLE가 \(\infty\) 로 가도 prior이 “유한한 영역으로 끌어당김” → posterior mode finite.
• Posterior mean/median도 합리적 범위.
• 95% credible interval도 well-defined.

빈도주의 penalty 해석: Laplace approximation의 prior mode는 \(-\log\) prior penalty:

\[ -\log \text{Cauchy}(0, 2.5) \propto \log(1 + \beta^2 / 6.25) \]

이는 로그 penalty — Tikhonov (ridge, \(\beta^2\)) 와 LASSO (\(|\beta|\)) 사이의 중간. 큰 \(\beta\) 에서는 매우 약한 penalty (big effects allowed), 작은 \(\beta\) 에서는 부드러운 shrinkage.

Firth (1993) 의 bias-reduced logistic과 관련 — 두 방법 모두 penalty 기반 separation 해결.

4.10 다른 Prior 옵션

• Normal(0, \(s\)): Tail 얇아서 큰 효과 과도 축소. 권장 안 함.
• Laplace(0, \(s\)): LASSO. Sparsity 유발하지만 분리 문제에서 불안정.
• Student-\(t_\nu\), \(\nu \in [3, 7]\): Cauchy보다 덜 관대, 대부분 상황에서 OK.
• Horseshoe (Ch.14 § 14.6): 큰 효과는 그대로, 작은 효과는 0 근처. Sparse GLM에 적합.

Gelman의 default는 Cauchy이지만 \(t_7\) 도 합리적 대안.

5 Python 완전 구현

5.1 예제 1: 로지스틱 + Weakly Informative Prior (분리 상황)

import numpy as np
import pymc as pm
import arviz as az

rng = np.random.default_rng(1964)

# simulate 1964-style separation: all black respondents -> Democrat
n_total = 1000
female = rng.binomial(1, 0.5, n_total)
black = rng.binomial(1, 0.10, n_total)
income = rng.standard_normal(n_total) * 0.5  # standardized

# base rates
beta0_true = -0.3
beta_female_true = 0.2
beta_income_true = 0.4

# introduce separation: black => always Democrat (y = 0)
logit_p = (beta0_true + beta_female_true * female
           + beta_income_true * income)  # no black coefficient
logit_p[black == 1] = -50.0  # effectively guaranteed y=0
p = 1 / (1 + np.exp(-logit_p))
y = rng.binomial(1, p)

X = np.column_stack([female, black, income])


# Model 1: Flat prior (will fail)
with pm.Model() as flat_model:
    b0 = pm.Normal("b0", 0, 1e6)
    b = pm.Normal("b", 0, 1e6, shape=3)
    eta = b0 + X @ b
    pm.Bernoulli("y", p=pm.math.sigmoid(eta), observed=y)
    try:
        trace_flat = pm.sample(1000, tune=500, target_accept=0.95, chains=2)
    except Exception as e:
        print(f"Flat prior failed or produced extreme values: {e}")


# Model 2: Weakly informative Cauchy prior (Gelman et al. 2008)
with pm.Model() as cauchy_model:
    b0 = pm.Cauchy("b0", alpha=0, beta=10)
    b = pm.Cauchy("b", alpha=0, beta=2.5, shape=3)
    eta = b0 + X @ b
    pm.Bernoulli("y", p=pm.math.sigmoid(eta), observed=y)
    trace_cauchy = pm.sample(2000, tune=1000, target_accept=0.97, chains=4)


print("\n=== Cauchy prior result ===")
print(az.summary(trace_cauchy, var_names=["b0", "b"]))
print(f"\nCoefficient names: [female, black, income]")
print(f"With separation, 'black' coefficient should be pulled to finite range")

예상 출력:

=== Cauchy prior result ===
      mean    sd  hdi_3%  hdi_97%  r_hat
b0   -0.31  0.08   -0.46    -0.17   1.00
b[0]  0.19  0.14   -0.08     0.45   1.00
b[1] -5.84  2.68  -11.20    -1.73   1.00
b[2]  0.43  0.07    0.30     0.58   1.00

해석:

• 처음 3개 계수 (female, income) 는 참값 근처.
• “black” 계수는 \(-5.84\) (95% CI [-11.2, -1.7]) — flat prior였으면 \(-\infty\) 였을 것. Cauchy가 유한한 큰 음수로 끌어와 해석 가능.
• MCMC 수렴 안정적 (\(\hat{R} = 1.00\)).

5.2 예제 2: Poisson + IWLS Manual Implementation

import numpy as np


def iwls_poisson(X, y, offset=None, max_iter=20, tol=1e-6):
    """IWLS for Poisson regression, demonstrating equation (16.4)."""
    n, k = X.shape
    if offset is None:
        offset = np.zeros(n)

    # initialize: intercept only
    beta = np.zeros(k)
    beta[0] = np.log(y.mean() + 1e-3)

    for it in range(max_iter):
        eta_hat = X @ beta + offset
        mu_hat = np.exp(eta_hat)

        # Poisson: L'(eta) = y - mu, L''(eta) = -mu
        # pseudo-data and variance (eq. 16.4)
        z = eta_hat + (y - mu_hat) / mu_hat - offset
        weights = mu_hat  # 1/sigma_i^2 = mu_i for Poisson

        # weighted regression: beta = (X' W X)^{-1} X' W z
        W = np.diag(weights)
        XtWX = X.T @ W @ X
        XtWz = X.T @ W @ z
        beta_new = np.linalg.solve(XtWX, XtWz)

        if np.max(np.abs(beta_new - beta)) < tol:
            beta = beta_new
            break
        beta = beta_new

    # asymptotic covariance
    V_beta = np.linalg.inv(XtWX)
    return beta, V_beta, it + 1


# simulate Poisson regression
rng = np.random.default_rng(42)
n = 300
X = np.column_stack([np.ones(n), rng.standard_normal(n), rng.binomial(1, 0.3, n)])
beta_true = np.array([1.0, 0.5, -0.7])
mu = np.exp(X @ beta_true)
y = rng.poisson(mu)

beta_hat, V_beta, iters = iwls_poisson(X, y)

print(f"IWLS converged in {iters} iterations")
print(f"beta_hat = {beta_hat.round(3)}")
print(f"SE       = {np.sqrt(np.diag(V_beta)).round(3)}")
print(f"True     = {beta_true}")

예상 출력:

IWLS converged in 6 iterations
beta_hat = [1.017 0.493 -0.719]
SE       = [0.058 0.037 0.075]
True     = [ 1.   0.5 -0.7]

식 (16.4) 에 따라 정확히 MLE 도출. 수렴 6회만에 완료 (뉴턴-랩슨의 이차 수렴).

6 § 16.1~16.3 실전 체크리스트

§ 16.1 — Likelihood 선택

• 반응변수 \(y\) 의 유형: continuous / count / binary / binomial / proportion.
• 적절한 분포 (Normal·Poisson·Binomial·Gamma) 선택.
• Canonical link가 기본; 다른 link 쓸 특별 이유 있는가?
• 과분산 의심되면 negative binomial 또는 random effects.

§ 16.2 — 계산

• Offset 필요? (exposure·time·area)
• 표준 데이터 전처리: 연속 \(x\) 를 SD 0.5로 표준화.
• IWLS 수렴? (\(\hat\beta\) 변화 < 1e-6)
• 큰 \(n\) 에서는 정규 근사로 충분, 작은 \(n\) 이면 MCMC.
• 잠재 변수 parameterization 적절한가 (probit의 경우)?
• 계층 모형이면 Non-centered reparameterization.

§ 16.3 — Prior

• 분리 가능성 점검: 카테고리별 \(y\) 빈도 (0 또는 \(n\) 인 cell 있는가)?
• Flat prior 사용 시 MLE 수렴 확인.
• 분리 의심되면 Cauchy(0, 2.5) 권장.
• 상수항 prior는 더 넓게 (Cauchy(0, 10)).
• arm::bayesglm (R) 또는 PyMC에서 구현.
• 대안: \(t_7\), horseshoe, normal(0, 2.5) 등.
• Posterior mode·mean·median 비교로 prior 영향 점검.

7 관련 주제

선행 지식

• Ch.16 Overview — GLM
• Ch.14 § 14.7~14.10 — Weighted Regression·Augmented Data
• Ch.13 § 13.1~13.3 — Mode Finding·Laplace Approximation
• Ch.3 § 3.7 — Bioassay (로지스틱 초기 예제)

후속 주제

• § 16.4~16.5 — Overdispersed Poisson·MRP (예정)
• § 16.6~16.9 — Multinomial·Loglinear·연습 (예정)
• Ch.17 Robust Inference — \(t\) 오차, robit

관련 개념 (cross-category)

• Logistic Regression Fundamentals
• Poisson Regression in Practice
• Firth Penalized Logistic (Separation)
• Polya-Gamma Augmentation

8 참고문헌

• Gelman, A., Carlin, J. B., Stern, H. S., Dunson, D. B., Vehtari, A., & Rubin, D. B. (2013). Bayesian Data Analysis (3rd ed.), Ch.16 § 16.1~16.3. CRC Press.
• Nelder, J. A., & Wedderburn, R. W. M. (1972). Generalized Linear Models. JRSS A, 135, 370-384.
• McCullagh, P., & Nelder, J. A. (1989). Generalized Linear Models (2nd ed.). Chapman & Hall.
• Gelman, A., Jakulin, A., Pittau, M. G., & Su, Y.-S. (2008). A Weakly Informative Default Prior Distribution for Logistic and Other Regression Models. Annals of Applied Statistics, 2(4), 1360-1383.
• Albert, J. H., & Chib, S. (1993). Bayesian Analysis of Binary and Polychotomous Response Data. JASA, 88, 669-679.
• Polson, N. G., Scott, J. G., & Windle, J. (2013). Bayesian Inference for Logistic Models Using Pólya-Gamma Latent Variables. JASA, 108, 1339-1349.
• Firth, D. (1993). Bias Reduction of Maximum Likelihood Estimates. Biometrika, 80, 27-38.