Fligner-Killeen 검정 (scipy)¶

(1) \(g_3\) 는 \(g_1\) 을 \(20\) 만큼 옮긴 것이다. 중앙값도 \(12.5\) 에서 \(32.5\) 로 함께 옮겨 가므로

로 편차가 글자 그대로 같다. 따라서 \(\bar a_3 = \bar a_1\) 이고, 세 집단 문제가 두 수 \(\bar a_1\) 과 \(\bar a_2\) 로 줄어든다.

편차값이 네 개뿐인 까닭. 세 집단 모두 연속한 정수 근처에 모여 있고 \(n_i = 8\) 이 짝수라 중앙값이 두 가운데 값의 중점, 곧 \(\ast.5\) 꼴이 된다. 그래서 모든 편차가 \(0.5,\ 1.5,\ 2.5,\ 3.5\) 네 값 위에만 놓인다. 집단별 구성은

이고 합치면 \(8, 8, 6, 2\) 다. 동점이 대량으로 생기므로 평균순위를 쓴다. \(0.5\) 가 순위 \(1\)–\(8\) 을 차지하니 평균순위가 \(4.5\), \(1.5\) 가 \(9\)–\(16\) 이니 \(12.5\), \(2.5\) 가 \(17\)–\(22\) 이니 \(19.5\), \(3.5\) 가 \(23\)–\(24\) 이니 \(23.5\) 다. 점수도 따라서 네 개뿐이다.

평균순위 \(12.5\) 가 사분위점을 맞히는 까닭. \(N = 24\) 이므로 점수함수의 인자가

이고 \(R = 12.5\) 를 넣으면 \(0.5 + 0.25 = 0.75\) 가 정확히 된다. 그래서 그 점수는 표준정규의 상위 사분위점 \(\Phi^{-1}(0.75) = 0.674490\) 이다. \(N = 24\) 와 \(R = 12.5\) 가 맞아떨어진 우연이지만, 손으로 검산할 고정점이 하나 생긴 셈이다.

닫힌 꼴. \(\bar a_1 = \bar a_3 = A\), \(\bar a_2 = B\), \(n_i = 8\), \(N = 24\), \(k = 3\) 이라 쓰자. 전체평균이 \(\bar a = (2A+B)/3\) 이므로

이고

이다. 그러므로

로 세 수만 알면 끝난다. 여기서 \(\hat V\) 는 24 개 점수 전체의 표본분산(\(\mathrm{ddof}=1\))이다.

(2) 확인한다.

import numpy as np
from scipy import stats
from scipy.stats import fligner

# 평균은 다르지만 퍼짐은 비슷한 세 집단. 기각되지 않아야 한다.
g1 = np.array([12, 15, 14, 10, 13, 14, 12, 11], dtype=float)
g2 = np.array([22, 25, 20, 18, 24, 23, 19, 21], dtype=float)
g3 = np.array([32, 35, 34, 30, 33, 34, 32, 31], dtype=float)

X2, p = fligner(g1, g2, g3, center='median')
print(f"Fligner-Killeen X2 = {X2:.6f}, p-value = {p:.6f}")

# (1) 1단계: 중앙값 편차. g3 = g1 + 20 이므로 Z3 = Z1 이어야 한다.
Z = [np.abs(g - np.median(g)) for g in (g1, g2, g3)]
print(f"\n중앙값 = {[np.median(g) for g in (g1, g2, g3)]}")
for i, z in enumerate(Z, start=1):
    print(f"  Z{i} 정렬 = {np.sort(z)}")
print(f"  Z3 == Z1 인가: {np.array_equal(Z[0], Z[2])}")

# 2단계: 동점이 대량으로 생긴다. 서로 다른 편차값이 넷뿐이다.
allZ = np.concatenate(Z)
N = len(allZ)
vals, cnts = np.unique(allZ, return_counts=True)
R = stats.rankdata(allZ)
print(f"\nN = {N},  서로 다른 편차값 = {vals}  (개수 {cnts})")
print(f"{'편차':>7}{'개수':>6}{'차지하는 순위':>16}{'평균순위':>10}{'점수':>10}")
lo = 0
for v, c in zip(vals, cnts):
    avg = lo + (c + 1) / 2
    sc = stats.norm.ppf(0.5 + avg / (2 * (N + 1)))
    print(f"{v:>7.1f}{c:>6}{f'{lo + 1}-{lo + c}':>16}{avg:>10.1f}{sc:>10.4f}")
    lo += c
print(f"  검산: 평균순위 12.5 -> 0.5 + 12.5/50 = {0.5 + 12.5 / 50:.2f},  "
      f"Phi^-1(0.75) = {stats.norm.ppf(0.75):.6f}")

# 3단계: 집단별 점수 평균과 닫힌 꼴
a = stats.norm.ppf(0.5 + R / (2 * (N + 1)))
abar_i = np.array([a[0:8].mean(), a[8:16].mean(), a[16:24].mean()])
abar, Vhat = a.mean(), a.var(ddof=1)
print(f"\nabar_i = {np.round(abar_i, 6)}   (1 과 3 이 같은가: {abar_i[0] == abar_i[2]})")
print(f"abar = {abar:.6f},  Vhat = {Vhat:.6f}")
A, B = abar_i[0], abar_i[1]
closed = 16 * (A - B) ** 2 / (3 * Vhat)
print(f"닫힌 꼴 16(A-B)^2/(3 Vhat) = 16*({A - B:.6f})^2/(3*{Vhat:.6f}) = {closed:.6f}")
print(f"  scipy = {X2:.6f},  차 = {abs(closed - X2):.2e}")

# (2) 기각하려면 A 와 B 가 얼마나 벌어져야 하는가
crit = stats.chi2(2).ppf(0.95)
need = np.sqrt(3 * Vhat * crit / 16)
print(f"\nchi2(0.95, 2) = {crit:.4f}  ->  |A - B| > sqrt(3 Vhat * {crit:.4f} / 16) = {need:.6f}")
print(f"  실제 |A - B| = {abs(A - B):.6f}  ->  {need / abs(A - B):.3f} 배 더 벌어져야 한다")

# 집단 2 만 가진 3.5 두 점이 A-B 에 얼마나 기여하는가
print(f"\n집단 2 만 가진 편차 3.5 의 점수 = {stats.norm.ppf(0.5 + 23.5 / 50):.4f} (두 개)")
print(f"  그 두 점이 B 에 보태는 몫 = 2*{stats.norm.ppf(0.5 + 23.5 / 50):.4f}/8 = "
      f"{2 * stats.norm.ppf(0.5 + 23.5 / 50) / 8:.4f}")
print(f"  B = {B:.4f} 이므로 그 몫은 B 의 {2 * stats.norm.ppf(0.5 + 23.5 / 50) / 8 / B:.1%}")
print(f"표본분산: {[round(float(np.var(g, ddof=1)), 4) for g in (g1, g2, g3)]}")

출력:

Fligner-Killeen X2 = 2.550715, p-value = 0.279331

중앙값 = [12.5, 21.5, 32.5]
  Z1 정렬 = [0.5 0.5 0.5 1.5 1.5 1.5 2.5 2.5]
  Z2 정렬 = [0.5 0.5 1.5 1.5 2.5 2.5 3.5 3.5]
  Z3 정렬 = [0.5 0.5 0.5 1.5 1.5 1.5 2.5 2.5]
  Z3 == Z1 인가: True

N = 24,  서로 다른 편차값 = [0.5 1.5 2.5 3.5]  (개수 [8 8 6 2])
     편차    개수         차지하는 순위      평균순위        점수
    0.5     8             1-8       4.5    0.2275
    1.5     8            9-16      12.5    0.6745
    2.5     6           17-22      19.5    1.2265
    3.5     2           23-24      23.5    1.8808
  검산: 평균순위 12.5 -> 0.5 + 12.5/50 = 0.75,  Phi^-1(0.75) = 0.674490

abar_i = [0.644895 1.002339 0.644895]   (1 과 3 이 같은가: True)
abar = 0.764043,  Vhat = 0.267149
닫힌 꼴 16(A-B)^2/(3 Vhat) = 16*(-0.357444)^2/(3*0.267149) = 2.550715
  scipy = 2.550715,  차 = 4.44e-16

chi2(0.95, 2) = 5.9915  ->  |A - B| > sqrt(3 Vhat * 5.9915 / 16) = 0.547827
  실제 |A - B| = 0.357444  ->  1.533 배 더 벌어져야 한다

집단 2 만 가진 편차 3.5 의 점수 = 1.8808 (두 개)
  그 두 점이 B 에 보태는 몫 = 2*1.8808/8 = 0.4702
  B = 1.0023 이므로 그 몫은 B 의 46.9%
표본분산: [2.8393, 6.0, 2.8393]

(1)의 모든 주장이 맞는다. \(Z_3\) 과 \(Z_1\) 이 같고(\(\texttt{True}\)), 편차값이 \(0.5, 1.5, 2.5, 3.5\) 네 개이고 개수가 \(8, 8, 6, 2\) 다. 평균순위 \(4.5, 12.5, 19.5, 23.5\) 와 점수 \(0.2275, 0.6745, 1.2265, 1.8808\) 도 그대로다. \(R = 12.5\) 의 점수가 \(\Phi^{-1}(0.75) = 0.674490\) 임이 검산 줄에서 확인된다.

닫힌 꼴이 SciPy 와 기계 정밀도까지 같다. 차가 \(4.44\times10^{-16}\) 이다. 손으로 따라가면

이다. 세 집단, 24 개 관측값, 동점 열여덟 쌍짜리 계산이 결국 세 수로 접힌다.

(2) \(1.53\) 배 더 벌어져야 한다. 닫힌 꼴을 뒤집으면 기각 조건이

인데 실제는 \(0.3574\) 다. 지금 간격의 \(1.53\) 배가 필요하다. 표본분산이 \(2.84,\ 6.00,\ 2.84\) 로 두 배 넘게 차이 나는데도 그 정도가 모자란다.

왜 모자라는지는 마지막 덩어리가 말해 준다. 집단 2 를 다른 두 집단과 구별해 주는 것은 오직 편차 \(3.5\) 두 점이다. 그 둘이 받는 점수 \(1.8808\) 은 \(N = 24\) 에서 가능한 최대 점수이고, 그 두 점이 \(B = 1.0023\) 의 \(46.9\%\) 를 혼자 떠받친다. 곧 집단 2 가 "더 퍼져 있다"는 증거 전체가 관측값 두 개에 걸려 있다. 나머지 여섯 개는 다른 집단과 같은 점수를 받는다.

순위 검정에서 이것은 구조적인 제약이다. 점수에 상한 \(1.8808\) 이 있으므로, 편차 \(3.5\) 를 \(35\) 로 바꾸어도 그 두 점의 점수는 한 자리도 변하지 않는다. 분산 차이를 더 크게 만들어도 \(X^2\) 가 커지지 않는다는 뜻이고, 플리그너–킬린이 작은 표본에서 검정력을 잃는 까닭이 바로 여기다. 더 많은 관측값이 상위 순위를 차지해야 비로소 \(A\) 와 \(B\) 가 벌어진다.

집단중앙값: \(\tilde{x}_1 = 5\), \(\tilde{x}_2 = 5\).

절대편차: 집단 1은 \(|3-5|=2\), \(|7-5|=2\), \(|5-5|=0\). 집단 2는 \(|1-5|=4\), \(|10-5|=5\), \(|6-5|=1\), \(|4-5|=1\).

합친 편차를 정렬하면 \(0, 1, 1, 2, 2, 4, 5\)이고 (동점은 평균순위로) 순위는 \(R = 1, 2.5, 2.5, 4.5, 4.5, 6, 7\)이다.

정규분위수 점수 \(a = \Phi^{-1}((1 + R/8)/2)\):

집단 1 점수: \(0.7764, 0.7764, 0.1573\), 평균 \(\bar{a}_1 = 0.5701\). 집단 2 점수: \(1.1503, 1.5341, 0.4023, 0.4023\), 평균 \(\bar{a}_2 = 0.8722\). 전체평균: \(\bar{a} = (3 \times 0.5701 + 4 \times 0.8722)/7 = 0.7427\).

점수 전체의 분산은 \(\hat{V} = \frac{1}{6}\sum (a - \bar{a})^2 = 0.2282\)이므로

import numpy as np
from scipy import stats

x1 = np.array([3, 7, 5.])
x2 = np.array([1, 10, 6, 4.])
print(stats.fligner(x1, x2, center='median'))

출력:

FlignerResult(statistic=0.6859648006823972, pvalue=0.40754030224732085)

\(k-1 = 1\) 자유도에서 \(X^2 = 0.686\)이고 \(\chi^2_{0.95,1} = 3.841\)이므로 등분산을 기각하지 못한다(\(p = 0.408\)).

표본분산이 \(4\) 대 \(14\)로 3.5배 차이인데도 그렇다. \(N = 7\)에서 순위는 \(\binom{7}{3} = 35\)가지 배열밖에 구별하지 못하므로, 가장 극단적인 배열이라도 \(p\)값이 \(1/35 = 0.029\)보다 작아질 수 없다. \(\square\)

import numpy as np
from scipy.stats import bartlett, levene, fligner

rng = np.random.default_rng(42)
rej = {"bartlett": 0, "brown_forsythe": 0, "fligner": 0}
n_sims = 2000

for _ in range(n_sims):
    g1 = rng.standard_cauchy(40)
    g2 = rng.standard_cauchy(40)
    g3 = rng.standard_cauchy(40)

    _, p_b = bartlett(g1, g2, g3)
    _, p_bf = levene(g1, g2, g3, center='median')
    _, p_fk = fligner(g1, g2, g3, center='median')

    if p_b < 0.05: rej["bartlett"] += 1
    if p_bf < 0.05: rej["brown_forsythe"] += 1
    if p_fk < 0.05: rej["fligner"] += 1

for name, count in rej.items():
    print(f"{name:20s} FPR = {count/n_sims:.4f}")

출력:

bartlett             FPR = 0.9530
brown_forsythe       FPR = 0.0215
fligner              FPR = 0.0435

Bartlett 검정의 거짓 양성률이 \(0.953\)이다. 사실상 항상 기각한다. 명목값의 19배로, 이 책 전체에서 관찰한 가장 극단적인 크기 붕괴이다.

Fligner-Killeen이 \(0.044\)로 명목값에 가장 가깝다. Brown-Forsythe도 통제되지만 \(0.022\)로 명목값의 절반, 곧 지나치게 보수적이다. 편차의 평균이 코시의 극단값에 여전히 영향을 받아 집단내 변동이 부풀려지고, F 통계량의 분모가 커져 검정이 둔감해진 것이다. 순위 변환은 이 영향을 원천 차단한다.

정규분위수 점수는 두 가지 목적을 갖는다.

1. 척도의 정규화. 원래의 순위는 균등분포를 따르지 정규분포가 아니다. \(\Phi^{-1}\) 변환이 이를 근사적으로 표준정규 값으로 바꾸므로, 검정통계량에 표준 카이제곱 점근이론이 더 정확하게 적용된다.

2. 극단값 압축. 원래의 순위는 모든 관측값에 같은 간격을 준다. 극단적 이상점(순위 \(N\))과 순위 \(N-1\) 사이의 "거리"가 다른 인접 쌍과 같다. 정규점수는 꼬리를 압축한다. 최대 점수와 두 번째 점수의 차이가 내부 인접 점수들보다는 크지만, 실제 자료값이 보일 수 있는 극단성만큼 크지는 않다. 순서 정보를 보존하면서 이상점에 대한 저항을 제공한다.

원래의 순위를 쓰면. 검정은 여전히 비모수적이고 타당하겠지만 카이제곱 근사가 부정확해지고 정규 대립가설에 대한 검정력이 낮아진다.

더 근본적인 문제. 15.5절 로버스트 검정 연습문제 1에서 보았듯, 분모 \(\hat V\) 없이 원래 순위를 쓴 \(H = \sum_i n_i(\bar R_i - \bar R)^2\) 형태는 척도가 \(N^2\)에 비례하여 커지므로 카이제곱분포를 아예 따르지 못한다. 정규점수 변환은 이 척도 문제도 함께 해결한다(최대 점수가 \(\sqrt{2\ln N}\) 규모로만 커진다).

정규점수 판은 바탕 분포가 정규일 때 순위 기반 검정 가운데 점근적으로 최적이다. \(\square\)

import numpy as np
from scipy.stats import fligner

rng = np.random.default_rng(42)
port_a = rng.normal(0.01, 0.04, 60)
port_b = rng.normal(0.01, 0.06, 60)
port_c = rng.normal(0.01, 0.08, 60)

X2, p = fligner(port_a, port_b, port_c, center='median')
print(f"Fligner-Killeen: X2 = {X2:.4f}, p = {p:.6g}")
print(f"Sample SDs: A={port_a.std(ddof=1):.4f}, "
      f"B={port_b.std(ddof=1):.4f}, C={port_c.std(ddof=1):.4f}")

출력:

Fligner-Killeen: X2 = 33.2275, p = 6.09168e-08
Sample SDs: A=0.0314, B=0.0462, C=0.0816

\(p = 6.1 \times 10^{-8}\)로 등분산을 압도적으로 기각한다.

참 표준편차 4%, 6%, 8%는 의미 있게 다른 위험 수준을 나타낸다. 포트폴리오당 월별 관측값 60개면 이 차이를 탐지하기에 충분하다.

실무에서 투자자가 변동성 차이에 관심을 갖는 것은 그것이 위험조정수익률에 영향을 주기 때문이다. 통계적으로 유의한 Fligner-Killeen 결과는 이 포트폴리오들이 실제로 다른 수준의 위험을 지닌다는 것을 확인해 주며, 차별화된 자산배분 전략을 뒷받침한다.

\(H_0\) 아래에서 모든 관측값이 같은 분포에서 오므로 편차 \(z_{ij} = |x_{ij} - \tilde{x}_i|\)가 집단에 걸쳐 동일하게 분포한다. 따라서 정규분위수 점수 \(a_{ij}\)도 집단에 걸쳐 동일하게 분포한다.

\(\bar{a}_{i\cdot}\)을 집단 \(i\)의 평균 점수, \(\bar{a}_{\cdot\cdot}\)을 전체평균이라 하자. \(H_0\) 아래에서 큰 \(n\)에 대해 \(\bar{a}_{i\cdot}\)들은 근사적으로 독립이며

여기서 \(\mu_a\)와 \(\sigma_a^2\)은 점수분포의 평균과 분산이다. 검정통계량은

다변량 중심극한정리에 의해 벡터 \(\sqrt{n}(\bar{a}_{1\cdot} - \mu_a, \ldots, \bar{a}_{k\cdot} - \mu_a)\)가 \(\mathcal{N}(\mathbf{0}, \sigma_a^2 I_k)\)로 수렴한다. \(\bar{a}_{\cdot\cdot}\)을 중심으로 중심화하면 이차형식이 상수벡터에 직교하는 \((k-1)\)차원 부분공간으로의 사영을 포함하게 되어

다변량 정규벡터의 이차형식에서 계수 \(k-1\)인 멱등행렬을 쓰면 \(\chi^2(k-1)\) 분포를 갖는다는 표준 결과에서 따라 나온다.

한 가지 미묘한 점: 점수분포가 자료에 의존하지 않는다. 순위가 항상 \(\{1, \ldots, N\}\)의 순열이므로 점수 집합 \(\{a_{(1)}, \ldots, a_{(N)}\}\)은 자료와 무관하게 고정되어 있다. 바뀌는 것은 어느 점수가 어느 집단에 배정되는가뿐이다.

이 때문에 \(\hat{\sigma}_a^2\)이 사실상 알려진 상수이며, 분모의 변동에서 오는 오차가 없다. F 분포를 참조하는 Levene 계열과 달리 자유도를 추정할 필요가 없는 이유이다.

유한표본 정확성. 이 구조 덕분에 \(N\)이 작아도 근사가 비교적 정확하지만, \(\bar{a}_{i\cdot}\)의 정규근사는 여전히 필요하다. \(n_i\)가 아주 작으면(연습문제 1의 \(n_1 = 3\)처럼) 정확한 순열 \(p\)값을 쓰는 편이 낫다. \(\square\)