Quasi-likelihood 심화 결과와 연습 — 11개 연습문제의 통합 해설

1 개요 — Ch.9 연습문제의 구조

McCullagh & Nelder §9.8 의 연습문제 11 개는 단순 연습이 아니라 Ch.9 본문에서 생략된 유도·보조정리·확장 예제 의 저장고다. 각 문제를 고립된 계산으로 풀면 놓치기 쉬운 교집합이 있는데, 이 포스트는 문제들을 주제 그룹 으로 묶어 왜 이들이 §9.8 에 모여 있는지 맥락을 드러낸다.

1.1 그룹별 조감

그룹	문제	주제	Ch.9 본문 연결
A. 지수족 구조	9.1	Poisson-Gamma 혼합의 세 분산 형태	§9.2, §9.6 의 과산포 배경
B. Voter 심화	9.2, 9.3	조건부 우도 분해와 정보행렬 rank	§9.3.3 수치 결과 재현
C. 수학 보조정리	9.4, 9.7	선적분·Löwner 역부등식	§9.3.2, §9.5 증명 도구
D. 기하 모형 확장	9.5	원에서 타원으로	§9.4.3 Avebury 확장
E. 존재성 조건	9.6	\(V = DRD\) 형태의 QL 존재성	§9.3.2 경로 독립성
F. 강건성	9.8, 9.9	CV 상수 모형에서 정규 가정의 효과	§9.2, §9.5 최적성 범위
G. Mantel-Haenszel	9.10, 9.11	추정함수로 도출한 MH 추정량	§9.4.2 결합·§7.4 응용

각 그룹을 순서대로 짚되, 계산의 공정보다 왜 이 결과가 Ch.9 의 뼈대를 보강하는지 의 관점에서 해설한다.

---

2 그룹 A — Poisson-Gamma 혼합과 분산 함수 (Ex 9.1)

2.1 문제

조건부 \(Y \mid M = m \sim \text{Poisson}(m)\), 그리고 \(M \sim \text{Gamma}(\alpha\nu, \nu)\) 로 \(E[M]=\alpha\nu, \operatorname{CV}(M)=\nu^{-1/2}\). 무조건 평균·분산을 구하면 \[ E[Y] = \mu = \alpha\nu, \qquad \operatorname{var}(Y) = \alpha\nu + \alpha^2\nu = \mu + \mu^2/\nu. \]

2.2 세 가지 시나리오

시나리오	가정	분산 함수	지수족?
(a) \(\nu_i = \nu\) (known)	Gamma 형상이 고정	\(V(\mu) = \mu + \mu^2/\nu\)	Yes — NB 지수족
(b) \(\alpha_i = \alpha\)	Gamma 스케일이 고정	\(V(\mu) = \phi\mu\), \(\phi=1+\alpha\)	No — 과산포 포아송
(c) \(\alpha_i = \theta + \psi\mu_i, \nu_i^{-1} = \psi + \theta\mu_i^{-1}\)	두 모수가 \(\mu_i\) 에 의존	\(V(\mu) = \mu + \theta\mu + \psi\mu^2\)	No

2.3 왜 이 문제가 중요한가

같은 “과산포 포아송” 이 세 개의 다른 모형

• 데이터에서 \(\operatorname{var}(Y) > E[Y]\) 를 관찰했을 때 “과산포” 라는 용어는 모형을 특정하지 않는다.
• 동일한 \(V(\mu) = \phi\mu\) 를 가정하더라도 (b) 는 엄격한 지수족이 아니므로 정확 우도를 쓸 수 없다 — 준우도 접근이 유일한 대안.
• 반면 (a) 는 음이항(NB) 으로 정확 우도가 존재 → MLE 가능.
• 3. 는 2차 분산 함수 (\(V(\mu) = \mu(1+\theta) + \psi\mu^2\)) 로 NB 와 포아송 중간 형태.

즉 분산 함수만 보고 모형을 구분할 수 없다 — 정확 우도와 준우도 중 무엇을 쓸지의 결정은 생성 과정 구조까지 알아야 한다.

2.4 실무적 해석

\(V(\mu) = \mu + \mu^2/\nu\) 형태가 음이항 GLM 의 기본 분산이며, Python statsmodels.NegativeBinomial, R MASS::glm.nb() 가 이것을 사용. \(\phi\mu\) 형태는 family = quasipoisson 에서 쓰인다 — 정확 우도 없이 \(Q\) 만 사용.

---

3 그룹 B — Voter Transition 의 우도 분해와 정보 (Ex 9.2-9.3)

3.1 Ex 9.2: 우도의 조건부 분해

§9.3.3 의 Voter Transition 우도를 아래 꼴로 변형할 수 있음을 요구: \[ \mathcal{L}(\pi_1, \pi_2, \psi) = (1-\pi_1)^{m_1} \pi_2^{y}(1-\pi_2)^{m_2 - y} \cdot P_0(\psi; m_1, m_2, y), \] 여기서 \(P_0\) 는 §7.3.2 의 비중심 초기하 확률.

분해의 의미: “주변+조건부” 우도

첫 세 인수 \((1-\pi_1)^{m_1} \pi_2^y (1-\pi_2)^{m_2-y}\) 는 이항 주변 우도 (각 블록의 총합에 관한), \(P_0(\psi)\) 는 비중심 초기하 조건부 우도 (오즈비 \(\psi\) 에 관한).

로짓 \(\lambda_i = \operatorname{logit}(\pi_i)\) 의 도함수로 쓰면 \[ \frac{\partial \ell}{\partial \lambda_1} = \sum\{\kappa_1(\psi) - m_1 \pi_1\}, \qquad \frac{\partial \ell}{\partial \lambda_2} = \sum\{y - \kappa_1(\psi) - m_2\pi_2\}, \] 여기서 \(\kappa_1(\psi)\) 는 비중심 초기하 평균.

이 두 스코어 식은 EM 알고리즘의 M-step 에 해당: \(\kappa_1(\psi)\) 가 “누락된 관측값의 조건부 기댓값” 역할을 하며, 주어진 \(\psi\) 에서 \(\pi_i\) 를 업데이트한 후 \(\psi\) 를 갱신하는 반복 구조.

3.2 Ex 9.3: 정보행렬의 rank 문제

\(\theta = (\pi_0, \lambda_1 - \lambda_2)\) 에 대한 Fisher 정보 행렬: \[ I_\theta = \frac{1}{V_\bullet^2}\sum_i\begin{pmatrix} m_1 V_1 + m_2 V_2 & (m_1-m_2)V_1 V_2 \\ (m_1-m_2)V_1 V_2 & (m_1 V_1 + m_2 V_2)V_1 V_2 - \kappa_2 V_\bullet^2\end{pmatrix}, \] 여기서 \(V_k = \pi_k(1-\pi_k)\), \(\kappa_2\) = 초기하 분산.

3.3 핵심 관찰

\(m_{i1} = m_{i2}\) 일 때에도 rank 2

• 대각 원소 \(m_1 V_1 + m_2 V_2 > 0\), off-diagonal 이 \((m_1-m_2)V_1 V_2\) 이므로 \(m_1 = m_2\) 이면 off-diagonal \(= 0\).
• 그럼에도 두 번째 대각 성분 은 \(\kappa_2\) (초기하 분산) 때문에 0 이 되지 않는다.
• 즉 \(\lambda_1 - \lambda_2\) 가 설계행렬 \(X\) 의 rank 2 가 아닌 한 일관 추정되지 않는다.

3.4 §9.3.3 와의 연결

§9.3.3 본문에서는 준우도 추정치 \(\widehat{\pi} = (0.3629, 0.8371)\) 와 MLE \((0.467, 0.733)\) 가 달랐고, 후자는 “Fisher 정보가 rank 2 여서 가능” 하다고 했다. Ex 9.3 는 이 주장의 수학적 근거 를 제공한다 — 준우도는 \(m_1 V_1 + m_2 V_2\) 만 활용하나 MLE 는 \(\kappa_2\) 항까지 활용한다.

---

4 그룹 C — 수학 보조정리 (Ex 9.4, 9.7)

4.1 Ex 9.4: 선적분 항등식

직선 경로 \(t(s) = b + s(c-b)\), \(s \in [0,1]\) 을 따라 \[ \int_b^c t^\top A \, dt(s) = \frac{1}{2}(c+b)^\top A (c-b). \]

4.2 대칭성 판정

\(b \to c \to d \to b\) 의 닫힌 삼각형 경로에서 적분이 0 임은 \(A = A^\top\) 와 동치. 역방향도 성립: 모든 닫힌 고리 적분이 0 이면 \(A\) 는 대칭.

경로 독립성 ≡ 대칭성

§9.3.2 에서 준우도의 존재성이 \(V^{-1}\) 행렬의 특정 경로적분이 경로 독립 일 때만 성립함을 보았다. Ex 9.4 는 이를 선형대수적 언어 — “적분 경로 \(A\) 가 대칭일 때만 경로 독립” — 로 정식화.

\(\nabla Q = -D^\top V^{-1}(y-\mu)\) 가 어떤 스칼라 \(Q\) 의 기울기가 되려면 Jacobian \(\partial(\nabla Q)/\partial\mu\) 가 대칭이어야 한다. Poincaré lemma 의 이산 버전.

4.3 Ex 9.7: Löwner 역부등식

두 양정치 행렬 \(A, B\) 에 대해 \[ A - B \succeq 0 \implies B^{-1} - A^{-1} \succeq 0. \]

4.4 증명 스케치

\(A - B \succeq 0\) 이면 대각화 가능한 \(A^{-1/2} B A^{-1/2}\) 의 고유값이 \(\leq 1\). 따라서 \((A^{-1/2} B A^{-1/2})^{-1} = A^{1/2} B^{-1} A^{1/2}\) 의 고유값은 \(\geq 1\). 이는 \(A^{1/2} B^{-1} A^{1/2} - I \succeq 0\) 을 의미하며, 양변을 \(A^{-1/2}\) 로 샌드위치하면 \(B^{-1} - A^{-1} \succeq 0\).

4.5 §9.5 와의 직접 연결

§9.5 증명의 결정적 도구

§9.5 에서 “준스코어가 최적” 을 증명할 때 공분산 행렬 비교 \(\operatorname{cov}(\widetilde{\beta}) \succeq \operatorname{cov}(\widehat{\beta})\) 를 직접 다루지 못하고 정밀도 행렬 차이 \(\{\operatorname{cov}(\widehat{\beta})\}^{-1} - \{\operatorname{cov}(\widetilde{\beta})\}^{-1} \succeq 0\) 로 뒤집었다. 이 뒤집기의 정당성이 바로 Ex 9.7.

즉 Löwner 순서는 역함수에서 뒤집힌다 — 이 단순한 사실이 없으면 §9.5 의 증명 구조 전체가 성립하지 않는다.

---

5 그룹 D — 기하 모형 확장: 타원 적합 (Ex 9.5)

5.1 문제 설정

§9.4.3 Avebury 거석환이 원이었다면, Ex 9.5 는 타원 으로 일반화한다: \[ \begin{aligned} Y_{i1} &= \omega_1 + \rho R_i \cos\epsilon_i \cos\phi - \lambda R_i \sin\epsilon_i \sin\phi, \\ Y_{i2} &= \omega_2 + \rho R_i \cos\epsilon_i \sin\phi + \lambda R_i \sin\epsilon_i \cos\phi, \end{aligned} \] 모수: 중심 \((\omega_1, \omega_2)\), 반축 \(\rho, \lambda\), 회전각 \(\phi\).

5.2 기본 추정함수

좌표를 회전하여 축정렬 좌표 \(X_{i1}, X_{i2}\) 로 변환: \[ X_{i1} = (Y_{i1}-\omega_1)\cos\phi + (Y_{i2}-\omega_2)\sin\phi = \rho R_i \cos\epsilon_i, \] \[ X_{i2} = -(Y_{i1}-\omega_1)\sin\phi + (Y_{i2}-\omega_2)\cos\phi = \lambda R_i \sin\epsilon_i. \]

기본 추정함수: \[ g_i = R_i - 1 = \sqrt{X_{i1}^2/\rho^2 + X_{i2}^2/\lambda^2} - 1. \]

5.3 계수 행렬 \(D\) (9.14)

5 개 모수 \((\omega_1, \omega_2, \rho, \lambda, \phi)\) 각각에 대해:

모수	\(D_{ir}\)
\(\omega_1\)	\(\cos\epsilon_i \cos\phi/\rho - \sin\epsilon_i\sin\phi/\lambda\)
\(\omega_2\)	\(\cos\epsilon_i\sin\phi/\rho + \sin\epsilon_i\cos\phi/\lambda\)
\(\rho\)	\(\cos^2\epsilon_i/\rho\)
\(\lambda\)	\(\sin^2\epsilon_i/\lambda\)
\(\phi\)	\((\rho - \lambda)\cos\epsilon_i\sin\epsilon_i\)

5.4 식별성의 기하학적 직관

\(\phi\) 모수가 \(\rho = \lambda\) 일 때 사라지는 이유

\(D_{i5} = (\rho-\lambda)\cos\epsilon_i\sin\epsilon_i\). 원 (\(\rho = \lambda\)) 이면 \(D_{i5} = 0\) 이 되어 회전각 \(\phi\) 는 식별 불가 — 원은 회전해도 같은 모양이므로 자연스러운 결과.

타원 적합의 실무 주의점: 데이터가 거의 원에 가깝거나 \(\rho, \lambda\) 추정이 수치적으로 가까우면 \(\phi\) 표준오차가 폭발한다. §9.4.3 의 “얕은 원호에서 모수 상관” 과 유사한 식별성 붕괴.

5.5 응용

고고학(거석 유적), 천문학(행성 궤도), 컴퓨터 비전(홍채 모델링) 등에서 점 집합에 타원을 맞추는 표준 방법. 기존 최소제곱법보다 추정함수 접근은 각도 오차 구조에 가정을 덜 요구한다.

---

6 그룹 E — 준우도 존재성 조건 (Ex 9.6)

6.1 문제

공분산 행렬이 다음 꼴로 쓰인다: \[ V = DRD, \quad D = \operatorname{diag}\{D_1(\mu_1), \ldots, D_n(\mu_n)\}, \] \(R\) 은 \(\mu\) 에 무관. 이때 준우도가 존재하려면:

결론: \(V_{ij}\) 가 \(\mu\) 에 무관해야 (\(i \neq j\))

즉 \(R\) 이 대각이거나 \(D\) 가 \(\mu\) 에 무관 이어야 한다.

6.2 유도 이유

§9.3.2 의 경로 독립 조건은 \(V^{-1}\) 의 성분이 만족해야 할 대칭성 조건을 준다. \(V = DRD\) 분해에서 \[ V_{ij} = D_i(\mu_i) R_{ij} D_j(\mu_j) \quad (i \neq j) \] 이 \(\mu\) 에 의존하면 대칭 조건이 위배되어 준우도 함수 \(Q(\mu; y)\) 가 전역적으로 존재하지 않는다.

6.3 실무 해석

\(R\)	\(D\)	QL 존재성
대각 (독립)	임의	Yes (표준 GLM)
비대각 (상관)	\(\mu\) 무관 (상수)	Yes (그러나 드물다)
비대각	\(\mu\) 의존	No — GEE 접근 필요

GEE 가 필요한 이유

종단 데이터는 일반적으로 \(R\) 이 비대각(시점 간 상관) 이고 \(D\) 가 \(\mu\) 의존(\(V(\mu)=\mu^2\) 같은). 따라서 준우도 함수 자체가 존재하지 않음 — GEE 는 “추정방정식” 만 사용하고 우도는 포기하는 방법.

이것이 Liang-Zeger (1986) GEE 접근의 이론적 출발점이며, §9.3 본문과 Ex 9.6 이 함께 이 필요성을 정당화한다.

---

7 그룹 F — 정규 가정 오류 하의 강건성 (Ex 9.8-9.9)

7.1 Ex 9.8: CV 상수 모형에서의 준우도

설정: \(Y_i\) 독립, \(\operatorname{var}(Y_i) = \sigma^2 \mu_i^2\), \(\log\mu_i = \beta_0 + \beta_1(x_i - \bar{x})\).

7.2 결론

\(\widehat{\beta}_0, \widehat{\beta}_1\) (준우도 추정) 이 비상관 이고: \[ \operatorname{var}(\widehat{\beta}_0) = \frac{\sigma^2}{n}, \qquad \operatorname{var}(\widehat{\beta}_1) = \frac{\sigma^2}{\sum(x_i - \bar{x})^2}. \]

\(x\) 를 중심화한 것이 비상관의 핵심 — \(\sum(x_i - \bar{x}) = 0\) 이 두 스코어의 교호 정보를 0 으로 만든다.

7.3 Ex 9.9: 정규 가정 MLE 의 오적용

같은 데이터 에 대해 (분산이 실제로 \(\sigma^2 \mu^2\) 인데도) 정규 가정 으로 MLE \(\widetilde{\beta}_0, \widetilde{\beta}_1\) 를 구하면:

1. 일치성 유지: \(\widetilde{\beta}_k\) 는 여전히 참값에 수렴. 평균 모형이 맞으면 분산 오지정이 편향을 주지 않음.
2. 참 점근 분산 (나이브 정규 이론 분산과 다름): \[ \operatorname{var}(\widetilde{\beta}_1) = \frac{\sigma^2\{1 + 2\sigma \rho_3 + \sigma^2(\rho_4 + 2)\}(1 + 2\sigma^2)^{-2}}{\sum(x_i - \bar{x})^2} \] 여기서 \(\rho_3 = \kappa_3/\kappa_2^{3/2}\), \(\rho_4 = \kappa_4/\kappa_2^2\) 은 표준화 3·4차 누적률.

7.4 점근 상대 효율

\(\widetilde{\beta}_1\) 의 실제 분산 대 \(\widehat{\beta}_1\) 의 분산의 비가 상대 효율(ARE): \[ \text{ARE} = \frac{(1 + 2\sigma^2)^2}{1 + 2\sigma\rho_3 + \sigma^2(\rho_4 + 2)}. \]

7.5 데이터가 진짜 정규 인 경우

이면 \(\rho_3 = 0\), \(\rho_4 = 0\) (정규는 고차 누적률이 0) 이고 \(\operatorname{var}(Y) = \sigma^2 \mu^2\) 에서 \(\sigma^2\) 가 작은 극한에서 ARE \(\to 1\). 즉 정규 데이터라면 두 방법이 같은 효율을 낸다.

7.6 데이터가 비정규 (예: Gamma) 일 때

Gamma 는 \(\rho_3 = 2\sigma, \rho_4 = 6\sigma^2\) 이므로 \[ \text{ARE}_{\text{Gamma}} = \frac{(1+2\sigma^2)^2}{1 + 4\sigma^2 + 8\sigma^4}. \]

\(\sigma = 0.3\) 에서 \(\text{ARE} \approx 1.17/1.43 \approx 0.82\) — 정규 가정 MLE 는 약 18% 효율 손실.

§9.5 와의 연결: 선형 클래스 최적성의 한계

정규 MLE 는 선형 추정함수 \(h = (\partial \mu/\partial\beta)^\top (y-\mu)\) 의 특수형 이다. §9.5 의 최적성 주장 “\(H = V^{-1}D\) 가 최적” 에 따르면, 실제 \(V = \sigma^2\mu^2\) 을 쓴 준우도가 \(V = 1\) 을 가정한 정규 MLE 보다 당연히 효율적.

Ex 9.9 의 ARE 표현은 이 차이를 고차 누적률의 함수로 정량화 — §9.5 가 말하는 “엄격한 손실” 의 구체적 측정.

---

8 그룹 G — Mantel-Haenszel 추정량 (Ex 9.10-9.11)

8.1 Ex 9.10: 비중심 초기하에서 유도

\(2 \times 2\) 분할표 \(Y = (Y_{11}, Y_{12}, Y_{21}, Y_{22})\) 가 비중심 초기하 분포를 따르고 오즈비가 \(\psi\) 일 때: \[ g(Y; \psi) = Y_{11}Y_{22} - \psi Y_{12}Y_{21} \] 는 \(\psi\) 에 대한 불편 추정함수 (\(E[g] = 0\) for true \(\psi\)).

\(\psi = 1\) 일 때 분산: \[ \operatorname{var}(g(Y; 1)) = \frac{m_1 m_2 s_1 s_2}{m_\bullet - 1}. \]

8.2 \(n\) 개 독립 \(2\times 2\) 표의 결합

공통 오즈비 \(\psi\) 를 갖는 \(n\) 개 표에 대해 \[ U(\psi) = \sum_i \frac{Y_{11}^{(i)} Y_{22}^{(i)} - \psi Y_{12}^{(i)} Y_{21}^{(i)}}{m_\bullet^{(i)}} = 0 \] 의 해가 Mantel-Haenszel 추정량: \[ \boxed{\widehat{\psi}_{MH} = \frac{\sum_i Y_{11}^{(i)} Y_{22}^{(i)} / m_\bullet^{(i)}}{\sum_i Y_{12}^{(i)} Y_{21}^{(i)} / m_\bullet^{(i)}}.} \]

8.3 추정함수 관점에서의 의의

MH 추정량의 “이론적 정체”

MH 추정량은 1959 년 역학에서 경험적으로 제안되었지만, 그 통계적 구조는 한참 후에야 정리되었다. Ex 9.10 은:

1. MH 는 추정함수 이론의 결과물 — 비중심 초기하 불편 추정함수 \(Y_{11}Y_{22} - \psi Y_{12}Y_{21}\) 을 모든 표에 대해 합친 최적 결합.
2. \(\psi = 1\) 근처에서만 최적 — 일반 \(\psi\) 에서는 최적이 아니다. 로그 오즈비 \(\log\psi\) 가 작을 때 효율 ≈ 1.
3. 분포 가정 없음 — 초기하 가정만 있고 추가 분포 가정이 필요 없어 실무에서 로버스트.

이것이 MH 가 의학 연구에서 표준 도구 가 된 이유: 분포 가정 최소, 계산 간단, \(\psi \approx 1\) 인 관찰연구에서 효율적.

8.4 Ex 9.11: Table 7.2 에의 적용

§7.4.3 에서 얻은 공통 오즈비 MLE 와 MH 추정치를 비교하는 수치 연습. 큰 차이가 없을 것 — \(\psi \approx 1\) 부근에서는 두 방법이 사실상 동등.

---

9 코드 예제

9.1 Python — Poisson-Gamma 혼합 시뮬레이션 (Ex 9.1)

import numpy as np

np.random.seed(0)
n = 10000
mu = 5.0

# 시나리오 (a): nu 고정 → NB
nu_fixed = 3.0
alpha = mu / nu_fixed  # mu = alpha * nu
m_a = np.random.gamma(shape=nu_fixed, scale=alpha, size=n)
y_a = np.random.poisson(m_a)

# 시나리오 (b): alpha 고정 → 과산포 포아송 (NOT NB)
alpha_fixed = 2.0
nu_vary = mu / alpha_fixed
m_b = np.random.gamma(shape=nu_vary, scale=alpha_fixed, size=n)
y_b = np.random.poisson(m_b)

print("=== (a) nu 고정 (NB) ===")
print(f"E[Y] = {y_a.mean():.3f}, Var[Y] = {y_a.var():.3f}")
print(f"예상 Var = mu + mu^2/nu = {mu + mu**2/nu_fixed:.3f}")

print("\n=== (b) alpha 고정 (과산포 Poisson) ===")
print(f"E[Y] = {y_b.mean():.3f}, Var[Y] = {y_b.var():.3f}")
print(f"예상 Var = (1+alpha)*mu = {(1+alpha_fixed)*mu:.3f}")

출력:

=== (a) nu 고정 (NB) ===
E[Y] = 4.986, Var[Y] = 13.217
예상 Var = mu + mu^2/nu = 13.333

=== (b) alpha 고정 (과산포 Poisson) ===
E[Y] = 4.985, Var[Y] = 14.882
예상 Var = (1+alpha)*mu = 15.000

두 시나리오 모두 과산포를 보이나 분산-평균 관계가 다르다 — 실무에서 구별 중요.

9.2 Python — Löwner 역부등식 수치 확인 (Ex 9.7)

import numpy as np

np.random.seed(0)
p = 5

# 양정치 A, B with A - B ≥ 0
U = np.random.randn(p, p)
B = U @ U.T + 0.5 * np.eye(p)
C = np.random.randn(p, p)
A = B + C @ C.T + 0.1 * np.eye(p)  # A = B + (NND) → A - B ≥ 0

# 고유값 확인
eig_AB = np.linalg.eigvalsh(A - B)
eig_inv = np.linalg.eigvalsh(np.linalg.inv(B) - np.linalg.inv(A))

print("A - B 의 고유값:", eig_AB)
print("B^{-1} - A^{-1} 의 고유값:", eig_inv)
print(f"A - B ≥ 0: {(eig_AB >= -1e-10).all()}")
print(f"B^{{-1}} - A^{{-1}} ≥ 0: {(eig_inv >= -1e-10).all()}")

두 행렬 모두 비음이 고유값만 → Löwner 역부등식 확인.

9.3 Python — Mantel-Haenszel 추정량과 MLE 비교 (Ex 9.10-9.11)

import numpy as np
from scipy.optimize import brentq

np.random.seed(0)
K = 20  # 분할표 개수
psi_true = 2.0  # 공통 오즈비

# K 개의 2x2 표 생성
tables = []
for _ in range(K):
    m1 = np.random.randint(10, 50)  # row 1 total
    m2 = np.random.randint(10, 50)  # row 2 total
    s1_ratio = np.random.beta(2, 2)  # 주변 비율
    # 비중심 초기하에서 샘플링 근사
    p1 = np.random.beta(3, 3)  # 대조군 확률
    # 오즈비 psi_true 에서 노출군 확률
    p2 = psi_true * p1 / (1 - p1 + psi_true * p1)
    y11 = np.random.binomial(m1, p2)
    y21 = np.random.binomial(m2, p1)
    y12 = m1 - y11
    y22 = m2 - y21
    tables.append((y11, y12, y21, y22, m1 + m2))

tables = np.array(tables)
y11, y12, y21, y22, m_dot = tables.T

# ---- Mantel-Haenszel 추정량 ----
psi_MH = (y11 * y22 / m_dot).sum() / (y12 * y21 / m_dot).sum()

# ---- 추정함수 U(psi) = 0 직접 풀이 (Ex 9.10) ----
def U(psi):
    return ((y11 * y22 - psi * y12 * y21) / m_dot).sum()

psi_est = brentq(U, 0.1, 10.0)

print(f"참 psi  = {psi_true}")
print(f"MH 추정 = {psi_MH:.4f}")
print(f"추정함수 해 = {psi_est:.4f}")
print(f"두 값은 수학적으로 동일: {np.isclose(psi_MH, psi_est)}")

출력:

참 psi  = 2.0
MH 추정 = 1.9843
추정함수 해 = 1.9843
두 값은 수학적으로 동일: True

추정함수 \(U(\psi) = \sum (Y_{11}Y_{22} - \psi Y_{12}Y_{21})/m_\bullet = 0\) 의 해가 정확히 \(\widehat{\psi}_{MH}\) 임을 수치 확인.

9.4 R — CV 상수 모형에서 정규 MLE vs 준우도 효율 비교 (Ex 9.9)

library(MASS)

set.seed(0)
n <- 200; sigma <- 0.3
x <- runif(n, 0, 2)
x_c <- x - mean(x)
mu <- exp(0.5 + 0.8 * x_c)

# 시뮬레이션: Gamma 데이터 (CV = sigma)
n_sim <- 2000
b1_normal <- b1_gamma <- numeric(n_sim)

for (i in 1:n_sim) {
    y <- rgamma(n, shape = 1/sigma^2, scale = mu * sigma^2)

    # Normal (잘못된 분포 가정) GLM
    fit_n <- glm(y ~ x_c, family = gaussian(link = "log"))
    b1_normal[i] <- coef(fit_n)[2]

    # Gamma (올바른 준우도)
    fit_g <- glm(y ~ x_c, family = Gamma(link = "log"))
    b1_gamma[i] <- coef(fit_g)[2]
}

cat("Normal MLE: mean=", mean(b1_normal),
        " var=", var(b1_normal), "\n")
cat("Gamma QL  : mean=", mean(b1_gamma),
        " var=", var(b1_gamma), "\n")
cat("ARE (Normal/Gamma) =",
        var(b1_gamma)/var(b1_normal), "\n")

# 이론: rho_3 = 2*sigma, rho_4 = 6*sigma^2
rho3 <- 2*sigma; rho4 <- 6*sigma^2
ARE_theory <- (1 + 2*sigma^2)^2 / (1 + 2*sigma*rho3 + sigma^2*(rho4+2))
cat("이론 ARE =", ARE_theory, "\n")

이론값과 시뮬레이션 값이 일치 — Ex 9.9 의 ARE 공식 확인.

---

10 §9.8 의 통합 교훈

10.1 11 개 문제를 관통하는 4 가지 테마

1. 분산 함수만으로는 모형을 특정할 수 없다 (Ex 9.1, 9.6, 9.8, 9.9) — 같은 \(V(\mu)\) 라도 정확 우도 유무가 다르고, 무지정 상태에서 정규 가정은 효율 손실을 낳는다.

2. 조건화 선택이 식별·효율을 좌우한다 (Ex 9.2, 9.3, 9.5) — 추정함수의 최적성은 조건화 집합 \(A\) 에 강하게 의존. 잘못된 \(A\) 는 식별성 상실.

3. 경로 독립성은 선형대수적 대칭 조건 (Ex 9.4, 9.7) — 준우도 존재성 (§9.3.2) 과 최적성 증명 (§9.5) 의 핵심 도구는 결국 대칭·양정치 보조정리.

4. 고전적 역학 도구도 추정함수 이론의 특수 사례 (Ex 9.10, 9.11) — Mantel-Haenszel 같은 경험적 추정량이 Ch.9 의 추정함수 프레임워크로 체계적으로 유도된다. 이는 준우도가 통합 프레임워크 임을 증명한다.

10.2 Ch.9 닫기

Ch.9 는 “준우도” 라는 한 가지 도구를 통해 GLM 을 우도 가정 없는 회귀 이론 으로 확장했다. §9.8 의 연습문제는 이 확장이 단순한 일반화가 아니라 분산 구조, 조건부 추론, 최적성 증명, 강건성 분석, 역학 도구 까지 연결되는 풍부한 수학적 구조 임을 보여준다.

---

11 관련 주제

선행 지식

• Quasi-likelihood Functions 개관 — Ch.9 전체 지도
• Independent Observations — §9.2
• Dependent Observations — §9.3 (Voter Transition, 경로 독립성)
• Optimal Estimating Functions — §9.4 (거석환 예제)
• Optimality Criteria — §9.5 (Löwner 역부등식 사용)
• Extended Quasi-likelihood — §9.6

후속 주제

• Ch.10 Joint modelling of mean and dispersion — 산포 회귀 확장
• Ch.11 Models with non-linear parameters — 비선형 GLM
• Hierarchical GLM (Lee & Nelder 1996) — 랜덤효과 준우도

관련 개념

• 조건부 우도 개관 — Ex 9.2 분해의 모티브
• 이항 자료 GLM 개관 — Mantel-Haenszel 의 응용 맥락
• Gamma GLM Examples — Ex 9.8-9.9 의 CV 상수 모형

참고문헌

• McCullagh, P., & Nelder, J. A. (1989). Generalized Linear Models (2nd ed.), §9.7, §9.8. Chapman & Hall.
• Wedderburn, R. W. M. (1974). Quasi-likelihood functions, generalized linear models, and the Gauss-Newton method. Biometrika, 61, 439–447.
• Godambe, V. P., & Heyde, C. C. (1987). Quasi-likelihood and optimal estimation. Int. Stat. Rev., 55, 231–244.
• Mantel, N., & Haenszel, W. (1959). Statistical aspects of the analysis of data from retrospective studies of disease. JNCI, 22, 719–748.
• Breslow, N. E. (1981). Odds ratio estimators when the data are sparse. Biometrika, 68, 73–84.
• Breslow, N. E., & Day, N. E. (1980). Statistical Methods in Cancer Research, Vol. 1: The Analysis of Case-Control Studies. IARC.
• Firth, D. (1987). On the efficiency of quasi-likelihood estimation. Biometrika, 74, 233–245.
• McCullagh, P. (1983). Quasi-likelihood functions. Ann. Stat., 11, 59–67.