Oneway ANOVA: Inferences Regarding Linear Combinations of Means

대비를 통한 추론 — t 통계량 유도, 계획 비교, 직교 대비 (Casella §11.2.3)

일원분류 ANOVA에서 처리 평균들의 선형결합(linear combination)에 대한 추론 체계를 상세히 다룬다. 합동 분산 추정량의 분포, t 통계량의 유도, 개별 대비 검정과 신뢰구간, 계획 비교(planned comparisons), 직교 대비(orthogonal contrasts)까지 Casella & Berger §11.2.3을 중심으로 설명한다.

Statistics
저자

Kwangmin Kim

공개

2026년 04월 05일

1 개요

이전 포스트(Oneway ANOVA: 고전적 가설과 F 검정)에서 ANOVA 귀무가설이 모든 대비(contrast)가 0이라는 조건과 동치임을 보였고, 이로부터 F 통계량을 유도했다.

F 통계량은 “어떤 집단들 사이에 차이가 있는가?”라는 전체적 질문에만 답한다. 실무에서 더 유용한 것은 “어떤 집단이 얼마나 다른가?”라는 구체적 질문이며, 이에 대한 답은 선형결합과 대비에 대한 추론을 통해 얻는다.

Casella & Berger(2002, §11.2.3)는 이 추론 체계를 단계적으로 구성한다: 집단 평균의 분포 → 선형결합의 분포 → 합동 분산 추정 → t 통계량 → 검정과 신뢰구간.


2 기초: 집단 평균의 분포

2.1 개별 집단 평균

ANOVA 가정 하에서 \(Y_{ij} \sim N(\theta_i, \sigma^2)\) 이므로, \(j = 1, \ldots, n_i\) 의 평균은:

\[ \bar{Y}_{i\cdot} = \frac{1}{n_i} \sum_{j=1}^{n_i} Y_{ij} \sim N\!\left(\theta_i,\, \frac{\sigma^2}{n_i}\right), \quad i = 1, \ldots, k \]

이는 정규분포의 선형변환 성질에서 바로 나온다. 집단 \(i\) 의 표본 크기가 클수록 \(\bar{Y}_{i\cdot}\) 의 분산이 작아지고, 추정이 정밀해진다.

또한 서로 다른 집단의 평균들은 독립이다:

\[ \bar{Y}_{1\cdot}, \bar{Y}_{2\cdot}, \ldots, \bar{Y}_{k\cdot} \text{ 는 독립} \]

이는 ANOVA 가정의 비상관 조건(오차 간 독립)에서 따른다.

2.2 선형결합의 분포

임의의 상수 벡터 \(a = (a_1, \ldots, a_k)\) 에 대해, 선형결합 \(\sum_{i=1}^k a_i \bar{Y}_{i\cdot}\) 은 정규분포를 따른다 (독립 정규분포들의 선형결합은 정규분포):

\[ \sum_{i=1}^{k} a_i \bar{Y}_{i\cdot} \sim N\!\left(\sum_{i=1}^{k} a_i \theta_i,\; \sigma^2 \sum_{i=1}^{k} \frac{a_i^2}{n_i}\right) \]

기댓값: \(E\!\left(\sum_i a_i \bar{Y}_{i\cdot}\right) = \sum_i a_i E(\bar{Y}_{i\cdot}) = \sum_i a_i \theta_i\)

분산: \(\text{Var}\!\left(\sum_i a_i \bar{Y}_{i\cdot}\right) = \sum_i a_i^2 \text{Var}(\bar{Y}_{i\cdot}) = \sigma^2 \sum_i \frac{a_i^2}{n_i}\) (독립이므로 공분산 항 없음)

따라서 표준화하면:

\[ \frac{\sum_{i=1}^{k} a_i \bar{Y}_{i\cdot} - \sum_{i=1}^{k} a_i \theta_i}{\sigma \sqrt{\sum_{i=1}^{k} a_i^2 / n_i}} \sim N(0, 1) \]

이 결과는 정확하지만, 실용적이지 않다. \(\sigma\) 를 알지 못하기 때문이다.


3 합동 분산 추정량 \(S_p^2\)

3.1 개별 집단 분산 추정

집단 \(i\) 의 표본분산은:

\[ S_i^2 = \frac{1}{n_i - 1} \sum_{j=1}^{n_i} (Y_{ij} - \bar{Y}_{i\cdot})^2, \quad i = 1, \ldots, k \]

\(S_i^2\)\(\sigma^2\) 의 불편추정량이며, 정규성 가정에 의해:

\[ \frac{(n_i - 1) S_i^2}{\sigma^2} \sim \chi^2_{n_i - 1} \]

또한 \(\bar{Y}_{i\cdot}\)\(S_i^2\) 는 독립이다 (정규분포의 성질).

3.2 합동 추정량의 구성

ANOVA의 등분산 가정 하에서 모든 집단이 동일한 \(\sigma^2\) 를 가지므로, \(k\) 개의 추정량을 합산하여 하나의 정밀한 추정량을 구성한다 (Casella & Berger, 2002, eq. 11.2.5):

\[ S_p^2 = \frac{1}{N-k} \sum_{i=1}^{k} (n_i - 1) S_i^2 = \frac{1}{N-k} \sum_{i=1}^{k} \sum_{j=1}^{n_i} (Y_{ij} - \bar{Y}_{i\cdot})^2 \]

여기서 \(N = \sum_{i=1}^k n_i\) 이고, \(N - k = \sum_{i=1}^k (n_i - 1)\) 이다.

왜 이 가중 평균을 사용하는가? 자유도로 가중하면 작은 집단보다 큰 집단이 더 많이 반영되어 전체적으로 더 효율적인 추정이 된다. 집단 크기가 모두 같으면 단순 평균과 같다.

3.3 \(S_p^2\) 의 분포

\(S_i^2\) 들이 독립이므로 Lemma 5.3.2(독립 카이제곱의 합은 카이제곱):

\[ \frac{(N-k) S_p^2}{\sigma^2} = \sum_{i=1}^{k} \frac{(n_i-1)S_i^2}{\sigma^2} \sim \chi^2_{N-k} \]

자유도는 \(N - k\) 이다. 직관적으로, 각 집단에서 평균 1개를 추정하므로 \(k\) 개의 자유도를 잃는다.

3.4 \(S_p^2\) 와 집단 평균의 독립성

\(S_p^2\)\(\bar{Y}_{1\cdot}, \ldots, \bar{Y}_{k\cdot}\) 각각과 독립이다. 이는 정규분포에서 표본 평균과 표본분산의 독립성(크레이머-라오 정리 관련)에서 따른다.

이 독립성이 있어야만 t 비율을 구성할 수 있다.


4 t 통계량 유도 (Casella eq. 11.2.6)

\(\sigma\)\(S_p\) 로 대체하면 t 분포가 나타난다:

\[ \frac{\sum_{i=1}^{k} a_i \bar{Y}_{i\cdot} - \sum_{i=1}^{k} a_i \theta_i}{\sqrt{S_p^2 \sum_{i=1}^{k} a_i^2 / n_i}} \sim t_{N-k} \]

유도: 표준정규 \(Z = N(0,1)\) 와 독립 카이제곱 \(V = \chi^2_{N-k}\) 가 있을 때, \(Z / \sqrt{V/(N-k)} \sim t_{N-k}\) 임을 이용한다.

여기서: - 분자의 표준화된 선형결합 \(\sim N(0,1)\) - \(S_p^2 = \sigma^2 V / (N-k)\) 이고 \(V \sim \chi^2_{N-k}\) - 분자와 \(S_p^2\) 는 독립

따라서 \(S_p\) 로 나누면 \(t_{N-k}\) 가 나온다.

왜 이 형태인가? 추정량(분자)을 추정 표준오차(분모)로 나누면 표준화된다. 분자는 \(H_0\) 하에서 \(N(0,1)\) 이고, 분모는 독립인 \(\chi^2\) 추정에 기반한다. 정규를 독립 카이제곱 제곱근으로 나눈 것이 정확히 \(t\)-분포의 정의이다. \(\sigma\) 를 모르기 때문에 \(S_p\) 로 대체하는 것이 정규에서 \(t\) 로 바뀌는 이유이다.

왜 자유도가 \(N-k\) 인가?

\(t_{N-k}\) 의 자유도는 분모의 카이제곱 분포의 자유도이다. \(S_p^2\) 를 계산할 때 \(k\) 개의 집단 평균을 각각 추정하므로 \(k\) 개의 자유도를 소모한다. 따라서 \(N\) 개의 관측값에서 \(k\) 개의 자유도를 뺀 \(N - k\) 가 된다.

자유도가 클수록 t 분포가 표준정규에 가까워지므로, 집단 수가 적고 집단 내 표본이 많을수록 검정력이 높다.

이 t 통계량은 선형결합 \(\sum a_i \theta_i\) 에 대한 pivotal quantity로, 검정과 신뢰구간을 동시에 구성하는 기반이 된다.


5 가설검정 (Casella eq. 11.2.7)

5.1 검정 절차

\(H_0 : \sum_{i=1}^k a_i \theta_i = c\)\(H_1 : \sum_{i=1}^k a_i \theta_i \neq c\) 를 검정할 때, (\(c = 0\) 이 가장 흔한 경우):

\[ T = \frac{\sum_{i=1}^{k} a_i \bar{Y}_{i\cdot} - c}{\sqrt{S_p^2 \sum_{i=1}^{k} a_i^2 / n_i}} \]

수준 \(\alpha\) 검정의 기각역:

\[ |T| > t_{N-k,\, \alpha/2} \]

단측 검정의 경우: - \(H_1 : \sum a_i \theta_i > c\) 이면 \(T > t_{N-k,\, \alpha}\) 일 때 기각 - \(H_1 : \sum a_i \theta_i < c\) 이면 \(T < -t_{N-k,\, \alpha}\) 일 때 기각

5.2 p-값 계산

\[ p = 2 \cdot P(t_{N-k} > |T_{\text{obs}}|) \quad \text{(양측검정)} \]

5.3 검정의 오류 관리

오류 유형 정의 제어 방법
제1종 오류 \(H_0\) 참인데 기각 \(\alpha\) 로 직접 제어
제2종 오류 \(H_1\) 참인데 기각 못 함 표본 크기, 대비 선택으로 검정력 확보
비교군이 많을 때 FWER 여러 검정의 누적 오류율 상승 Bonferroni, Scheffé 등으로 제어

6 신뢰구간 (Casella eq. 11.2.8)

\(\sum_{i=1}^k a_i \theta_i\) 에 대한 \(1-\alpha\) 신뢰구간은:

\[ \sum_{i=1}^{k} a_i \bar{Y}_{i\cdot} \;\pm\; t_{N-k,\, \alpha/2} \sqrt{S_p^2 \sum_{i=1}^{k} \frac{a_i^2}{n_i}} \]

이는 t 통계량을 pivot으로 사용해 역산한 것이다.

6.1 신뢰구간의 폭 결정 요인

\[ \text{폭} = 2 t_{N-k,\, \alpha/2} \underbrace{\sqrt{S_p^2}}_{\text{오차 크기}} \underbrace{\sqrt{\sum_i a_i^2 / n_i}}_{\text{대비 구조}} \]

요인 폭을 좁히려면
\(t_{N-k,\alpha/2}\) 신뢰수준을 낮추거나, \(N-k\) 증가 (표본 늘리기)
\(S_p^2\) 실험 오차를 줄임 (블로킹, 무작위화 등)
\(\sum a_i^2/n_i\) 비교 집단의 표본 크기를 늘리고, 불필요한 집단 제외
두 집단 비교 시 t 검정과의 차이

\(a = (1, -1, 0, \ldots, 0)\) 으로 \(\theta_1 - \theta_2\) 를 추론하면, ANOVA 기반 검정의 기각역은:

\[ \left| \frac{\bar{Y}_{1\cdot} - \bar{Y}_{2\cdot}}{\sqrt{S_p^2(1/n_1 + 1/n_2)}} \right| > t_{N-k,\, \alpha/2} \]

이는 2-표본 t 검정과 형태가 동일하지만, 두 가지 차이가 있다:

  1. \(S_p^2\) 의 범위: 2-표본 t 검정은 두 집단의 정보만 사용하지만, ANOVA 기반 검정은 모든 \(k\) 개 집단의 정보를 사용하여 \(S_p^2\) 를 추정한다.
  2. 자유도: 2-표본 t 검정의 자유도는 \(n_1 + n_2 - 2\) 이지만, ANOVA 기반 검정의 자유도는 \(N - k \geq n_1 + n_2 - 2\) 이다.

결과적으로 ANOVA 기반 검정이 자유도가 더 크고 검정력이 높다. 단, 등분산 가정이 충족될 때만 이 이점이 유효하다.


7 계획 비교 (Planned Comparisons)

7.1 선험적 대비 선택의 중요성

대비를 데이터를 보기 전에 결정하면 다음 이점이 있다:

  1. 과학적 질문을 명확히 진술하게 된다
  2. 다중 비교 문제가 발생하지 않거나, 비교 수가 제한되어 Bonferroni 보정이 효율적이다
  3. 검정력이 사후(post-hoc) 방법보다 높다

7.2 자주 사용하는 대비 패턴

비교 목적 대비 벡터 \(a\) 대비 \(\sum a_i \theta_i\)
집단 1 vs. 집단 2 \((1, -1, 0, 0)\) \(\theta_1 - \theta_2\)
처리 vs. 대조군 \((1/3, 1/3, 1/3, -1)\) \(\bar{\theta}_{\text{trt}} - \theta_{\text{ctrl}}\)
처리군 내 비교 \((1, -1/2, -1/2, 0)\) \(\theta_1 - (\theta_2+\theta_3)/2\)
선형 추세 (등간격) \((-3, -1, 1, 3)\) 선형 trend component
이차 추세 (등간격) \((1, -1, -1, 1)\) 이차 trend component

7.3 선형·이차 추세 대비 (등간격 설계)

\(k\) 집단이 등간격으로 배치된 실험(예: 용량 반응 실험)에서 추세 대비를 사용하면 단조적 효과나 역U형 효과를 감지할 수 있다.

용량-반응 실험을 생각해보자. 4개 집단이 저용량, 중용량, 고용량, 초고용량이라면, 단순히 ‘모든 쌍’을 비교하는 것보다 ’용량이 증가하면 반응도 증가하는가?’(선형 추세)를 먼저 묻는 것이 더 날카롭다. 직교 추세 대비는 이 질문을 수학적으로 구현한다: 선형 대비는 ‘전반적 기울기’, 이차 대비는 ‘기울기의 변화(가속/감속)’, 삼차 대비는 ’변곡점’을 검출한다.

\(k = 4\) 인 경우 표준화된 직교 대비 계수:

차수 \(a_1\) \(a_2\) \(a_3\) \(a_4\) \(\sum a_i^2\)
선형 \(-3\) \(-1\) \(1\) \(3\) \(20\)
이차 \(1\) \(-1\) \(-1\) \(1\) \(4\)
삼차 \(-1\) \(3\) \(-3\) \(1\) \(20\)

8 직교 대비 (Orthogonal Contrasts)

8.1 직교성의 정의

두 대비 \(a^{(l)} = (a_1^{(l)}, \ldots, a_k^{(l)})\)\(a^{(l')} = (a_1^{(l')}, \ldots, a_k^{(l')})\)직교(orthogonal)하다는 것은:

\[ \sum_{i=1}^{k} \frac{a_i^{(l)} a_i^{(l')}}{n_i} = 0 \quad \text{for } l \neq l' \]

균형 설계(\(n_1 = \cdots = n_k = n\))에서는 \(\sum_i a_i^{(l)} a_i^{(l')} = 0\) 으로 단순화된다.

기하학적 직관: 직교 대비는 수직인 벡터와 같다. 두 대비가 직교하면 하나의 대비가 설명하는 변동과 다른 대비가 설명하는 변동이 전혀 겹치지 않는다. 따라서 SSB 를 직교 대비들로 분해하면 각 성분의 제곱합이 겹침 없이 정확히 더해진다: \(\text{SSB} = \text{SS}_{\text{대비1}} + \text{SS}_{\text{대비2}} + \cdots\) . 이것은 3차원 벡터를 \(x, y, z\) 성분으로 분해하면 \(\|\mathbf{v}\|^2 = v_x^2 + v_y^2 + v_z^2\) 인 것과 같은 원리이다.

8.2 직교 대비의 성질

  1. 독립성: 정규성 가정 하에서 직교 대비들은 독립이다 (Lemma 5.3.3)
  2. SSB의 분해: \(k-1\) 개의 직교 대비로 SSB를 \((k-1)\) 개의 독립적인 제곱합으로 분해할 수 있다
  3. F 검정과의 관계: \(F\) 통계량은 \(k-1\) 개 직교 대비 검정의 평균과 같다

8.3 직교 대비의 SSB 분해

각 대비의 대비 제곱합(contrast sum of squares):

\[ SS_{a^{(l)}} = \frac{\left(\sum_{i=1}^k a_i^{(l)} \bar{Y}_{i\cdot}\right)^2}{\sum_{i=1}^k (a_i^{(l)})^2 / n_i} \]

\(k-1\) 개의 직교 대비로 SSB가 분해된다:

\[ \mathrm{SSB} = SS_{a^{(1)}} + SS_{a^{(2)}} + \cdots + SS_{a^{(k-1)}} \]

\(SS_{a^{(l)}}\) 는 자유도 1인 카이제곱 분포를 따르고 서로 독립이다. 따라서 각 대비에 대해 개별 F 검정(자유도 \(1, N-k\))을 수행할 수 있다.

직교 대비의 실용적 의미

\(k-1\) 개의 직교 대비를 미리 설계하면: - SSB를 \(k-1\) 개의 독립적인 변동 성분으로 완전 분해 - 각 성분이 특정 과학적 질문에 대응 - 개별 검정의 제1종 오류가 독립적으로 \(\alpha\) 유지 (단, 전체 FWER은 Bonferroni 필요)

예: \(k=4\) 집단이 처리 3개 + 대조군 1개이면, 대조군 vs. 처리 평균, 처리 내 선형 추세, 처리 내 이차 추세라는 3개의 직교 대비로 처리 효과를 완전히 분해할 수 있다.


9 예시: 어류 독소 실험의 구체적 대비

Casella & Berger(2002, Example 11.2.6 연장)의 어류 독소 데이터:

\[ \text{T1: } (28, 23, 14, 27), \quad \text{T2: } (33, 36, 34, 29, 31, 34) \] \[ \text{T3: } (18, 21, 20, 22, 24), \quad \text{Ctrl: } (11, 14, 11, 16) \]

\(k=4\), \(N=19\), \(S_p^2 = 12.72\) (자유도 15).

9.1 대비 1: T2 vs. T1 (가장 높은 vs. 가장 낮은 독소)

\(a = (-1, 1, 0, 0)\):

\[ \sum a_i \bar{y}_{i\cdot} = -23.0 + 32.83 = 9.83 \]

\[ SE = \sqrt{12.72 \cdot (1/4 + 1/6)} = \sqrt{12.72 \times 0.4167} = 2.302 \]

\[ t = 9.83 / 2.302 = 4.27, \quad p < 0.001 \]

T2와 T1의 평균 생존 시간은 유의하게 다르다.

9.2 대비 2: 독소 평균 vs. 대조군

\(a = (1/3, 1/3, 1/3, -1)\):

\[ \sum a_i \bar{y}_{i\cdot} = \frac{23.0 + 32.83 + 21.0}{3} - 13.0 = 25.61 - 13.0 = 12.61 \]

\[ SE = \sqrt{12.72 \cdot (1/(9 \cdot 4) + 1/(9 \cdot 6) + 1/(9 \cdot 5) + 1/4)} = \sqrt{12.72 \cdot 0.3569} = 2.130 \]

\[ t = 12.61 / 2.130 = 5.92, \quad p < 0.0001 \]

독소들이 대조군에 비해 생존 시간을 유의하게 연장시킨다는 강한 증거가 있다.


10 코드 예시

10.1 Step 1: 순수 Python — 선형결합 추론 전체 구현

import numpy as np
from scipy import stats

# 데이터 (Casella & Berger, Example 11.2.1)
groups = {
    "T1":   np.array([28, 23, 14, 27]),
    "T2":   np.array([33, 36, 34, 29, 31, 34]),
    "T3":   np.array([18, 21, 20, 22, 24]),
    "Ctrl": np.array([11, 14, 11, 16]),
}
k    = len(groups)
ns   = np.array([len(v) for v in groups.values()])
N    = ns.sum()
means = np.array([v.mean() for v in groups.values()])

# 합동 분산 S_p^2
SSW = sum(((v - m) ** 2).sum() for v, m in zip(groups.values(), means))
Sp2 = SSW / (N - k)
print(f"S_p^2 = {Sp2:.4f}  (자유도 {N-k})")
# => S_p^2 = 12.7222  (자유도 15)

# ------- 선형결합 추론 함수 -------
def linear_combo_inference(a, c=0.0, alpha=0.05):
    """
    선형결합 sum(a_i * theta_i) = c 에 대한 검정과 신뢰구간.
    a: 계수 벡터, c: 귀무가설 값 (default=0)
    """
    a = np.array(a, dtype=float)
    # 추정값
    est = (a * means).sum()
    # 표준오차 (eq. 11.2.6 분모)
    se  = np.sqrt(Sp2 * (a**2 / ns).sum())
    # t 통계량
    t_stat = (est - c) / se
    # 임계값 및 p-값 (양측)
    t_crit = stats.t.ppf(1 - alpha/2, df=N - k)
    p_val  = 2 * stats.t.sf(abs(t_stat), df=N - k)
    # 신뢰구간 (eq. 11.2.8)
    ci_lo  = est - t_crit * se
    ci_hi  = est + t_crit * se
    return {
        "추정값": est, "SE": se, "t": t_stat, "자유도": N-k,
        "p-value": p_val, "95% CI": (ci_lo, ci_hi),
        "기각": abs(t_stat) > t_crit,
    }

# 대비 1: T2 vs T1
print("\n[대비] T2 - T1")
r1 = linear_combo_inference([-1, 1, 0, 0])
for k_, v in r1.items():
    print(f"  {k_}: {v:.4f}" if isinstance(v, float) else f"  {k_}: {v}")

# 대비 2: 독소 평균 vs 대조군
print("\n[대비] 독소 평균 - Ctrl")
r2 = linear_combo_inference([1/3, 1/3, 1/3, -1])
for k_, v in r2.items():
    print(f"  {k_}: {v:.4f}" if isinstance(v, float) else f"  {k_}: {v}")

# 대비 3: T3 vs Ctrl (단측)
print("\n[대비 단측] T3 - Ctrl > 0")
r3 = linear_combo_inference([0, 0, 1, -1])
t_obs = r3["t"]
p_one = stats.t.sf(t_obs, df=N - k)    # P(T > t_obs)
print(f"  t = {t_obs:.4f},  p(단측) = {p_one:.4f}")

10.2 Step 2: 직교 대비 분해 — SSB 분해 확인

import numpy as np

# k=4 등간격 설계를 위한 표준 직교 대비 계수
# (균형 설계 n_1=...=n_4=n 가정하여 단순화)
ortho_contrasts = {
    "선형":  np.array([-3, -1,  1,  3]),
    "이차":  np.array([ 1, -1, -1,  1]),
    "삼차":  np.array([-1,  3, -3,  1]),
}

# 직교성 확인: sum(a^(l) * a^(l')) = 0 for l != l'
ks = list(ortho_contrasts.keys())
for i in range(len(ks)):
    for j in range(i+1, len(ks)):
        a1 = ortho_contrasts[ks[i]]
        a2 = ortho_contrasts[ks[j]]
        inner = np.dot(a1, a2)
        print(f"직교성 <{ks[i]}, {ks[j]}> = {inner}  ({'직교 O' if inner==0 else '직교 X'})")

# 각 직교 대비의 대비 제곱합
# SS_a = (sum a_i * ybar_i)^2 / (sum a_i^2 / n_i)
def contrast_ss(a, means, ns):
    """대비 제곱합 (contrast sum of squares)."""
    a = np.array(a, dtype=float)
    num = (a * means).sum() ** 2
    den = (a**2 / ns).sum()
    return num / den

print("\n--- 직교 대비 SSB 분해 ---")
total_css = 0
for name, a in ortho_contrasts.items():
    css = contrast_ss(a, means, ns)
    total_css += css
    print(f"SS_{name} = {css:.4f}  (df=1)")

print(f"\n합계 = {total_css:.4f}")
print(f"SSB  = {SSW:.4f}")   # 비교용 (실제로는 groups로 계산된 SSB)
# sum of orthogonal contrast SS == SSB (균형 설계에서 정확히 성립)

10.3 Step 3: statsmodels로 계획 비교 수행

import pandas as pd
import statsmodels.formula.api as smf
from scipy import stats

# 데이터프레임
rows = []
for name, obs in groups.items():
    for val in obs:
        rows.append({"group": name, "y": float(val)})
df = pd.DataFrame(rows)

# 모형 적합
model = smf.ols("y ~ C(group, Treatment('Ctrl'))", data=df).fit()
print(model.summary())

# 사용자 정의 대비 (T2 - T1)
from statsmodels.sandbox.stats.multicomp import multipletests
contrasts_list = [
    ([0, 0, 0, 0, 1, -1], "T3 vs T2"),  # 계수는 모형 디자인에 따라 다름
]

# 직접 대비 계산 (위에서 정의한 linear_combo_inference 사용)
contrasts_to_test = {
    "T2 vs T1":           [-1,  1,  0,  0],
    "T3 vs T1":           [-1,  0,  1,  0],
    "T3 vs T2":           [ 0, -1,  1,  0],
    "독소평균 vs Ctrl":    [ 1/3, 1/3, 1/3, -1],
}

alpha_bonf = 0.05 / len(contrasts_to_test)  # Bonferroni 보정
print(f"\nBonferroni 보정 유의수준: {alpha_bonf:.4f}")
print("-" * 60)
for label, a in contrasts_to_test.items():
    r = linear_combo_inference(a)
    sig = "***" if r["p-value"] < alpha_bonf else (
          "."  if r["p-value"] < 0.05 else "ns")
    print(f"{label:25s}  추정={r['추정값']:7.3f}  "
          f"t={r['t']:6.3f}  p={r['p-value']:.4f}  {sig}")

11 검정력 분석 (Power Analysis)

특정 대비 \(\Delta = \sum a_i \theta_i\) 가 존재할 때 검정의 검정력은 비심 t 분포로 계산한다.

비중심 모수 \(\lambda\)신호 대 잡음 비(signal-to-noise ratio)이다. 분자는 참된 대비 크기 \(\Delta\) 이고, 분모는 추정의 불확실성 \(\sigma\sqrt{\sum a_i^2/n_i}\) 이다. \(\lambda\) 가 클수록 \(H_0\) 에서 멀리 떨어져 있어 기각이 쉽다. 검정력은 본질적으로 \(\lambda\) 의 증가함수이다.

비심도 모수(noncentrality parameter):

\[ \lambda = \frac{\Delta}{\sigma \sqrt{\sum_{i=1}^k a_i^2 / n_i}} \]

검정력:

\[ \text{Power} = P\!\left(|T| > t_{N-k,\,\alpha/2} \;\Big|\; \lambda\right) = P\!\left(|t_{N-k}(\lambda)| > t_{N-k,\,\alpha/2}\right) \]

from scipy import stats

def contrast_power(a, delta, sigma, ns, N_total, k, alpha=0.05):
    """
    delta: 검출하고자 하는 대비의 실제 크기
    sigma: 공통 표준편차 (sigma, not sigma^2)
    """
    a = np.array(a, dtype=float)
    df = N_total - k
    se_factor = np.sqrt((a**2 / ns).sum())
    lam = delta / (sigma * se_factor)   # 비심도 모수
    t_crit = stats.t.ppf(1 - alpha/2, df=df)
    # 비심 t 분포로 검정력 계산
    nc_dist = stats.nct(df=df, nc=lam)
    power = 1 - nc_dist.cdf(t_crit) + nc_dist.cdf(-t_crit)
    return {"lambda": lam, "power": power}

# 어류 독소 데이터에서 T2-T1 대비 검정력
# delta=9.83 (관측된 차이), sigma=sqrt(12.72)=3.57
a_ex = np.array([-1, 1, 0, 0])
result = contrast_power(
    a=a_ex, delta=9.83, sigma=np.sqrt(12.72),
    ns=ns, N_total=N, k=4
)
print(f"비심도 lambda = {result['lambda']:.4f}")
print(f"검정력 = {result['power']:.4f}")

# 표본 크기 증가에 따른 검정력 변화
for n_per_group in [4, 6, 10, 15, 20]:
    ns_eq = np.array([n_per_group] * 4)
    N_eq  = ns_eq.sum()
    r     = contrast_power(a_ex, delta=9.83, sigma=3.57, ns=ns_eq, N_total=N_eq, k=4)
    print(f"n={n_per_group:2d} (N={N_eq:2d}): lambda={r['lambda']:.3f}, power={r['power']:.4f}")

12 핵심 정리

ANOVA에서 선형결합 추론의 흐름:

ANOVA 가정: Y_ij ~ N(θ_i, σ²)
       ↓
집단 평균: Ȳ_i. ~ N(θ_i, σ²/n_i), 독립
       ↓
선형결합: Σa_i Ȳ_i. ~ N(Σa_i θ_i, σ² Σa_i²/n_i)
       ↓
σ 대체: S_p² = SSW/(N-k) ~ σ² χ²_{N-k}/(N-k), Ȳ와 독립
       ↓
t 통계량: (Σa_i Ȳ_i. - Σa_i θ_i) / √(S_p² Σa_i²/n_i) ~ t_{N-k}
       ↓
검정 (eq. 11.2.7): |T| > t_{N-k, α/2}
신뢰구간 (eq. 11.2.8): Σa_i Ȳ_i. ± t_{N-k, α/2} · SE
개념 핵심 수식 Casella 참조
선형결합 분포 \(\sum a_i \bar{Y}_{i\cdot} \sim N(\sum a_i \theta_i, \sigma^2 \sum a_i^2/n_i)\) §11.2.3
합동 분산 \(S_p^2 = \mathrm{SSW}/(N-k)\), \((N-k)S_p^2/\sigma^2 \sim \chi^2_{N-k}\) eq. 11.2.5
t 통계량 \(T \sim t_{N-k}\) eq. 11.2.6
검정 기각역 \(|T| > t_{N-k, \alpha/2}\) eq. 11.2.7
신뢰구간 \(\hat{\Delta} \pm t_{N-k,\alpha/2} \cdot SE\) eq. 11.2.8
직교 대비 \(\sum a_i^{(l)} a_i^{(l')}/n_i = 0\) eq. 11.2.18

13 관련 주제

선행 지식

후속 주제

관련 개념

Subscribe

Enjoy this blog? Get notified of new posts by email: