1 도입 — 왜 독립성을 풀어야 하는가
§9.2 의 독립 관측 하에서는 \(\mathbf{V}(\boldsymbol{\mu}) = \operatorname{diag}\{V_i(\mu_i)\}\) 라는 강한 대각 구조가 모든 대수를 단순하게 만들었다. 실무의 많은 자료는 이 가정을 만족하지 않는다.
- 종단 자료(longitudinal): 같은 피험자의 반복 측정은 서로 양의 상관을 보인다. 시간 간격이 가까울수록 상관이 크다.
- 군집 샘플링(cluster): 같은 학급 학생의 성적, 같은 가구 구성원의 질병 유무는 군집 내 상관을 갖는다.
- 공간 상관(spatial): 인접 지역의 측정치는 서로 연관된다.
- 주변 도수(marginal totals): 표의 행·열 합만 관측되고 내부 셀이 숨겨진 경우, 관측 간 정보가 간접적으로 연결된다.
이 장은 \(\operatorname{cov}(\mathbf{Y}) = \sigma^2 \mathbf{V}(\boldsymbol{\mu})\) 에서 \(\mathbf{V}\) 가 대각이 아닐 수 있는 일반 상황을 다룬다. 결과는 두 가지다.
- 추정방정식은 §9.2 와 같은 형태로 유지된다. 준-스코어의 세 성질이 그대로 성립하여, \(\widehat{\beta}\) 의 점근 이론은 거의 변하지 않는다.
- 그러나 스칼라 함수로서의 준우도 \(Q(\boldsymbol{\mu}; \mathbf{y})\) 는 일반적으로 존재하지 않는다. 적분이 경로에 의존할 수 있기 때문이다. 이를 피하려면 공분산에 추가 구조적 조건을 부과해야 한다.
이 두 결과는 현대 일반화추정방정식(Generalized Estimating Equations, GEE; Liang & Zeger, 1986)의 이론적 토대이다.
2 비대각 공분산 하의 준-스코어
2.1 준-스코어 함수
공분산 \(\operatorname{cov}(\mathbf{Y}) = \sigma^2 \mathbf{V}(\boldsymbol{\mu})\) 이 \(n \times n\) 대칭 양정치 행렬이라 하자. 독립 관측과 동일한 형태의
\[ U(\beta) = \frac{1}{\sigma^2}\mathbf{D}^T \mathbf{V}^{-1}(\mathbf{Y} - \boldsymbol{\mu}) \tag{9.5} \]
을 준-스코어로 취한다. 여기서 \(\mathbf{D}_{ir} = \partial \mu_i / \partial \beta_r\) 이다.
2.2 세 성질의 유지
독립 관측에서 개별 성분에 대해 증명한 세 성질은, 비대각 \(\mathbf{V}\) 에서도 행렬 수준에서 동일하게 성립한다.
\[ \begin{aligned} E\{U_r(\beta)\} &= 0,\\ \operatorname{cov}\{U(\beta)\} &= \mathbf{D}^T \mathbf{V}^{-1}\mathbf{D}/\sigma^2 = \mathbf{i}_\beta,\\ - E\!\left(\frac{\partial U_r(\beta)}{\partial \beta_s}\right) &= \mathbf{D}^T \mathbf{V}^{-1}\mathbf{D}/\sigma^2. \end{aligned} \tag{9.8} \]
\(E(U) = 0\) 증명. \(E(\mathbf{Y} - \boldsymbol{\mu}) = \mathbf{0}\) 과 선형성에서 자명. \(\mathbf{V}\) 의 구조는 무관하다.
\(\operatorname{cov}(U)\) 증명. \[ \operatorname{cov}(U) = \frac{1}{\sigma^4}\mathbf{D}^T\mathbf{V}^{-1}\operatorname{cov}(\mathbf{Y})\mathbf{V}^{-1}\mathbf{D} = \frac{1}{\sigma^4}\mathbf{D}^T\mathbf{V}^{-1}\cdot\sigma^2\mathbf{V}\cdot\mathbf{V}^{-1}\mathbf{D} = \frac{\mathbf{D}^T\mathbf{V}^{-1}\mathbf{D}}{\sigma^2}. \] 중간 항 \(\mathbf{V}^{-1}\mathbf{V}\mathbf{V}^{-1} = \mathbf{V}^{-1}\) 가 핵심이다. 이는 \(\mathbf{V}\) 가 대각이든 비대각이든 무관하게 성립한다.
\(-E(\partial U/\partial \beta)\) 증명. \(\partial(\mathbf{Y}-\boldsymbol{\mu})/\partial \beta^T = -\mathbf{D}\) 이고 \(\mathbf{D}, \mathbf{V}\) 자체의 \(\beta\) 미분은 \((\mathbf{Y}-\boldsymbol{\mu})\) 와 곱해져 기댓값이 0 이 된다. 따라서 주도 항만 남아
\[ E\!\left(\frac{\partial U}{\partial \beta^T}\right) = \frac{1}{\sigma^2}\mathbf{D}^T\mathbf{V}^{-1}\cdot(-\mathbf{D}) = -\mathbf{i}_\beta. \]
부호를 뒤집으면 주장이 나온다.
준-스코어의 세 성질은 \(\mathbf{V}\) 의 대각성과 무관하게 성립한다. 1차 점근 이론(일치성, 정규성, \(\mathbf{i}_\beta^{-1}\) 공분산)은 독립 관측과 동일하게 따라 나온다.
2.3 추정과 점근 분포
추정방정식 \(U(\widehat{\beta}) = \mathbf{0}\), 즉
\[ \widehat{\mathbf{D}}^T \widehat{\mathbf{V}}^{-1}(\mathbf{Y} - \widehat{\boldsymbol{\mu}}) = \mathbf{0} \]
을 Newton-Raphson 또는 Fisher 스코어링으로 푼다. 점근 분포는
\[ \widehat{\beta} \approx N\!\left(\beta,\; \sigma^2 (\mathbf{D}^T \mathbf{V}^{-1}\mathbf{D})^{-1}\right). \]
일치성과 점근 정규성 증명은 §9.2 와 거의 같지만, 공분산 구조로 인해 고유값이 무한대로 발산하는 조건(regularity) 이 더 섬세해진다.
2.4 GEE 로서의 해석 — 실무 연결
Liang & Zeger (1986) 의 GEE 는 정확히 이 식을 쓴다. 종단 자료에서 피험자 \(k\) 의 관측 벡터 \(\mathbf{Y}_k\) 에 대해 작업 공분산(working covariance) \(\mathbf{V}_k(\boldsymbol{\mu}_k, \alpha)\) 을 가정한다(교환가능·AR(1)·비구조 등).
\[ U(\beta) = \sum_k \mathbf{D}_k^T \mathbf{V}_k^{-1}(\mathbf{Y}_k - \boldsymbol{\mu}_k). \]
놀라운 성질 — 일치성의 강건성. 작업 공분산이 틀려도 (진짜 상관 구조와 다르게 지정해도) \(\widehat{\beta}\) 는 여전히 일치 추정량이다. 이유는 \(E(\mathbf{Y} - \boldsymbol{\mu}) = \mathbf{0}\) 이 \(\mathbf{V}\) 에 무관하게 성립하여 \(E\{U(\beta_{\text{true}})\} = 0\) 이 보장되기 때문이다.
표준오차의 보정. 단, \(\operatorname{cov}(U)\) 는 틀린 \(\mathbf{V}\) 하에서 \(\mathbf{D}^T\mathbf{V}^{-1}\mathbf{D}/\sigma^2\) 이 아니므로 샌드위치 추정량
\[ \operatorname{cov}(\widehat{\beta})_{\text{robust}} = (\mathbf{D}^T\mathbf{V}^{-1}\mathbf{D})^{-1}\cdot \mathbf{D}^T\mathbf{V}^{-1}\widehat{\operatorname{cov}}(\mathbf{Y})\mathbf{V}^{-1}\mathbf{D}\cdot(\mathbf{D}^T\mathbf{V}^{-1}\mathbf{D})^{-1} \]
이 필요하다. 중앙의 \(\widehat{\operatorname{cov}}(\mathbf{Y})\) 는 잔차로부터 경험적으로 추정한다.
3 경로 독립성 문제
3.1 문제의 등장
준-스코어 \(U(\beta)\) 가 어떤 스칼라 함수의 그래디언트인가가 비대각 \(\mathbf{V}\) 에서는 자명하지 않다. 수학적으로 \(U\) 가 그래디언트 벡터이려면 교환성이 필요하다.
\[ \frac{\partial U_r(\beta)}{\partial \beta_s} = \frac{\partial U_s(\beta)}{\partial \beta_r},\qquad \forall r, s. \]
§9.2 의 대각 \(\mathbf{V}\) 에서는 이 조건이 저절로 성립하여 \(Q = \sum_i Q_i(\mu_i; y_i)\) 라는 스칼라가 존재했다. 비대각이면 일반적으로 성립하지 않는다 — 기댓값 수준에서는 대칭이지만(두 표현 모두 \(-\mathbf{i}_\beta\)), 점별로는 다를 수 있다.
3.2 선적분의 경로 의존성
이를 조금 다른 각도에서 본다. \(\mu\) -공간의 곡선 \(t(s), s \in [s_0, s_1]\) 을 따라
\[ Q(\boldsymbol{\mu}; \mathbf{y}, t(\cdot)) = \sigma^{-2}\int_{t(s)=\mathbf{y}}^{t(s)=\boldsymbol{\mu}} (\mathbf{y} - t)^T\{\mathbf{V}(t)\}^{-1} dt(s) \]
을 생각한다. 만약 이 선적분 값이 경로 선택에 의존하면 \(Q(\boldsymbol{\mu}; \mathbf{y})\) 는 함수로서 잘 정의되지 않는다. 경로 독립이어야만 \(Q\) 가 점별(pointwise)로 유일한 값을 갖는 스칼라 함수가 된다.
왜 이것이 중요한가. 스칼라 \(Q\) 가 존재하지 않으면
- 로그우도 비(LRT) 유사 통계량을 정의할 수 없다. 모형 비교가 이탈도 차이로 불가능.
- score 통계량 \(-2\log\Lambda\) 대신 Wald 또는 score-squared** 형태를 써야 한다.
- 최대 우도 해석이 불가능. “우도를 최대화한다”는 표현이 의미를 잃는다.
추정방정식 자체는 여전히 유효하지만, 우도 기반 추론 도구가 제한된다는 실무적 결과로 이어진다.
3.3 경로 독립성 조건 (9.9)
\(\mathbf{W}(\boldsymbol{\mu}) = \mathbf{V}^{-1}(\boldsymbol{\mu})\) 로 쓰자. 선적분이 경로 독립일 필요충분조건은 세 지표에 대한 편미분 대칭이다.
\[ \frac{\partial W_{ij}}{\partial \mu_k} = \frac{\partial W_{ik}}{\partial \mu_j} = \frac{\partial W_{jk}}{\partial \mu_i}\quad \forall\, i, j, k. \]
이는 다시 다음과 동치이다.
선적분 \(\int_{\mathbf{y}}^{\boldsymbol{\mu}} (\mathbf{y}-t)^T \mathbf{V}^{-1}(t)\,dt\) 가 경로 독립일 필요충분조건은, 볼록 함수 \(b^*(\boldsymbol{\mu})\) 가 존재하여
\[ \mathbf{V}^{-1}(\boldsymbol{\mu}) = \nabla^2 b^*(\boldsymbol{\mu}) \equiv \frac{\partial^2 b^*(\boldsymbol{\mu})}{\partial \boldsymbol{\mu}\,\partial \boldsymbol{\mu}^T} \]
을 만족하는 것이다. 이 경우 정준모수 \(\boldsymbol{\theta}(\boldsymbol{\mu})\) 와 누율함수 \(b(\boldsymbol{\theta})\) 가 자연히 정의되어
\[ \boldsymbol{\theta} = \nabla b^*(\boldsymbol{\mu}),\qquad \boldsymbol{\mu} = \nabla b(\boldsymbol{\theta}),\qquad \mathbf{V}(\boldsymbol{\mu}) = \nabla^2 b(\boldsymbol{\theta}). \tag{9.9} \]
직관. 경로 독립성은 포텐셜 함수의 존재와 동치이다. 물리학의 보존장(conservative field)에서처럼, 임의의 폐경로를 따른 순환적분이 0 이어야 한다. 이것이 2차 미분이 대칭 \(W_{ij,k} = W_{ik,j}\) 라는 조건으로 나타난다. 이 조건이 바로 지수족 구조의 자연 복원이다.
so what. 지수족 분포에서는 \(\mathbf{V}(\boldsymbol{\mu}) = \nabla^2 b(\boldsymbol{\theta})\) 가 항상 성립한다. 즉 지수족은 경로 독립성을 자동으로 만족한다. 역으로, 경로 독립성을 요구하면 사실상 지수족 구조를 되찾게 된다. 준우도가 완전한 스칼라 함수이려면 지수족 수준의 구조가 필요하다는 의미이다.
3.4 \(\mathbf{V}^{-1}\) 의 분해 조건 (9.10)
경로 독립성을 만족시키는 비대각 \(\mathbf{V}\) 를 구성하는 한 방법은 다음 분해이다.
\[ \mathbf{V}^{-1}(\boldsymbol{\mu}) = \sum_{j=1}^{k} \mathbf{A}_j^T \mathbf{W}_j(\mathbf{A}_j \boldsymbol{\mu})\mathbf{A}_j \tag{9.10} \]
여기서 \(\mathbf{A}_j\) 는 \(\boldsymbol{\mu}\) 에 무관한 행렬이며, \(\mathbf{W}_j(\boldsymbol{\gamma}_j)\) 는 \(\boldsymbol{\gamma}_j = \mathbf{A}_j \boldsymbol{\mu}\) 의 대각 함수이다.
특수 경우 \(k = 1\): \(\mathbf{A}_1 = \mathbf{I}\) 로 잡으면 \(\mathbf{V}^{-1} = \mathbf{W}_1(\boldsymbol{\mu})\) 가 대각 — 이것이 바로 §9.2 의 독립 관측이다.
일반 \(k\): 공분산이 여러 개의 “축”을 따라 선형 결합된 대각 구조로 표현됨을 의미한다. 예를 들어 다항(multinomial) 공분산 \(\boldsymbol{\Sigma} = \operatorname{diag}(\boldsymbol{\pi}) - \boldsymbol{\pi}\boldsymbol{\pi}^T\) 는 이 분해를 허용한다.
실용적 함의. “경로 독립 공분산을 직접 설계하는 것은 어렵다”(McCullagh & Nelder, 1989, §9.3.2). 실무에서는 분해 (9.10) 을 만족하는 기존 공분산(대표적으로 다항)을 쓰거나, 스칼라 \(Q\) 를 포기하고 추정방정식만 쓴다.
3.5 직선 경로 선적분 (9.11)
경로 독립성이 성립하면 어떤 경로를 택해도 \(Q\) 값이 같다. 계산 편의상 직선 경로 \(t(s) = \mathbf{y} + (\boldsymbol{\mu} - \mathbf{y})s, s \in [0, 1]\) 을 택하면
\[ Q(\boldsymbol{\mu}; \mathbf{y}) = -(\mathbf{y}-\boldsymbol{\mu})^T\!\left\{\sigma^{-2}\int_0^1 s\,\{\mathbf{V}(t(s))\}^{-1}ds\right\}(\mathbf{y}-\boldsymbol{\mu}). \tag{9.11} \]
유도 아이디어. \(t(s)\) 를 대입하면 피적분 벡터 \((\mathbf{y}-t(s)) = -s(\boldsymbol{\mu}-\mathbf{y}) = s(\mathbf{y}-\boldsymbol{\mu})\cdot(-1)\), 그리고 \(dt(s) = (\boldsymbol{\mu}-\mathbf{y})ds\) 이므로 이들을 곱하여 정리한다.
3.6 이차 근사 (9.12)
\(\mathbf{V}^{-1}(t)\) 가 경로에서 \(t\) 에 대해 거의 선형이면, 적분을 사다리꼴 공식으로 근사하여
\[ Q(\boldsymbol{\mu}; \mathbf{y}) \simeq -\frac{1}{3}(\mathbf{y}-\boldsymbol{\mu})^T\mathbf{V}^{-1}(\boldsymbol{\mu})(\mathbf{y}-\boldsymbol{\mu})/\sigma^2 - \frac{1}{6}(\mathbf{y}-\boldsymbol{\mu})^T\mathbf{V}^{-1}(\mathbf{y})(\mathbf{y}-\boldsymbol{\mu})/\sigma^2. \tag{9.12} \]
실용. \(\mathbf{V}^{-1}\) 의 closed-form 이 없을 때도, \(\boldsymbol{\mu}\) 와 \(\mathbf{y}\) 에서의 값만 쓰면 간단히 근사 준우도를 계산할 수 있다.
3.7 직관 — 편미분 대칭이 왜 볼록 함수의 존재를 보장하는가
식 (9.9) 의 대칭 조건 \(\partial (V^{-1})_{ij}/\partial \mu_k = \partial (V^{-1})_{ik}/\partial \mu_j\) 는 미적분에서 closed form (완전 미분 형태) 조건이다. 세 변수 \(\boldsymbol\mu = (\mu_1, \mu_2, \mu_3)\) 공간에서 \(\mathbf V^{-1}\) 을 벡터장으로 보면, 이 대칭은 “벡터장이 회전 없는 (curl-free) 장” 이라는 뜻이고, 그 결과 어떤 스칼라 함수 \(b^*(\boldsymbol\mu)\) 가 존재해 \(\mathbf V^{-1} = \nabla^2 b^*\) 를 만족한다.
기하학적으로는 \(\boldsymbol\mu\) - 공간에서 시작점과 끝점만 같으면 어떤 경로를 따라가도 선적분 값이 같다 는 뜻. 반대로 대칭이 깨지면 열린 곡선의 선적분처럼 경로에 따라 값이 달라져, 준-우도 \(Q\) 자체가 잘 정의되지 않는다. 이것이 §9.2 의 대각 \(\mathbf V\) 가 자동으로 경로 독립이었던 이유와, §9.3 에서 비대각 \(\mathbf V\) 일 때 조건 (9.9) 을 따로 검증해야 하는 이유이다.
3.8 정준모수 표현
경로 독립성 하에 지수족 구조가 복원되면 다음 친숙한 표현이 성립한다.
\[ Q(\boldsymbol{\mu}; \mathbf{y}) = \sigma^{-2}\{\mathbf{y}^T \boldsymbol{\theta} - b(\boldsymbol{\theta}) - b^*(\mathbf{y})\}. \]
확인. 정규((\(b(\theta) = \theta^2/2\))), 포아송((\(b(\theta) = e^\theta\))), 감마((\(b(\theta) = -\log(-\theta)\)))에서 이 공식은 각각의 로그우도의 \(\mu\) - 의존 부분을 준다. 즉 §9.2 의 Table 9.1 에서 얻은 결과들이 이 틀로 통합된다.
4 예제 — Voter Transition Probabilities
이 예제는 \(\operatorname{var}(Y)\) 가 \(E(Y)\) 만의 함수가 아닌, 따라서 §9.2 준우도가 직접 적용되지 않는 상황을 보여준다.
4.1 자료와 모형
두 연속 선거에서 각 선거구의 정당 \(C, L\) 득표 총수만 관측된다. 각 선거구의 전환 표는 다음과 같이 부분적으로 가려져 있다.
| 2차 선거 \(C\) | 2차 선거 \(L\) | 합계 | |
|---|---|---|---|
| 1차 \(C\) | \(X_1\) (비관측) | \(m_1 - X_1\) | \(m_1\) |
| 1차 \(L\) | \(X_2\) (비관측) | \(m_2 - X_2\) | \(m_2\) |
| 2차 합 | \(Y = X_1 + X_2\) (관측) | \(m_\bullet - Y\) | \(m_\bullet\) |
관측값은 2차 선거의 \(C\) 당 총 득표수 \(Y\) 뿐이다. 관심 모수는
- \(\pi_1 = P(\text{1차 }C\text{ 투표자가 2차에서도 }C)\)
- \(\pi_2 = P(\text{1차 }L\text{ 투표자가 2차에서는 }C\text{ 로 전환})\)
4.2 모멘트 구조
이항 모형 \(X_1 \sim B(m_1, \pi_1), X_2 \sim B(m_2, \pi_2)\), \(X_1 \perp X_2\) 가정 하에
\[ E(Y) = m_1\pi_1 + m_2\pi_2 = \mu,\qquad \operatorname{var}(Y) = m_1\pi_1(1-\pi_1) + m_2\pi_2(1-\pi_2). \]
핵심 관찰. \(\operatorname{var}(Y)\) 는 \(\pi_1, \pi_2\) 각각에 복잡하게 의존하며, \(\mu = m_1\pi_1 + m_2\pi_2\) 만으로는 결정되지 않는다. 따라서 §9.2 의 방식으로 스칼라 준우도를 정의할 수 없다.
그럼에도 식 (9.5) 의 준-스코어는 구성할 수 있다. \(n\) 개 선거구에서
\[ E(\mathbf{Y}) = \mathbf{M}\boldsymbol{\pi},\qquad \operatorname{cov}(\mathbf{Y}) = \operatorname{diag}\{m_{i1}\pi_1(1-\pi_1) + m_{i2}\pi_2(1-\pi_2)\} = \mathbf{V}(\boldsymbol{\pi}) \]
이고, \(\mathbf{M}\) 은 \(n \times 2\) 행렬 \((m_{i1}, m_{i2})\) 이다. 준-스코어는
\[ U(\boldsymbol{\pi}) = \mathbf{M}^T \mathbf{V}^{-1}(\boldsymbol{\pi})(\mathbf{Y} - \mathbf{M}\boldsymbol{\pi}). \tag{9.13} \]
4.3 경로 의존성의 구체적 확인
두 성분을 전개하면
\[ \begin{aligned} U_1(\boldsymbol{\pi}) &= \sum_i m_{i1}(y_i - m_{i1}\pi_1 - m_{i2}\pi_2)/V_i(\boldsymbol{\pi}),\\ U_2(\boldsymbol{\pi}) &= \sum_i m_{i2}(y_i - m_{i1}\pi_1 - m_{i2}\pi_2)/V_i(\boldsymbol{\pi}). \end{aligned} \]
직접 계산하면
\[ \frac{\partial U_1}{\partial \pi_2} \neq \frac{\partial U_2}{\partial \pi_1}. \]
왜 비대칭이 발생하는가. 분산 \(V_i(\boldsymbol\pi) = \pi_1(1-\pi_1) + \pi_2(1-\pi_2)\) 가 \(\pi_1, \pi_2\) 모두에 의존한다. 따라서 \(U_1\) 을 \(\pi_2\) 로 미분하면 분자의 \(\pi_2\) 항과 분모 \(V_i\) 의 \(\pi_2\) 항이 함께 떨어져 \(U_2\) 의 \(\pi_1\) 미분과 구조가 달라진다. 만약 \(V_i = \pi_1(1-\pi_1)\) 처럼 \(\pi_2\) 무관했다면 대칭이 복구된다 — 이것이 §9.2 의 함수적 독립 조건 이 깨지는 실제 얼굴이다.
결과. \(U\) 는 어떤 스칼라 \(Q(\boldsymbol{\pi})\) 의 그래디언트가 아니다. 경로 독립성이 성립하지 않으므로 스칼라 준우도는 없다. 그럼에도 추정방정식 \(U(\widehat{\boldsymbol{\pi}}) = 0\) 은 풀 수 있다.
4.4 수치 결과
McCullagh & Nelder 가 제시한 3 개 선거구 미니 자료.
| 선거구 | \(Y\) | \(m_1\) | \(m_2\) |
|---|---|---|---|
| 1 | 7 | 5 | 5 |
| 2 | 5 | 6 | 4 |
| 3 | 6 | 4 | 6 |
반복 수렴 결과
- 준우도 추정: \(\widehat{\boldsymbol{\pi}} = (0.3629, 0.8371)\)
- 적합값: \(\widehat{\boldsymbol{\mu}} = \mathbf{M}\widehat{\boldsymbol{\pi}} = (6.000, 5.526, 6.474)^T\)
- 정보 행렬: \(\mathbf{i}_{\boldsymbol{\pi}} = \begin{pmatrix} 41.41 & 39.79 \\ 39.79 & 42.54 \end{pmatrix}\)
- 표준오차: \(\operatorname{SE}(\widehat{\pi}_1) = 0.489, \operatorname{SE}(\widehat{\pi}_2) = 0.482\)
- 상관: \(\operatorname{corr}(\widehat{\pi}_1, \widehat{\pi}_2) = -0.948\)
해석. 두 추정량의 상관이 \(-0.948\) 로 매우 강하다. 이는 합 \(\widehat{\pi}_1 + \widehat{\pi}_2\) 는 꽤 잘 추정되지만, 차이 \(\widehat{\pi}_1 - \widehat{\pi}_2\) 나 비 \(\widehat{\pi}_1/\widehat{\pi}_2\) 같은 대비(contrast)는 정보가 거의 없다는 뜻이다. 실제로 \(\operatorname{SE}(\widehat{\pi}_1 + \widehat{\pi}_2) = 0.1565\).
4.5 MLE 와의 비교
실제 이항 로그우도를 최대화하면(EM 알고리즘 이용) \(\widehat{\boldsymbol{\pi}}_{\text{ML}} = (0.2, 1.0)\) — 모수 공간 경계. Fisher 정보 행렬
\[ \mathbf{I}_{\boldsymbol{\pi}} = \begin{pmatrix} 62.18 & 4.57 \\ 4.57 & 102.24 \end{pmatrix},\quad \mathbf{I}_{\boldsymbol{\pi}}^{-1} = \begin{pmatrix} 0.0161 & -0.0007 \\ -0.0007 & 0.0098 \end{pmatrix} \]
즉 MLE 표준오차는 약 \(0.127, 0.099\). 특히 오프-대각 상관이 작다 — MLE 의 두 추정량이 거의 독립.
왜 이렇게 다른가. 이항 로그우도에는 고차 모멘트 정보가 들어있다. 예를 들어 \(\operatorname{var}(Y_i) = m_{i1}\pi_1(1-\pi_1) + m_{i2}\pi_2(1-\pi_2)\) 는 오즈비 \(\psi = \pi_1(1-\pi_2)/\{\pi_2(1-\pi_1)\}\) 의 정보를 담는다. 준우도는 \(\mathbf{V}\) 를 경유할 뿐 이 정보를 효과적으로 쓰지 못한다.
4.6 스케일 불변성 비교
모든 관측을 100 배로 하면:
- 준우도: 같은 \(\widehat{\boldsymbol{\pi}}\), 표준오차는 \(1/10\) 배 — 자료의 선형 함수.
- MLE: \(\widehat{\boldsymbol{\pi}}_{\text{ML}} = (0.467, 0.733)\) 으로 변함 — 비선형.
교훈. 준-스코어 (9.5) 는 \((\mathbf{y}-\boldsymbol{\mu})\) 에 선형이다. 따라서 1단계 해는 자료의 선형 결합이 되어 스케일에 대해 homogeneous 한 방식으로 바뀐다. MLE 는 비선형 정보까지 이용하므로 규모에 비선형적으로 반응한다.
실무 교훈. \(\operatorname{var}(Y)\) 가 \(E(Y)\) 의 함수로 환원되지 않는 문제에서는 준우도의 효율 손실이 클 수 있다. 특히 대비(contrast)의 추정에서 손실이 두드러진다. 이 경우 완전한 우도 접근(EM, 직접 적합)이 더 효율적이다.
5 응용 분야 — 현대 GEE 로의 확장
§9.3 의 틀은 이후 30 년간 일반화추정방정식(GEE)으로 확장되어 종단·군집 자료의 표준 분석 도구가 되었다.
5.1 종단 자료
피험자 \(k\) 의 시점 \(t = 1, \ldots, T_k\) 반복 측정 \(\mathbf{Y}_k\). 피험자 간 독립, 피험자 내 상관.
\[ U(\beta) = \sum_{k=1}^K \mathbf{D}_k^T \mathbf{V}_k^{-1}(\mathbf{Y}_k - \boldsymbol{\mu}_k). \]
작업 상관(working correlation) 구조의 예.
| 구조 | \(\operatorname{corr}(Y_{kt}, Y_{kt'})\) | 용도 |
|---|---|---|
| Independence | 0 | baseline |
| Exchangeable | \(\rho\) | 군집 |
| AR(1) | \(\rho^{|t-t'|}\) | 시계열형 |
| Unstructured | \(\rho_{tt'}\) | T 작을 때 |
Liang-Zeger 의 중요 결과. 작업 상관이 틀려도 \(\widehat{\beta}\) 는 일치. 표준오차는 샌드위치 보정.
5.2 군집 샘플링
학교-학급-학생 같은 계층 자료. 군집 내 관측은 양의 상관. 각 군집을 독립 단위로 보고 합.
5.3 공간 자료
지리적 인접에 따른 상관. \(\mathbf{V}_{ij} = \sigma^2 \rho(d_{ij})\) 형태로 거리의 함수.
5.4 다항 응답
\(\operatorname{cov}(\mathbf{Y}) = \operatorname{diag}(\boldsymbol{\pi}) - \boldsymbol{\pi}\boldsymbol{\pi}^T\) 는 분해 (9.10) 을 만족하므로 스칼라 준우도가 존재한다. 다항 로지스틱 회귀의 준우도 버전(Dirichlet-multinomial 포함).
6 코드 예시
6.1 Python statsmodels.GEE
import numpy as np
import pandas as pd
import statsmodels.api as sm
from statsmodels.genmod.generalized_estimating_equations import GEE
from statsmodels.genmod.cov_struct import Exchangeable, Autoregressive, Independence
from statsmodels.genmod.families import Binomial, Poisson, Gaussian
# 종단 자료: id 가 피험자, time 이 시점
# y: 반응, X: 공변량
df = pd.DataFrame({
"id": np.repeat(range(100), 5),
"time": np.tile(range(5), 100),
"x": np.random.randn(500),
"y": np.random.binomial(1, 0.5, 500),
})
# 교환가능 작업 상관 하의 로지스틱 GEE
model = GEE.from_formula(
"y ~ x + time",
groups="id",
data=df,
family=Binomial(),
cov_struct=Exchangeable(),
)
fit = model.fit()
print(fit.summary())
print("robust SE:", fit.bse) # 샌드위치 표준오차
print("naive SE:", fit.standard_errors(cov_type="naive")) # V 기반
print("공분산 모수 alpha:", fit.cov_struct.dep_params)출력 해석.
Coef: 준우도 해 \(\widehat{\beta}\).robust SE: 작업 상관 오규정에 강건한 샌드위치 추정.naive SE: \(\mathbf{V}\) 가 정확하다고 가정한 표준오차.- 두 값이 크게 다르면 작업 상관 가정이 의심된다.
6.2 R geepack::geeglm
library(geepack)
fit_ex <- geeglm(y ~ x + time,
id = id,
data = df,
family = binomial(link = "logit"),
corstr = "exchangeable")
summary(fit_ex)
# 작업 상관 구조 비교
fit_ar <- update(fit_ex, corstr = "ar1")
fit_un <- update(fit_ex, corstr = "unstructured")
fit_in <- update(fit_ex, corstr = "independence")
QIC(fit_ex, fit_ar, fit_un, fit_in) # Pan 의 QIC 로 상관 구조 비교포인트.
geeglm의 표준오차는 기본이 샌드위치.- QIC(Quasi-Information Criterion) 로 작업 상관 구조를 선택한다.
- AIC/BIC 는 우도 기반이므로 GEE 에는 적용되지 않는다.
6.3 Voter Transition — 직접 구현
import numpy as np
# 3개 선거구 자료
Y = np.array([7., 5., 6.])
m1 = np.array([5., 6., 4.])
m2 = np.array([5., 4., 6.])
M = np.column_stack([m1, m2]) # (3, 2)
def variance(pi):
p1, p2 = pi
return m1 * p1 * (1 - p1) + m2 * p2 * (1 - p2) # (3,) 대각 V
def U(pi):
V = variance(pi)
return M.T @ ((Y - M @ pi) / V)
def info(pi):
V = variance(pi)
return M.T @ (M / V[:, None]) # D'V^{-1}D
# Fisher 스코어링
pi = np.array([0.5, 0.5])
for it in range(50):
g = U(pi)
H = info(pi)
step = np.linalg.solve(H, g)
pi_new = pi + step
pi_new = np.clip(pi_new, 1e-6, 1 - 1e-6)
if np.max(np.abs(pi_new - pi)) < 1e-10:
pi = pi_new
break
pi = pi_new
print("pi_hat =", pi) # ~ (0.363, 0.837)
print("info =\n", info(pi))
cov = np.linalg.inv(info(pi))
print("SE =", np.sqrt(np.diag(cov))) # ~ (0.489, 0.482)
print("corr =", cov[0, 1] / np.sqrt(cov[0, 0] * cov[1, 1])) # ~ -0.9486.4 경로 의존성 검증
두 편미분이 다른지 수치적으로 확인한다.
def U1(pi):
return U(pi)[0]
def U2(pi):
return U(pi)[1]
eps = 1e-6
pi = np.array([0.363, 0.837])
dU1_dpi2 = (U1(pi + np.array([0, eps])) - U1(pi - np.array([0, eps]))) / (2 * eps)
dU2_dpi1 = (U2(pi + np.array([eps, 0])) - U2(pi - np.array([eps, 0]))) / (2 * eps)
print(f"dU1/dpi2 = {dU1_dpi2:.4f}")
print(f"dU2/dpi1 = {dU2_dpi1:.4f}")
# 두 값이 달라서 경로 의존이 확인된다7 요약 — §9.2 와 §9.3 의 비교
| 항목 | §9.2 독립 | §9.3 종속 |
|---|---|---|
| 공분산 | 대각, \(V_i(\mu_i)\) | 비대각 일반 |
| 준-스코어의 3 성질 | 성립 | 성립 (식 9.8) |
| 스칼라 \(Q(\boldsymbol{\mu}; \mathbf{y})\) | 항상 존재 | 경로 독립 시만 |
| 경로 독립성 조건 | 자동 성립 | \(\mathbf{V}^{-1} = \nabla^2 b^*\) (식 9.9) |
| 추정방정식 | \(\mathbf{D}^T\mathbf{V}^{-1}(\mathbf{y}-\boldsymbol{\mu}) = \mathbf{0}\) | 동일 |
| 일치성 | 2차 모멘트 정확 시 | \(E(\mathbf{Y})\) 맞으면 \(\mathbf{V}\) 오규정에도 성립 |
| 점근 공분산 | \(\sigma^2(\mathbf{D}^T\mathbf{V}^{-1}\mathbf{D})^{-1}\) | 동일(정확 \(\mathbf{V}\)) / 샌드위치(오규정) |
| LRT 유사 검정 | 가능(이탈도 차이) | 스칼라 \(Q\) 없으면 score/Wald 만 |
| 대표 응용 | 과산포 GLM | GEE, 종단·군집 자료 |
한 문장 결론. §9.3 의 이론적 핵심은 “추정은 여전히 가능하지만, 우도 기반 추론 도구(LRT)는 경로 독립성 없이는 정의되지 않는다”는 것이다. 이 점이 현대 GEE 가 Wald 와 score 검정을 선호하는 이유의 밑에 있다.